 Verificar datos devueltos
Como podria verificar que la linea 52 devuelve los valores correctos:
Código:
# -*- coding: utf-8 -*-
from __future__ import unicode_literals
import numpy as np
from numpy import ma
from matplotlib import scale as mscale
from matplotlib import transforms as mtransforms
from matplotlib.ticker import Formatter, FixedLocator, FormatStrFormatter
class EscalaAcumuladaGumbel(mscale.ScaleBase):
name = 'gumbel'
def __init__(self, axis, **kwargs):
mscale.ScaleBase.__init__(self)
umbral = kwargs.pop("umbral", 99.99)
if not (0.0 < umbral < 100.0):
raise ValueError("umbral debe ser mayor que 0.0 o menor que 100.0")
self.umbral = umbral
def get_transform(self):
return self.TransformacionAcumuladaGumbel(self.umbral)
def set_default_locators_and_formatters(self, axis):
class Porcentaje(Formatter):
def __call__(self, x, pos=None):
return "%d%%" % x
axis.set_major_locator(FixedLocator(np.arange(0.0, 100.0, 3.0)))
axis.set_major_formatter(Porcentaje())
axis.set_minor_formatter(Porcentaje())
def limit_range_for_scale(self, vmin, vmax, minpos):
return max(vmin, -self.umbral), min(vmax, self.umbral)
class TransformacionAcumuladaGumbel(mtransforms.Transform):
input_dims = 1
output_dims = 1
is_separable = True
def __init__(self, umbral):
mtransforms.Transform.__init__(self)
self.umbral = umbral
def transform(self, a):
masked = ma.masked_where((a < -self.umbral) | (a > self.umbral), a)
if masked.mask.any():
return 19.4243825345 + (8.74843400474*(-ma.log(-ma.log(a/100.0))))
else:
if not np.all(a):
a = np.empty(a.size)
a.fill(99.99)
return 19.4243825345 + (8.74843400474*(-np.log(-np.log(a/100.0))))
def inverted(self):
return EscalaAcumuladaGumbel.TransformacionInversaAcumuladaGumbel(self.umbral)
class TransformacionInversaAcumuladaGumbel(mtransforms.Transform):
input_dims = 1
output_dims = 1
is_separable = True
def __init__(self, umbral):
mtransforms.Transform.__init__(self)
self.umbral = umbral
def transform(self, a):
return np.exp(-np.exp(-a))
def inverted(self):
return EscalaAcumuladaGumbel.TransformacionAcumuladaGumbel(self.umbral)
mscale.register_scale(EscalaAcumuladaGumbel)
if __name__ == '__main__':
import matplotlib.pyplot as plt
Datos = np.array([ 2690.0, 2700.0, 2700.667, 2701.333, 2702.0, 3196.0, 2372.0, 2395.0, 2128.0, 2727.0, 2431.0, 2850.0, 2216.0, 2057.0, 2269.0, 2208.0, 2628.0, 2729.0, 2588.0, 3448.0, 2508.0, 3081.0, 2417.0, 2770.0, 2283.0, 2455.0, 1963.0, 2786.0, 2885.0, 2357.0, 3422.0, 2423.0, 2148.0, 1305.0, 2472.0, 2186.0, 2720.0, 2430.0, 2304.0, 2556.0, 2625.0, 2164.0, 2585.0, ])
DatosOrdenados = np.sort(Datos)
ProbabilidadAsignada = (np.array(range(1,np.size(DatosOrdenados)+ 1))*100.0)/float(np.size(DatosOrdenados)+ 1)
beta = (np.std(DatosOrdenados, ddof=1)*np.sqrt(6.0))/np.pi
mu = np.mean(DatosOrdenados) - (0.577215664901532860606*beta)
y1 = ( DatosOrdenados - mu) / beta
x1 = (np.exp(-np.exp(-y1)))*100.0
plt.plot(x1, DatosOrdenados,'-', lw=2)
plt.plot(ProbabilidadAsignada,DatosOrdenados,'ro')
plt.gca().set_xscale('gumbel')
plt.xlabel('F(z)')
plt.xticks(rotation='vertical', fontsize=7)
plt.ylabel('Precipitación')
plt.gca().yaxis.set_major_formatter(FormatStrFormatter('%d mm./año'))
plt.yticks(fontsize=7)
plt.title('Papel de probabilidad de Gumbel')
plt.grid(True)
plt.show()
La escala del eje x no parace ser la correcta, comparando con otro código;
Código:
import numpy as np
a = np.array([ 0.01, 3.0, 6.0, 9.0, 12.0, 15.0, 18.0, 21.0, 24.0, 27.0, 30.0, 33.0, 36.0, 39.0, 42.0, 45.0, 48.0, 51.0, 54.0, 57.0, 60.0, 63.0, 66.0, 69.0, 72.0, 75.0, 78.0, 81.0, 84.0, 87.0, 90.0, 93.0, 96.0, 99.99])
escala_eje_x = 19.4243825345 + (8.74843400474*(-np.log(-np.log(a/100.0))))
print escala_eje_x
Si 'a' = 3.0:
'escala_eje_x' = 8.448
Pero en el gráfico midiendo con regla 'a' = 3.0:
'escala_eje_x' = 31
También comparando con otros papeles de probabilidad las escalas son diferentes. |