diff --git a/pylot/core/pick/compare.py b/pylot/core/pick/compare.py index 37534660..0adaa388 100644 --- a/pylot/core/pick/compare.py +++ b/pylot/core/pick/compare.py @@ -350,13 +350,15 @@ class PDFstatistics(object): self.stations = {} self.p_std = {} self.s_std = {} - self.makeDirlist() - self.getData() - self.getStatistics() - self.arraylen = 0 - self.Theta015 = [] - self.Theta25 = [] - self.Theta485 = [] + self.theta015 = [] + self.theta1 = [] + self.theta2 = [] + self.dirlist = list() + self.evtdict = {} + #self.makeDirlist() + #self.getData() + #self.getStatistics() + #self.showData() @@ -367,8 +369,10 @@ class PDFstatistics(object): self.evtdict[rd] = glob.glob1((os.path.join(self.directory, rd)), '*.xml') def getData(self): + arraylen = 0 for dir in self.dirlist: for evt in self.evtdict[dir]: + print evt self.stations[evt] = [] self.p_std[evt] = [] self.s_std[evt] = [] @@ -377,16 +381,24 @@ class PDFstatistics(object): # print station, pdfs try: 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: p_std = np.nan try: 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: s_std = np.nan self.stations[evt].append(station) self.p_std[evt].append(p_std) self.s_std[evt].append(s_std) - self.arraylen += 1 + arraylen += 1 + self.arraylen = arraylen + self.makeArray() def makeArray(self): @@ -426,5 +438,11 @@ class PDFstatistics(object): if __name__ == "__main__": - rootdir = '/home/sebastianp/Data/Reassessment/Insheim/' - Insheim = PDFstatistics(rootdir) \ No newline at end of file + rootdir = '/home/sebastianp/Data/Reassessment/Insheim/2012.10/' + Insheim = PDFstatistics(rootdir) + Insheim.dirlist = [rootdir] + Insheim.evtdict[rootdir] = ['e0002.294.12.xml'] + Insheim.getData() + print Insheim.theta015 + print Insheim.theta1 + print Insheim.theta2 diff --git a/pylot/core/util/pdf.py b/pylot/core/util/pdf.py index 9e1b0bee..443bf072 100644 --- a/pylot/core/util/pdf.py +++ b/pylot/core/util/pdf.py @@ -325,14 +325,15 @@ class ProbabilityDensityFunction(object): raise ValueError('value out of bounds: {0}'.format(value)) return self.prob_limits((value, self.axis[-1])) - def prob_limits(self, limits): - lim = np.arange(limits[0], limits[1], self.incr) + def prob_limits(self, limits, oversampling=1.): + sampling = self.incr / oversampling + lim = np.arange(limits[0], limits[1], sampling) data = [self.data(t) for t in lim] min_est, max_est = 0., 0. for n in range(len(data) - 1): min_est += min(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): if not (self.axis[0] <= value <= self.axis[-1]): @@ -352,7 +353,7 @@ class ProbabilityDensityFunction(object): m = (r + l) / 2 diff = prob_value - self.prob_lt_val(m) - while abs(diff) > eps: + while abs(diff) > eps and ((r - l) > self.incr): if diff > 0: l = m else: @@ -367,6 +368,13 @@ class ProbabilityDensityFunction(object): 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): import matplotlib.pyplot as plt