undid earlier changes in PDFStatistics

This commit is contained in:
sebastianp 2016-07-13 14:15:44 +02:00
parent 6fd1f9e2cc
commit c7d8d0ecf5
2 changed files with 40 additions and 14 deletions

View File

@ -350,13 +350,15 @@ class PDFstatistics(object):
self.stations = {} self.stations = {}
self.p_std = {} self.p_std = {}
self.s_std = {} self.s_std = {}
self.makeDirlist() self.theta015 = []
self.getData() self.theta1 = []
self.getStatistics() self.theta2 = []
self.arraylen = 0 self.dirlist = list()
self.Theta015 = [] self.evtdict = {}
self.Theta25 = [] #self.makeDirlist()
self.Theta485 = [] #self.getData()
#self.getStatistics()
#self.showData() #self.showData()
@ -367,8 +369,10 @@ class PDFstatistics(object):
self.evtdict[rd] = glob.glob1((os.path.join(self.directory, rd)), '*.xml') self.evtdict[rd] = glob.glob1((os.path.join(self.directory, rd)), '*.xml')
def getData(self): def getData(self):
arraylen = 0
for dir in self.dirlist: for dir in self.dirlist:
for evt in self.evtdict[dir]: for evt in self.evtdict[dir]:
print evt
self.stations[evt] = [] self.stations[evt] = []
self.p_std[evt] = [] self.p_std[evt] = []
self.s_std[evt] = [] self.s_std[evt] = []
@ -377,16 +381,24 @@ class PDFstatistics(object):
# print station, pdfs # print station, pdfs
try: try:
p_std = pdfs['P'].standard_deviation() p_std = pdfs['P'].standard_deviation()
self.theta015.append(pdfs['P'].qtile_dist_quot(0.015))
self.theta1.append(pdfs['P'].qtile_dist_quot(0.1))
self.theta2.append(pdfs['P'].qtile_dist_quot(0.2))
except KeyError: except KeyError:
p_std = np.nan p_std = np.nan
try: try:
s_std = pdfs['S'].standard_deviation() s_std = pdfs['S'].standard_deviation()
self.theta015.append(pdfs['S'].qtile_dist_quot(0.015))
self.theta1.append(pdfs['S'].qtile_dist_quot(0.1))
self.theta2.append(pdfs['S'].qtile_dist_quot(0.2))
except KeyError: except KeyError:
s_std = np.nan s_std = np.nan
self.stations[evt].append(station) self.stations[evt].append(station)
self.p_std[evt].append(p_std) self.p_std[evt].append(p_std)
self.s_std[evt].append(s_std) self.s_std[evt].append(s_std)
self.arraylen += 1 arraylen += 1
self.arraylen = arraylen
self.makeArray()
def makeArray(self): def makeArray(self):
@ -426,5 +438,11 @@ class PDFstatistics(object):
if __name__ == "__main__": if __name__ == "__main__":
rootdir = '/home/sebastianp/Data/Reassessment/Insheim/' rootdir = '/home/sebastianp/Data/Reassessment/Insheim/2012.10/'
Insheim = PDFstatistics(rootdir) Insheim = PDFstatistics(rootdir)
Insheim.dirlist = [rootdir]
Insheim.evtdict[rootdir] = ['e0002.294.12.xml']
Insheim.getData()
print Insheim.theta015
print Insheim.theta1
print Insheim.theta2

View File

@ -325,14 +325,15 @@ class ProbabilityDensityFunction(object):
raise ValueError('value out of bounds: {0}'.format(value)) raise ValueError('value out of bounds: {0}'.format(value))
return self.prob_limits((value, self.axis[-1])) return self.prob_limits((value, self.axis[-1]))
def prob_limits(self, limits): def prob_limits(self, limits, oversampling=1.):
lim = np.arange(limits[0], limits[1], self.incr) sampling = self.incr / oversampling
lim = np.arange(limits[0], limits[1], sampling)
data = [self.data(t) for t in lim] data = [self.data(t) for t in lim]
min_est, max_est = 0., 0. min_est, max_est = 0., 0.
for n in range(len(data) - 1): for n in range(len(data) - 1):
min_est += min(data[n], data[n + 1]) min_est += min(data[n], data[n + 1])
max_est += max(data[n], data[n + 1]) max_est += max(data[n], data[n + 1])
return (min_est + max_est) / 2. * self.incr return (min_est + max_est) / 2. * sampling
def prob_val(self, value): def prob_val(self, value):
if not (self.axis[0] <= value <= self.axis[-1]): if not (self.axis[0] <= value <= self.axis[-1]):
@ -352,7 +353,7 @@ class ProbabilityDensityFunction(object):
m = (r + l) / 2 m = (r + l) / 2
diff = prob_value - self.prob_lt_val(m) diff = prob_value - self.prob_lt_val(m)
while abs(diff) > eps: while abs(diff) > eps and ((r - l) > self.incr):
if diff > 0: if diff > 0:
l = m l = m
else: else:
@ -367,6 +368,13 @@ class ProbabilityDensityFunction(object):
return qu - ql return qu - ql
def qtile_dist_quot(self,x):
if x <= 0 or x >= 0.5:
raise ValueError('Value out of range.')
return self.quantile_distance(0.5-x)/self.quantile_distance(x)
def plot(self, label=None): def plot(self, label=None):
import matplotlib.pyplot as plt import matplotlib.pyplot as plt