summaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
authorAaron LI <aaronly.me@outlook.com>2016-06-25 16:14:35 +0800
committerAaron LI <aaronly.me@outlook.com>2016-06-25 16:14:35 +0800
commite856955776a2fc072068cf44be7b1b5f999e1eed (patch)
tree5c7c1cdc5b720c14048afbe993525d7e4649446e
parentdf786c667ebb052155aa9b99b22d556ff42e5d50 (diff)
downloadcexcess-e856955776a2fc072068cf44be7b1b5f999e1eed.tar.bz2
calc_pei.py: Use InterpolatedUnivariateSpline
-rwxr-xr-xcalc_pei.py15
1 files changed, 8 insertions, 7 deletions
diff --git a/calc_pei.py b/calc_pei.py
index 961c19c..bb5879d 100755
--- a/calc_pei.py
+++ b/calc_pei.py
@@ -2,9 +2,11 @@
#
# Aaron LI
# Created: 2016-04-29
-# Updated: 2016-05-18
+# Updated: 2016-06-25
#
# Change log:
+# 2016-06-25:
+# * Use 'InterpolatedUnivariateSpline' instead of 'interp1d'
# 2016-05-18:
# * Roughly implement the PEI uncertainty estimation
# * Fix/Update PEI Y positions determination
@@ -34,8 +36,8 @@ import json
from collections import OrderedDict
import numpy as np
-import scipy.interpolate
-import scipy.integrate
+import scipy.interpolate as interpolate
+import scipy.integrate as integrate
import matplotlib.pyplot as plt
from matplotlib.backends.backend_agg import FigureCanvasAgg as FigureCanvas
@@ -85,9 +87,8 @@ def calc_pei(data, r500, interp_np=101, pei_pos=None):
# data points within the PEI range
mask_pei = np.logical_and(x >= pei_xmin, x <= pei_xmax)
y_pei = y[mask_pei]
- # interpolate the power spectrum
- f_interp = scipy.interpolate.interp1d(x, y, kind="linear",
- assume_sorted=True)
+ # interpolate the power spectrum by fitting a smoothing spline
+ f_interp = interpolate.InterpolatedUnivariateSpline(x, y)
x_interp = np.linspace(pei_xmin, pei_xmax, num=interp_np)
y_interp = f_interp(x_interp)
# determine the Y positions of PEI rectangle
@@ -111,7 +112,7 @@ def calc_pei(data, r500, interp_np=101, pei_pos=None):
y_interp[y_interp < pei_ymin] = pei_ymin
# calculate the PEI
area_total = (pei_xmax - pei_xmin) * (pei_ymax - pei_ymin)
- area_below = scipy.integrate.trapz((y_interp-pei_ymin), x_interp)
+ area_below = integrate.trapz((y_interp-pei_ymin), x_interp)
pei_value = area_below / area_total
results = {
"pei_xmin": pei_xmin,