From 8eff4b5d0de7428f48b3d2c80119eaaa4a5f8aae Mon Sep 17 00:00:00 2001 From: Dengda98 Date: Sat, 8 Aug 2026 12:09:19 +0800 Subject: [PATCH] REFAC: rename static x/y API to north/east with xarr/yarr compatibility --- .../Advanced/integ_converg/run_ptam/run.py | 6 +- docs/source/Gallery/ex15/plot.py | 2 +- docs/source/Gallery/ex15/plot_all.py | 2 +- docs/source/Gallery/ex16/plot.py | 10 +- docs/source/Tutorial/static/run/run.py | 42 +++--- docs/source/Tutorial/static/run_upar/run.py | 12 +- .../include/grt/static/static_postprocess.h | 6 +- .../src/static/grt_static_greenfn.c | 103 ++++++------- .../src/static/grt_static_rotation.c | 44 +++--- .../src/static/grt_static_strain.c | 44 +++--- .../src/static/grt_static_stress.c | 44 +++--- pygrt/pymod.py | 84 ++++++----- pygrt/utils.py | 136 +++++++++++++----- test/_compare_c_py/compare.py | 6 +- test/_compare_c_py/compare_func.py | 4 +- test/_compare_c_py/compare_staticXY.py | 14 +- test/static_greenfn/test_static_greenfn.py | 26 ++-- test/static_syn/test_static_syn.py | 6 +- test/static_tensors/test_static_tensors.py | 6 +- 19 files changed, 341 insertions(+), 256 deletions(-) diff --git a/docs/source/Advanced/integ_converg/run_ptam/run.py b/docs/source/Advanced/integ_converg/run_ptam/run.py index 520711c5..731f8b91 100755 --- a/docs/source/Advanced/integ_converg/run_ptam/run.py +++ b/docs/source/Advanced/integ_converg/run_ptam/run.py @@ -44,9 +44,9 @@ deprcv = 0.0 pymod = pygrt.PyModel1D(modarr, depsrc=depsrc, deprcv=deprcv) -xarr = np.array([2.0]) -yarr = np.array([2.0]) -static_grn = pymod.compute_static_grn(xarr, yarr, converg_method='PTAM', statsfile=f"static_pygrtstats_{depsrc}_{deprcv}", k0=3, use_kmax_ref=True) +norths = np.array([2.0]) +easts = np.array([2.0]) +static_grn = pymod.compute_static_grn(norths, easts, converg_method='PTAM', statsfile=f"static_pygrtstats_{depsrc}_{deprcv}", k0=3, use_kmax_ref=True) ir = 0 statsdata1, statsdata2, ptamdata, dist = pygrt.utils.read_statsfile_ptam(f"static_pygrtstats_{depsrc}_{deprcv}/PTAM_{ir:04d}_*/PTAM") diff --git a/docs/source/Gallery/ex15/plot.py b/docs/source/Gallery/ex15/plot.py index 00d93cd2..4ae1ae3c 100644 --- a/docs/source/Gallery/ex15/plot.py +++ b/docs/source/Gallery/ex15/plot.py @@ -16,7 +16,7 @@ pygrnLst, _, _ = pymod._get_grn_spectra(distarr, nt=1, dt=100, zeta=1.0, keepAllFreq=True) # 静态解 -static_grn = pymod.compute_static_grn(xarr=[0.0], yarr=distarr) +static_grn = pymod.compute_static_grn(norths=[0.0], easts=distarr) # 绘制 coef = 1e-20 * 1e25 # 1e25 为地震矩 diff --git a/docs/source/Gallery/ex15/plot_all.py b/docs/source/Gallery/ex15/plot_all.py index 26e70305..4a990add 100644 --- a/docs/source/Gallery/ex15/plot_all.py +++ b/docs/source/Gallery/ex15/plot_all.py @@ -16,7 +16,7 @@ pygrnLst, _, _ = pymod._get_grn_spectra(distarr, nt=1, dt=500, zeta=1.0, keepAllFreq=True) # 静态解 -static_grn = pymod.compute_static_grn(xarr=[0.0], yarr=distarr) +static_grn = pymod.compute_static_grn(norths=[0.0], easts=distarr) # 绘制零频结果 fig, axs = plt.subplots(2, 3, figsize=(12, 8), gridspec_kw=dict(hspace=0.3, wspace=0.3)) # diff --git a/docs/source/Gallery/ex16/plot.py b/docs/source/Gallery/ex16/plot.py index 422efff9..2c0aecf1 100644 --- a/docs/source/Gallery/ex16/plot.py +++ b/docs/source/Gallery/ex16/plot.py @@ -4,7 +4,7 @@ from matplotlib.figure import Figure import sys -yarr = np.arange(0.01, 10.01, 0.01) +easts = np.arange(0.01, 10.01, 0.01) rs = np.array([10]) # 震中距数组,km nt = 801 # 总点数,不要求2的幂次 dt = 0.01 # 采样时间间隔(s) @@ -32,7 +32,7 @@ pymod1 = pygrt.PyModel1D(modarr, depsrc, deprcv, topbound=bound1, botbound=bound2) st1 = pymod1.compute_grn(distarr=rs, nt=nt, dt=dt, keepAllFreq=True)[0] pygrt.utils.stream_integral(st1) -static1 = pymod1.compute_static_grn(xarr=[0.0], yarr=yarr) +static1 = pymod1.compute_static_grn(norths=[0.0], easts=easts) # ============================================================= # 设置上下翻转模型 @@ -54,7 +54,7 @@ pymod2 = pygrt.PyModel1D(modarr2, depsrc2, deprcv2, topbound=bound2, botbound=bound1) # 整理好的模型对象 st2 = pymod2.compute_grn(distarr=rs, nt=nt, dt=dt, keepAllFreq=True)[0] -static2 = pymod2.compute_static_grn(xarr=[0.0], yarr=yarr) +static2 = pymod2.compute_static_grn(norths=[0.0], easts=easts) pygrt.utils.stream_integral(st2) @@ -87,8 +87,8 @@ ax.set_ymargin(0.3) ax = axs2[i] - ax.plot(yarr, static1[chLst[i]][0], **prop1) - ax.plot(yarr, static2[chLst[i]][0] * sgn, **prop2) + ax.plot(easts, static1[chLst[i]][0], **prop1) + ax.plot(easts, static2[chLst[i]][0] * sgn, **prop2) ax.text(0.96, 0.9, chLst[i], transform=ax.transAxes, ha='right', va='top', bbox=dict(fc='w')) ax.ticklabel_format(axis='y', style='sci', scilimits=(0,0)) diff --git a/docs/source/Tutorial/static/run/run.py b/docs/source/Tutorial/static/run/run.py index c6bf5dfa..33922045 100755 --- a/docs/source/Tutorial/static/run/run.py +++ b/docs/source/Tutorial/static/run/run.py @@ -7,14 +7,14 @@ pymod = pygrt.PyModel1D(modarr, depsrc=2.0, deprcv=0.0) -xarr = np.linspace(-3, 3, 41) -yarr = np.linspace(-2.5, 2.5, 33) +norths = np.linspace(-3, 3, 41) +easts = np.linspace(-2.5, 2.5, 33) # 可以设置 distarr 来指定震中距序列 # static_grn = pymod.compute_static_grn(distarr=np.arange(0,10+1e-8,0.1)) -# 也可以设置 xarr 和 yarr 来指定 XY 网格 -static_grn = pymod.compute_static_grn(xarr=xarr, yarr=yarr) +# 也可以设置 norths 和 easts 来指定 north/east 网格 +static_grn = pymod.compute_static_grn(norths=norths, easts=easts) print(static_grn.keys()) -# dict_keys(['_xarr', '_yarr', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'EXZ', 'VFZ', 'DDZ', 'HFZ', 'DSZ', 'SSZ', 'EXR', 'VFR', 'DDR', 'HFR', 'DSR', 'SSR', 'HFT', 'DST', 'SST']) +# dict_keys(['_norths', '_easts', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'EXZ', 'VFZ', 'DDZ', 'HFZ', 'DSZ', 'SSZ', 'EXR', 'VFR', 'DDR', 'HFR', 'DSR', 'SSR', 'HFT', 'DST', 'SST']) # END GRN # --------------------------------------------------------------------------------- @@ -27,11 +27,11 @@ def plot_static(static_syn:dict, out:Union[str,None]=None): fig, ax = plt.subplots(1, 1, figsize=(10,8)) # 设计对称色标 m = np.max(np.abs(static_syn[f'Z'])) * 1.2 - pcm = ax.pcolormesh(yarr, xarr, static_syn[f'Z'], cmap='bwr', vmin=-m, vmax=m) - ax.quiver(yarr, xarr, static_syn[f'E'], static_syn[f'N'], + pcm = ax.pcolormesh(easts, norths, static_syn[f'Z'], cmap='bwr', vmin=-m, vmax=m) + ax.quiver(easts, norths, static_syn[f'E'], static_syn[f'N'], angles='uv', pivot='mid') - ax.set_ylim([xarr[0], xarr[-1]]) - ax.set_xlim([yarr[0], yarr[-1]]) + ax.set_ylim([norths[0], norths[-1]]) + ax.set_xlim([easts[0], easts[-1]]) ax.set_aspect('equal') cbar = fig.colorbar(pcm, ax=ax, label='Z(cm)') cbar.formatter.set_powerlimits((0, 0)) @@ -46,7 +46,7 @@ def plot_static(static_syn:dict, out:Union[str,None]=None): # BEGIN SYN EX static_syn = pygrt.utils.gen_syn_from_gf_EX(static_grn, M0=1e24, ZNE=True) print(static_syn.keys()) -# dict_keys(['_xarr', '_yarr', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) +# dict_keys(['_norths', '_easts', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) plot_static(static_syn, "syn_ex.svg") # END SYN EX # --------------------------------------------------------------------------------- @@ -56,7 +56,7 @@ def plot_static(static_syn:dict, out:Union[str,None]=None): # BEGIN SYN SF static_syn = pygrt.utils.gen_syn_from_gf_SF(static_grn, S=1e16, fN=1, fE=-0.5, fZ=2, ZNE=True) print(static_syn.keys()) -# dict_keys(['_xarr', '_yarr', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) +# dict_keys(['_norths', '_easts', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) plot_static(static_syn, "syn_sf.svg") # END SYN SF # --------------------------------------------------------------------------------- @@ -66,7 +66,7 @@ def plot_static(static_syn:dict, out:Union[str,None]=None): # BEGIN SYN DC static_syn = pygrt.utils.gen_syn_from_gf_DC(static_grn, M0=1e24, strike=33, dip=50, rake=120, ZNE=True) print(static_syn.keys()) -# dict_keys(['_xarr', '_yarr', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) +# dict_keys(['_norths', '_easts', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) plot_static(static_syn, "syn_dc.svg") # END SYN DC # --------------------------------------------------------------------------------- @@ -75,7 +75,7 @@ def plot_static(static_syn:dict, out:Union[str,None]=None): # BEGIN SYN DC2 static_syn = pygrt.utils.gen_syn_from_gf_DC(static_grn, M0=1e24, strike=33, dip=90, rake=0, ZNE=True) print(static_syn.keys()) -# dict_keys(['_xarr', '_yarr', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) +# dict_keys(['_norths', '_easts', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) plot_static(static_syn, "syn_dc2.svg") # END SYN DC2 # --------------------------------------------------------------------------------- @@ -84,7 +84,7 @@ def plot_static(static_syn:dict, out:Union[str,None]=None): # BEGIN SYN TS static_syn = pygrt.utils.gen_syn_from_gf_TS(static_grn, M0=1e24, strike=33, dip=50, ZNE=True) print(static_syn.keys()) -# dict_keys(['_xarr', '_yarr', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) +# dict_keys(['_norths', '_easts', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) plot_static(static_syn, "syn_ts.svg") # END SYN TS # --------------------------------------------------------------------------------- @@ -94,7 +94,7 @@ def plot_static(static_syn:dict, out:Union[str,None]=None): # BEGIN SYN TS2 static_syn = pygrt.utils.gen_syn_from_gf_TS(static_grn, M0=1e24, strike=33, dip=90, ZNE=True) print(static_syn.keys()) -# dict_keys(['_xarr', '_yarr', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) +# dict_keys(['_norths', '_easts', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) plot_static(static_syn, "syn_ts2.svg") # END SYN TS2 # --------------------------------------------------------------------------------- @@ -104,7 +104,7 @@ def plot_static(static_syn:dict, out:Union[str,None]=None): # BEGIN SYN MT static_syn = pygrt.utils.gen_syn_from_gf_MT(static_grn, M0=1e24, MT=[0.1,-0.2,1.0,0.3,-0.5,-2.0], ZNE=True) print(static_syn.keys()) -# dict_keys(['_xarr', '_yarr', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) +# dict_keys(['_norths', '_easts', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) plot_static(static_syn, "syn_mt.svg") # END SYN MT # --------------------------------------------------------------------------------- @@ -113,7 +113,7 @@ def plot_static(static_syn:dict, out:Union[str,None]=None): # BEGIN SYN MT2 static_syn = pygrt.utils.gen_syn_from_gf_MT(static_grn, M0=1e24, MT=[0,-0.2,0,0,0,0], ZNE=True) print(static_syn.keys()) -# dict_keys(['_xarr', '_yarr', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) +# dict_keys(['_norths', '_easts', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) plot_static(static_syn, "syn_mt2.svg") # END SYN MT2 # --------------------------------------------------------------------------------- @@ -121,11 +121,11 @@ def plot_static(static_syn:dict, out:Union[str,None]=None): # --------------------------------------------------------------------------------- # BEGIN NEW XY -xarr2 = np.arange(-3, 3+1e-8, 0.2) -yarr2 = np.arange(-2.5, 2.5+1e-8, 0.25) -static_syn = pygrt.utils.gen_syn_from_gf_DC(static_grn, M0=1e24, strike=33, dip=90, rake=0, ZNE=True, xarr=xarr2, yarr=yarr2) +norths2 = np.arange(-3, 3+1e-8, 0.2) +easts2 = np.arange(-2.5, 2.5+1e-8, 0.25) +static_syn = pygrt.utils.gen_syn_from_gf_DC(static_grn, M0=1e24, strike=33, dip=90, rake=0, ZNE=True, norths=norths2, easts=easts2) print(static_syn.keys()) -# dict_keys(['_xarr', '_yarr', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) +# dict_keys(['_norths', '_easts', '_src_va', '_src_vb', '_src_rho', '_rcv_va', '_rcv_vb', '_rcv_rho', 'Z', 'N', 'E']) plot_static(static_syn, "synXY_dc2.svg") # END NEW XY # --------------------------------------------------------------------------------- diff --git a/docs/source/Tutorial/static/run_upar/run.py b/docs/source/Tutorial/static/run_upar/run.py index 90d94332..296e969a 100755 --- a/docs/source/Tutorial/static/run_upar/run.py +++ b/docs/source/Tutorial/static/run_upar/run.py @@ -5,8 +5,8 @@ def plot6(data:dict, title:str, out:str|None=None): chs = [k for k in data.keys() if k[0]!='_'] chs.sort(reverse=True) - xarr = data['_xarr'] - yarr = data['_yarr'] + norths = data['_norths'] + easts = data['_easts'] fig, axs = plt.subplots(len(chs)//3, 3, figsize=(10, len(chs))) axs = axs.ravel() @@ -25,7 +25,7 @@ def plot6(data:dict, title:str, out:str|None=None): vmin = -1 vmax = 1 - pcm = ax.pcolormesh(yarr, xarr, data[ch], shading='nearest', vmin=vmin, vmax=vmax, rasterized=True) + pcm = ax.pcolormesh(easts, norths, data[ch], shading='nearest', vmin=vmin, vmax=vmax, rasterized=True) ax.set_aspect('equal') ax.set_title(ch) cbar = fig.colorbar(pcm, ax=ax) @@ -42,10 +42,10 @@ def plot6(data:dict, title:str, out:str|None=None): pymod = pygrt.PyModel1D(modarr, depsrc=2.0, deprcv=0.0) -xarr = np.linspace(-3, 3, 41) -yarr = np.linspace(-2.5, 2.5, 33) +norths = np.linspace(-3, 3, 41) +easts = np.linspace(-2.5, 2.5, 33) # 传入calc_upar=True可计算空间导数 -static_grn = pymod.compute_static_grn(xarr, yarr, calc_upar=True) +static_grn = pymod.compute_static_grn(norths, easts, calc_upar=True) # 传入calc_upar=True可计算空间导数 # 传入ZNE=True返回ZNE分量 diff --git a/pygrt/C_extension/include/grt/static/static_postprocess.h b/pygrt/C_extension/include/grt/static/static_postprocess.h index c54b857b..68a2ed46 100644 --- a/pygrt/C_extension/include/grt/static/static_postprocess.h +++ b/pygrt/C_extension/include/grt/static/static_postprocess.h @@ -18,7 +18,7 @@ * 与 syn 中轴点处 (1/r)∂_θ 有限部分配套。 */ void grt_static_compute_stress( - size_t nx, size_t ny, const real_t *xs, const real_t *ys, + size_t nnorth, size_t neast, const real_t *norths, const real_t *easts, real_t *const u[GRT_CHANNEL_NUM], real_t *const upar[GRT_CHANNEL_NUM][GRT_CHANNEL_NUM], real_t *const res[GRT_CHANNEL_NUM][GRT_CHANNEL_NUM], @@ -29,7 +29,7 @@ void grt_static_compute_stress( * 数组布局同 grt_static_compute_stress()。 */ void grt_static_compute_strain( - size_t nx, size_t ny, const real_t *xs, const real_t *ys, + size_t nnorth, size_t neast, const real_t *norths, const real_t *easts, real_t *const u[GRT_CHANNEL_NUM], real_t *const upar[GRT_CHANNEL_NUM][GRT_CHANNEL_NUM], real_t *const res[GRT_CHANNEL_NUM][GRT_CHANNEL_NUM], bool rot2ZNE); @@ -39,7 +39,7 @@ void grt_static_compute_strain( * 数组布局同 grt_static_compute_stress()。仅写入 res 的非对角分量。 */ void grt_static_compute_rotation( - size_t nx, size_t ny, const real_t *xs, const real_t *ys, + size_t nnorth, size_t neast, const real_t *norths, const real_t *easts, real_t *const u[GRT_CHANNEL_NUM], real_t *const upar[GRT_CHANNEL_NUM][GRT_CHANNEL_NUM], real_t *const res[GRT_CHANNEL_NUM][GRT_CHANNEL_NUM], bool rot2ZNE); diff --git a/pygrt/C_extension/src/static/grt_static_greenfn.c b/pygrt/C_extension/src/static/grt_static_greenfn.c index 0b967213..a8d1659f 100644 --- a/pygrt/C_extension/src/static/grt_static_greenfn.c +++ b/pygrt/C_extension/src/static/grt_static_greenfn.c @@ -68,17 +68,17 @@ typedef struct { bool active; char *s_statsdir; ///< 保存目录,和当前目录同级 } S; - /** X 坐标 */ + /** -X: north 坐标 */ struct { bool active; - size_t nx; - real_t *xs; + size_t nnorth; + real_t *norths; } X; - /** Y 坐标 */ + /** -Y: east 坐标 */ struct { bool active; - size_t ny; - real_t *ys; + size_t neast; + real_t *easts; } Y; /** 输出 nc 文件名 */ struct { @@ -106,10 +106,10 @@ static void free_Ctrl(GRT_MODULE_CTRL *Ctrl){ GRT_SAFE_FREE_PTR(Ctrl->D.s_deprcv); // X - GRT_SAFE_FREE_PTR(Ctrl->X.xs); + GRT_SAFE_FREE_PTR(Ctrl->X.norths); // Y - GRT_SAFE_FREE_PTR(Ctrl->Y.ys); + GRT_SAFE_FREE_PTR(Ctrl->Y.easts); // O GRT_SAFE_FREE_PTR(Ctrl->O.s_outgrid); @@ -145,8 +145,11 @@ printf("\n" " [-K[+k][+f][+e]] [-S] [-e]\n" "\n" " There're two ways to define the \"epicentral distances\":\n" -" 1. set both -X and -Y to define a XY grid in advance.\n" -" 2. simply set -R, which equal to set the Y-coord with X=0.0.\n" +" 1. set both -X and -Y (north/east). The kernel still depends\n" +" only on r=hypot(north,east); results are stored with north/east\n" +" dims for compatibility.\n" +" 2. set -R for a 1D distance list. Stored as nnorth=1, norths=[0],\n" +" easts=R (recommended for building a GF library).\n" "\n\n" "Options:\n" "----------------------------------------------------------------\n" @@ -442,10 +445,10 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ GRTBadOptionError(X, "x1(%f) > x2(%f).", a1, a2); } - Ctrl->X.nx = floor((a2-a1)/delta) + 1; - Ctrl->X.xs = (real_t*)calloc(Ctrl->X.nx, sizeof(real_t)); - for(size_t i=0; iX.nx; ++i){ - Ctrl->X.xs[i] = a1 + delta*i; + Ctrl->X.nnorth = floor((a2-a1)/delta) + 1; + Ctrl->X.norths = (real_t*)calloc(Ctrl->X.nnorth, sizeof(real_t)); + for(size_t i=0; iX.nnorth; ++i){ + Ctrl->X.norths[i] = a1 + delta*i; } } break; @@ -465,25 +468,25 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ GRTBadOptionError(Y, "y1(%f) > y2(%f).", a1, a2); } - Ctrl->Y.ny = floor((a2-a1)/delta) + 1; - Ctrl->Y.ys = (real_t*)calloc(Ctrl->Y.ny, sizeof(real_t)); - for(size_t i=0; iY.ny; ++i){ - Ctrl->Y.ys[i] = a1 + delta*i; + Ctrl->Y.neast = floor((a2-a1)/delta) + 1; + Ctrl->Y.easts = (real_t*)calloc(Ctrl->Y.neast, sizeof(real_t)); + for(size_t i=0; iY.neast; ++i){ + Ctrl->Y.easts[i] = a1 + delta*i; } } break; - // -R 算是别名,相当于 -X0/0/1, -Yy1/y2/dy, 以此方式指定一维震中距序列 + // -R:一维震中距序列,存成 nnorth=1, norths=[0], easts=R // -R,[,...]|//| case 'R': Ctrl->X.active = Ctrl->Y.active = true; { real_t a1, a2, delta; - char **s_ys = NULL; + char **s_easts = NULL; // 如果输入仅由数字、小数点和间隔符组成,则直接读取 if(grt_string_composed_of(optarg, GRT_NUM_STR "eE+-" ".,")){ - s_ys = grt_string_split(optarg, ",", &Ctrl->Y.ny); + s_easts = grt_string_split(optarg, ",", &Ctrl->Y.neast); } // 尝试按照 // 读取 else if(3 == sscanf(optarg, "%lf/%lf/%lf", &a1, &a2, &delta)){ @@ -494,32 +497,32 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ GRTBadOptionError(R, "r1(%f) > r2(%f).", a1, a2); } - Ctrl->Y.ny = floor((a2-a1)/delta) + 1; - s_ys = (char **)calloc(Ctrl->Y.ny, sizeof(char*) * Ctrl->Y.ny); - for(size_t ir = 0; ir < Ctrl->Y.ny; ++ir){ - GRT_SAFE_ASPRINTF(&s_ys[ir], "%.*f", 8, a1 + delta*ir); + Ctrl->Y.neast = floor((a2-a1)/delta) + 1; + s_easts = (char **)calloc(Ctrl->Y.neast, sizeof(char*) * Ctrl->Y.neast); + for(size_t ir = 0; ir < Ctrl->Y.neast; ++ir){ + GRT_SAFE_ASPRINTF(&s_easts[ir], "%.*f", 8, a1 + delta*ir); } } // 否则从文件读取 else { FILE *fp = GRTCheckOpenFile(optarg, "r"); - s_ys = grt_string_from_file(fp, &Ctrl->Y.ny); + s_easts = grt_string_from_file(fp, &Ctrl->Y.neast); fclose(fp); } // 转为浮点数 - Ctrl->Y.ys = (real_t*)realloc(Ctrl->Y.ys, sizeof(real_t)*(Ctrl->Y.ny)); - for(size_t i=0; i < Ctrl->Y.ny; ++i){ - Ctrl->Y.ys[i] = atof(s_ys[i]); - if(Ctrl->Y.ys[i] < 0.0){ - GRTBadOptionError(R, "Can't set negative epicentral distance(%f).", Ctrl->Y.ys[i]); + Ctrl->Y.easts = (real_t*)realloc(Ctrl->Y.easts, sizeof(real_t)*(Ctrl->Y.neast)); + for(size_t i=0; i < Ctrl->Y.neast; ++i){ + Ctrl->Y.easts[i] = atof(s_easts[i]); + if(Ctrl->Y.easts[i] < 0.0){ + GRTBadOptionError(R, "Can't set negative epicentral distance(%f).", Ctrl->Y.easts[i]); } } - GRT_SAFE_FREE_PTR_ARRAY(s_ys, Ctrl->Y.ny); + GRT_SAFE_FREE_PTR_ARRAY(s_easts, Ctrl->Y.neast); - Ctrl->X.nx = 1; - Ctrl->X.xs = (real_t*)calloc(Ctrl->X.nx, sizeof(real_t)); - Ctrl->X.xs[0] = 0.0; + Ctrl->X.nnorth = 1; + Ctrl->X.norths = (real_t*)calloc(Ctrl->X.nnorth, sizeof(real_t)); + Ctrl->X.norths[0] = 0.0; } break; @@ -552,11 +555,11 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ GRTCheckOptionActive(Ctrl, O); // 设置震中距数组 - Ctrl->nr = Ctrl->X.nx*Ctrl->Y.ny; + Ctrl->nr = Ctrl->X.nnorth*Ctrl->Y.neast; Ctrl->rs = (real_t*)calloc(Ctrl->nr, sizeof(real_t)); - for(size_t ix=0; ixX.nx; ++ix){ - for(size_t iy=0; iyY.ny; ++iy){ - Ctrl->rs[iy + ix*Ctrl->Y.ny] = hypot(Ctrl->X.xs[ix], Ctrl->Y.ys[iy]); + for(size_t inorth=0; inorthX.nnorth; ++inorth){ + for(size_t ieast=0; ieastY.neast; ++ieast){ + Ctrl->rs[ieast + inorth*Ctrl->Y.neast] = hypot(Ctrl->X.norths[inorth], Ctrl->Y.easts[ieast]); } } @@ -677,10 +680,10 @@ int static_greenfn_main(int argc, char **argv){ // ================================================================================== // 将结果保存为 nc 格式 // ================================================================================== - int ncid, x_dimid, y_dimid; + int ncid, north_dimid, east_dimid; const int ndims = 2; int dimids[ndims]; - int x_varid, y_varid; + int north_varid, east_varid; intChnlGrid u_varids; intChnlGrid uiz_varids; intChnlGrid uir_varids; @@ -689,6 +692,8 @@ int static_greenfn_main(int argc, char **argv){ NC_CHECK(nc_create(Ctrl->O.s_outgrid, NC_CLOBBER, &ncid)); // 写入全局属性 + NC_CHECK(NC_FUNC_REAL(nc_put_att) (ncid, NC_GLOBAL, "depsrc", NC_REAL, 1, &mod1d->depsrc)); + NC_CHECK(NC_FUNC_REAL(nc_put_att) (ncid, NC_GLOBAL, "deprcv", NC_REAL, 1, &mod1d->deprcv)); NC_CHECK(NC_FUNC_REAL(nc_put_att) (ncid, NC_GLOBAL, "src_va", NC_REAL, 1, &src_va)); NC_CHECK(NC_FUNC_REAL(nc_put_att) (ncid, NC_GLOBAL, "src_vb", NC_REAL, 1, &src_vb)); NC_CHECK(NC_FUNC_REAL(nc_put_att) (ncid, NC_GLOBAL, "src_rho", NC_REAL, 1, &src_rho)); @@ -702,14 +707,14 @@ int static_greenfn_main(int argc, char **argv){ } // 定义维度 - NC_CHECK(nc_def_dim(ncid, "north", Ctrl->X.nx, &x_dimid)); - NC_CHECK(nc_def_dim(ncid, "east", Ctrl->Y.ny, &y_dimid)); - dimids[0] = x_dimid; - dimids[1] = y_dimid; + NC_CHECK(nc_def_dim(ncid, "north", Ctrl->X.nnorth, &north_dimid)); + NC_CHECK(nc_def_dim(ncid, "east", Ctrl->Y.neast, &east_dimid)); + dimids[0] = north_dimid; + dimids[1] = east_dimid; // 定义维度数组 - NC_CHECK(nc_def_var(ncid, "north", NC_REAL, 1, &x_dimid, &x_varid)); - NC_CHECK(nc_def_var(ncid, "east", NC_REAL, 1, &y_dimid, &y_varid)); + NC_CHECK(nc_def_var(ncid, "north", NC_REAL, 1, &north_dimid, &north_varid)); + NC_CHECK(nc_def_var(ncid, "east", NC_REAL, 1, &east_dimid, &east_varid)); // 定义不同震源不同分量的格林函数数组 GRT_LOOP_ChnlGrid(im, c){ @@ -735,8 +740,8 @@ int static_greenfn_main(int argc, char **argv){ NC_CHECK(nc_enddef(ncid)); // 写入数据 - NC_CHECK(NC_FUNC_REAL(nc_put_var) (ncid, x_varid, Ctrl->X.xs)); - NC_CHECK(NC_FUNC_REAL(nc_put_var) (ncid, y_varid, Ctrl->Y.ys)); + NC_CHECK(NC_FUNC_REAL(nc_put_var) (ncid, north_varid, Ctrl->X.norths)); + NC_CHECK(NC_FUNC_REAL(nc_put_var) (ncid, east_varid, Ctrl->Y.easts)); real_t *tmpdata = (real_t *)calloc(Ctrl->nr, sizeof(real_t)); GRT_LOOP_ChnlGrid(im, c){ int modr = GRT_SRC_M_ORDERS[im]; diff --git a/pygrt/C_extension/src/static/grt_static_rotation.c b/pygrt/C_extension/src/static/grt_static_rotation.c index 8d01a498..0485753a 100644 --- a/pygrt/C_extension/src/static/grt_static_rotation.c +++ b/pygrt/C_extension/src/static/grt_static_rotation.c @@ -50,17 +50,17 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ } void grt_static_compute_rotation( - size_t nx, size_t ny, const real_t *xs, const real_t *ys, + size_t nnorth, size_t neast, const real_t *norths, const real_t *easts, real_t *const u[GRT_CHANNEL_NUM], real_t *const upar[GRT_CHANNEL_NUM][GRT_CHANNEL_NUM], real_t *const res[GRT_CHANNEL_NUM][GRT_CHANNEL_NUM], bool rot2ZNE) { const char *chs = rot2ZNE ? GRT_ZNE_CODES : GRT_ZRT_CODES; - for(size_t ix=0; ix Tuple[object, object]: + """解析 north/east 坐标,兼容旧公开 API 名 xarr/yarr""" + if xarr is not None or yarr is not None: + _warn_xarr_yarr_deprecated(stacklevel=stacklevel) + if norths is not None or easts is not None: + raise ValueError("Use either norths/easts or xarr/yarr, not both.") + norths, easts = xarr, yarr + if norths is None: + norths = default_norths + if easts is None: + easts = default_easts + return norths, easts + + +def _pop_static_ne_from_kwargs( + kwargs:dict, + default_norths=None, + default_easts=None, + stacklevel:int=3, +) -> Tuple[object, object]: + """从 kwargs 取出 norths/easts(或旧名 xarr/yarr)""" + return _resolve_static_ne_coords( + kwargs.pop('norths', None), kwargs.pop('easts', None), + kwargs.pop('xarr', None), kwargs.pop('yarr', None), + default_norths=default_norths, default_easts=default_easts, + stacklevel=stacklevel, + ) + + +def _static_dict_norths_easts(d:dict, stacklevel:int=3) -> Tuple[object, object]: + """从静态结果字典读取 north/east,兼容旧键 _xarr/_yarr""" + if '_norths' in d and '_easts' in d: + return d['_norths'], d['_easts'] + if '_xarr' in d and '_yarr' in d: + _warn_xarr_yarr_deprecated(stacklevel=stacklevel) + return d['_xarr'], d['_yarr'] + raise KeyError("static dict missing '_norths'/'_easts' (or deprecated '_xarr'/'_yarr').") + from enum import Enum, unique @unique @@ -200,13 +253,16 @@ def _gen_syn_from_static_gf(grnDct:dict, calc_upar:bool, compute_type:GRT_SYN_TY # 为张裂计算 Vp/Vs VpVs_ratio = float(grnDct['_src_va'] / grnDct['_src_vb']) - xarr0 = np.ascontiguousarray(grnDct['_xarr'], dtype=NPCT_REAL_TYPE) - yarr0 = np.ascontiguousarray(grnDct['_yarr'], dtype=NPCT_REAL_TYPE) - xarr = np.ascontiguousarray(kwargs.get('xarr', grnDct['_xarr']), dtype=NPCT_REAL_TYPE) - yarr = np.ascontiguousarray(kwargs.get('yarr', grnDct['_yarr']), dtype=NPCT_REAL_TYPE) - nx0, ny0 = len(xarr0), len(yarr0) - nx, ny = len(xarr), len(yarr) - nr0, nr = nx0 * ny0, nx * ny + norths0, easts0 = _static_dict_norths_easts(grnDct, stacklevel=4) + norths0 = np.ascontiguousarray(norths0, dtype=NPCT_REAL_TYPE) + easts0 = np.ascontiguousarray(easts0, dtype=NPCT_REAL_TYPE) + norths, easts = _pop_static_ne_from_kwargs( + kwargs, default_norths=norths0, default_easts=easts0, stacklevel=4) + norths = np.ascontiguousarray(norths, dtype=NPCT_REAL_TYPE) + easts = np.ascontiguousarray(easts, dtype=NPCT_REAL_TYPE) + nnorth0, neast0 = len(norths0), len(easts0) + nnorth, neast = len(norths), len(easts) + nr0, nr = nnorth0 * neast0, nnorth * neast def _pack_gf(prefix:str=''): """打包为 C 侧 realChnlGrid[nr]:arr[震中距点][震源][分量]""" @@ -229,8 +285,8 @@ def _pack_gf(prefix:str=''): # 调用 C 函数 mchn = _set_source_mechanism(compute_type, **kwargs) C_grt_static_syn_from_gf( - nx0, npct.as_ctypes(xarr0), ny0, npct.as_ctypes(yarr0), - nx, npct.as_ctypes(xarr), ny, npct.as_ctypes(yarr), + nnorth0, npct.as_ctypes(norths0), neast0, npct.as_ctypes(easts0), + nnorth, npct.as_ctypes(norths), neast, npct.as_ctypes(easts), npct.as_ctypes(pygrn), npct.as_ctypes(pygrn_uiz) if calc_upar else None, npct.as_ctypes(pygrn_uir) if calc_upar else None, @@ -241,15 +297,15 @@ def _pack_gf(prefix:str=''): # ======================================================================== resDct = {k: deepcopy(v) for k, v in grnDct.items() if k.startswith('_')} - resDct['_xarr'] = xarr - resDct['_yarr'] = yarr + resDct['_norths'] = norths + resDct['_easts'] = easts chs = ZNEchs if ZNE else ZRTchs for i1, c1 in enumerate(chs): - resDct[c1] = syn[:, i1].reshape((nx, ny)) + resDct[c1] = syn[:, i1].reshape((nnorth, neast)) if calc_upar: for i2, c2 in enumerate(chs): - resDct[f'{c2.lower()}{c1}'] = syn_upar[:, i2, i1].reshape((nx, ny)) + resDct[f'{c2.lower()}{c1}'] = syn_upar[:, i2, i1].reshape((nnorth, neast)) return resDct @@ -292,8 +348,9 @@ def gen_syn_from_gf_DC(st:Union[Stream,dict], M0:float, strike:float, dip:float, :param az: azimuth, 0 <= az <= 360 (not used for static case) :param ZNE: whether output in 'ZNE'-coord, default is 'ZRT' :param calc_upar: whether calculate the spatial derivatives of displacements. - :param kwargs: For static rsults, you can set "xarr" and "yarr" to define a new XY grid, - the program will automatically find the nearest Green's functions for each node. + :param kwargs: For static results, set ``norths``/``easts`` (preferred) or + deprecated ``xarr``/``yarr`` to define a new north/east grid; + synthesis interpolates in epicentral distance. :return: - **stream** - :class:`obspy.Stream` @@ -318,8 +375,9 @@ def gen_syn_from_gf_TS(st:Union[Stream,dict], M0:float, strike:float, dip:float, :param az: azimuth, 0 <= az <= 360 (not used for static case) :param ZNE: whether output in 'ZNE'-coord, default is 'ZRT' :param calc_upar: whether calculate the spatial derivatives of displacements. - :param kwargs: For static rsults, you can set "xarr" and "yarr" to define a new XY grid, - the program will automatically find the nearest Green's functions for each node. + :param kwargs: For static results, set ``norths``/``easts`` (preferred) or + deprecated ``xarr``/``yarr`` to define a new north/east grid; + synthesis interpolates in epicentral distance. :return: - **stream** - :class:`obspy.Stream` @@ -346,8 +404,9 @@ def gen_syn_from_gf_SF(st:Union[Stream,dict], S:float, fN:float, fE:float, fZ:fl :param az: azimuth, 0 <= az <= 360 (not used for static case) :param ZNE: whether output in 'ZNE'-coord, default is 'ZRT' :param calc_upar: whether calculate the spatial derivatives of displacements. - :param kwargs: For static rsults, you can set "xarr" and "yarr" to define a new XY grid, - the program will automatically find the nearest Green's functions for each node. + :param kwargs: For static results, set ``norths``/``easts`` (preferred) or + deprecated ``xarr``/``yarr`` to define a new north/east grid; + synthesis interpolates in epicentral distance. :return: - **stream** - :class:`obspy.Stream` @@ -371,8 +430,9 @@ def gen_syn_from_gf_EX(st:Union[Stream,dict], M0:float, az:float=-999, ZNE=False :param az: azimuth, 0 <= az <= 360 (not used for static case) :param ZNE: whether output in 'ZNE'-coord, default is 'ZRT' :param calc_upar: whether calculate the spatial derivatives of displacements. - :param kwargs: For static rsults, you can set "xarr" and "yarr" to define a new XY grid, - the program will automatically find the nearest Green's functions for each node. + :param kwargs: For static results, set ``norths``/``easts`` (preferred) or + deprecated ``xarr``/``yarr`` to define a new north/east grid; + synthesis interpolates in epicentral distance. :return: - **stream** - :class:`obspy.Stream` @@ -397,8 +457,9 @@ def gen_syn_from_gf_MT(st:Union[Stream,dict], M0:float, MT:ArrayLike, az:float=- :param az: azimuth, 0 <= az <= 360 (not used for static case) :param ZNE: whether output in 'ZNE'-coord, default is 'ZRT' :param calc_upar: whether calculate the spatial derivatives of displacements. - :param kwargs: For static rsults, you can set "xarr" and "yarr" to define a new XY grid, - the program will automatically find the nearest Green's functions for each node. + :param kwargs: For static results, set ``norths``/``easts`` (preferred) or + deprecated ``xarr``/``yarr`` to define a new north/east grid; + synthesis interpolates in epicentral distance. :return: - **stream** - :class:`obspy.Stream` @@ -506,20 +567,21 @@ def _prepare_static_postprocess_arrays(syn:dict, chs:List[str]): RPtrs = PREAL * CHANNEL_NUM RMat = RPtrs * CHANNEL_NUM - xarr = np.ascontiguousarray(syn['_xarr'], dtype=np.float64) - yarr = np.ascontiguousarray(syn['_yarr'], dtype=np.float64) - if xarr.ndim != 1 or yarr.ndim != 1: - raise ValueError("'_xarr' and '_yarr' must be one-dimensional arrays.") + norths, easts = _static_dict_norths_easts(syn, stacklevel=4) + norths = np.ascontiguousarray(norths, dtype=np.float64) + easts = np.ascontiguousarray(easts, dtype=np.float64) + if norths.ndim != 1 or easts.ndim != 1: + raise ValueError("'_norths' and '_easts' must be one-dimensional arrays.") u = [np.ascontiguousarray(syn[c], dtype=np.float64) for c in chs] upar = [ [np.ascontiguousarray(syn[f"{d.lower()}{c}"], dtype=np.float64) for c in chs] for d in chs ] - expected_shape = (len(xarr), len(yarr)) + expected_shape = (len(norths), len(easts)) if any(arr.shape != expected_shape for arr in u) or any( arr.shape != expected_shape for row in upar for arr in row): - raise ValueError("Static displacement and derivative arrays must match '_xarr'/'_yarr'.") + raise ValueError("Static displacement and derivative arrays must match '_norths'/'_easts'.") u_ptrs = RPtrs(*(arr.ctypes.data_as(PREAL) for arr in u)) upar_ptrs = RMat(*( @@ -528,7 +590,7 @@ def _prepare_static_postprocess_arrays(syn:dict, chs:List[str]): res_ptrs = RMat(*( RPtrs(*(resarr[c2, c1].ctypes.data_as(PREAL) for c1 in range(CHANNEL_NUM))) for c2 in range(CHANNEL_NUM))) - return xarr, yarr, u, upar, u_ptrs, upar_ptrs, resarr, res_ptrs + return norths, easts, u, upar, u_ptrs, upar_ptrs, resarr, res_ptrs def _compute_static_strain_rotation(syn:dict, Type:str): @@ -557,7 +619,7 @@ def _compute_static_strain_rotation(syn:dict, Type:str): if f"nN" in syn.keys(): chs = ZNEchs - xarr, yarr, u, upar, u_ptrs, upar_ptrs, resarr, res_ptrs = \ + norths, easts, u, upar, u_ptrs, upar_ptrs, resarr, res_ptrs = \ _prepare_static_postprocess_arrays(syn, chs) # 结果字典 @@ -571,11 +633,11 @@ def _compute_static_strain_rotation(syn:dict, Type:str): if Type == 'strain': C_grt_static_compute_strain( - len(xarr), len(yarr), xarr.ctypes.data_as(PREAL), yarr.ctypes.data_as(PREAL), + len(norths), len(easts), norths.ctypes.data_as(PREAL), easts.ctypes.data_as(PREAL), u_ptrs, upar_ptrs, res_ptrs, chs == ZNEchs) else: C_grt_static_compute_rotation( - len(xarr), len(yarr), xarr.ctypes.data_as(PREAL), yarr.ctypes.data_as(PREAL), + len(norths), len(easts), norths.ctypes.data_as(PREAL), easts.ctypes.data_as(PREAL), u_ptrs, upar_ptrs, res_ptrs, chs == ZNEchs) for i1 in range(i1_end): @@ -687,7 +749,7 @@ def _compute_static_stress(syn:dict): chs = ZNEchs rot2ZNE = True - xarr, yarr, u, upar, u_ptrs, upar_ptrs, resarr, res_ptrs = \ + norths, easts, u, upar, u_ptrs, upar_ptrs, resarr, res_ptrs = \ _prepare_static_postprocess_arrays(syn, chs) va = syn['_rcv_va'] vb = syn['_rcv_vb'] @@ -705,7 +767,7 @@ def _compute_static_stress(syn:dict): resDct[k] = deepcopy(syn[k]) C_grt_static_compute_stress( - len(xarr), len(yarr), xarr.ctypes.data_as(PREAL), yarr.ctypes.data_as(PREAL), + len(norths), len(easts), norths.ctypes.data_as(PREAL), easts.ctypes.data_as(PREAL), u_ptrs, upar_ptrs, res_ptrs, rot2ZNE, mu, lam) for i1 in range(CHANNEL_NUM): diff --git a/test/_compare_c_py/compare.py b/test/_compare_c_py/compare.py index af74b842..ff24f8e1 100644 --- a/test/_compare_c_py/compare.py +++ b/test/_compare_c_py/compare.py @@ -103,10 +103,10 @@ #-------------------------- Static ----------------------------------------- # 为了方便测试,避免引入其他因素的误差,这里有意避开 0 -xarr = np.arange(-3.1, 3.2, 0.6) -yarr = np.arange(-4.1, 4.2, 0.8) +norths = np.arange(-3.1, 3.2, 0.6) +easts = np.arange(-4.1, 4.2, 0.8) -static_grn = pymod.compute_static_grn(xarr, yarr, calc_upar=True) +static_grn = pymod.compute_static_grn(norths, easts, calc_upar=True) AVGRERR2 = [] # plot_static(static_grn, "static/stgrn.nc") diff --git a/test/_compare_c_py/compare_func.py b/test/_compare_c_py/compare_func.py index fbef374e..894ca306 100644 --- a/test/_compare_c_py/compare_func.py +++ b/test/_compare_c_py/compare_func.py @@ -91,8 +91,8 @@ def plot_static(resDct:dict, c_prefix:str): n = len([k for k in resDct.keys() if k[0] != '_']) fig, axs = plt.subplots(n, 3, figsize=(8, 3*n)) with netcdf_file(c_prefix, mmap=False) as f: - norths = resDct['_xarr'] - easts = resDct['_yarr'] + norths = resDct['_norths'] + easts = resDct['_easts'] keys = f.variables keys.pop('north') diff --git a/test/_compare_c_py/compare_staticXY.py b/test/_compare_c_py/compare_staticXY.py index 9904d544..67cdd39b 100644 --- a/test/_compare_c_py/compare_staticXY.py +++ b/test/_compare_c_py/compare_staticXY.py @@ -14,8 +14,8 @@ #-------------------------- Static ----------------------------------------- # 为了方便测试,避免引入其他因素的误差,这里有意避开 0 -xarr = np.arange(-3.1, 3.2, 0.6) -yarr = np.arange(-4.1, 4.2, 0.8) +norths = np.arange(-3.1, 3.2, 0.6) +easts = np.arange(-4.1, 4.2, 0.8) S=1e24 fn=2 @@ -39,7 +39,7 @@ for ZNE in [False, True]: suffix = "-N" if ZNE else "" - static_syn = pygrt.utils.gen_syn_from_gf_EX(static_grn, S, ZNE=ZNE, calc_upar=True, xarr=xarr, yarr=yarr) + static_syn = pygrt.utils.gen_syn_from_gf_EX(static_grn, S, ZNE=ZNE, calc_upar=True, norths=norths, easts=easts) ststrain = pygrt.utils.compute_strain(static_syn) strotation = pygrt.utils.compute_rotation(static_syn) ststress = pygrt.utils.compute_stress(static_syn) @@ -48,7 +48,7 @@ update_dict(static_syn, strotation, "rotation_") AVGRERR2.append(static_compare3(static_syn, f"static/stsyn_ex{suffix}.nc")) - static_syn = pygrt.utils.gen_syn_from_gf_SF(static_grn, S, fn, fe, fz, ZNE=ZNE, calc_upar=True, xarr=xarr, yarr=yarr) + static_syn = pygrt.utils.gen_syn_from_gf_SF(static_grn, S, fn, fe, fz, ZNE=ZNE, calc_upar=True, norths=norths, easts=easts) ststrain = pygrt.utils.compute_strain(static_syn) strotation = pygrt.utils.compute_rotation(static_syn) ststress = pygrt.utils.compute_stress(static_syn) @@ -57,7 +57,7 @@ update_dict(static_syn, strotation, "rotation_") AVGRERR2.append(static_compare3(static_syn, f"static/stsyn_sf{suffix}.nc")) - static_syn = pygrt.utils.gen_syn_from_gf_DC(static_grn, S, stk, dip, rak, ZNE=ZNE, calc_upar=True, xarr=xarr, yarr=yarr) + static_syn = pygrt.utils.gen_syn_from_gf_DC(static_grn, S, stk, dip, rak, ZNE=ZNE, calc_upar=True, norths=norths, easts=easts) ststrain = pygrt.utils.compute_strain(static_syn) strotation = pygrt.utils.compute_rotation(static_syn) ststress = pygrt.utils.compute_stress(static_syn) @@ -66,7 +66,7 @@ update_dict(static_syn, strotation, "rotation_") AVGRERR2.append(static_compare3(static_syn, f"static/stsyn_dc{suffix}.nc")) - static_syn = pygrt.utils.gen_syn_from_gf_TS(static_grn, S, stk, dip, ZNE=ZNE, calc_upar=True, xarr=xarr, yarr=yarr) + static_syn = pygrt.utils.gen_syn_from_gf_TS(static_grn, S, stk, dip, ZNE=ZNE, calc_upar=True, norths=norths, easts=easts) ststrain = pygrt.utils.compute_strain(static_syn) strotation = pygrt.utils.compute_rotation(static_syn) ststress = pygrt.utils.compute_stress(static_syn) @@ -75,7 +75,7 @@ update_dict(static_syn, strotation, "rotation_") AVGRERR2.append(static_compare3(static_syn, f"static/stsyn_ts{suffix}.nc")) - static_syn = pygrt.utils.gen_syn_from_gf_MT(static_grn, S, [M11,M12,M13,M22,M23,M33], ZNE=ZNE, calc_upar=True, xarr=xarr, yarr=yarr) + static_syn = pygrt.utils.gen_syn_from_gf_MT(static_grn, S, [M11,M12,M13,M22,M23,M33], ZNE=ZNE, calc_upar=True, norths=norths, easts=easts) ststrain = pygrt.utils.compute_strain(static_syn) strotation = pygrt.utils.compute_rotation(static_syn) ststress = pygrt.utils.compute_stress(static_syn) diff --git a/test/static_greenfn/test_static_greenfn.py b/test/static_greenfn/test_static_greenfn.py index 99f8952e..af6c1e7a 100644 --- a/test/static_greenfn/test_static_greenfn.py +++ b/test/static_greenfn/test_static_greenfn.py @@ -3,29 +3,29 @@ depsrc=2 deprcv=0 -xarr = np.arange(-3., 3.1, 0.2) -yarr = np.arange(-2., 2.1, 0.2) +norths = np.arange(-3., 3.1, 0.2) +easts = np.arange(-2., 2.1, 0.2) modname="../milrow" modarr = np.loadtxt(modname) pymod = pygrt.PyModel1D(modarr, depsrc, deprcv) -stgrn = pymod.compute_static_grn(xarr, yarr) -stgrn = pymod.compute_static_grn(xarr, yarr, calc_upar=True) -stgrn = pymod.compute_static_grn(xarr, yarr, Length=20) +stgrn = pymod.compute_static_grn(norths, easts) +stgrn = pymod.compute_static_grn(norths, easts, calc_upar=True) +stgrn = pymod.compute_static_grn(norths, easts, Length=20) -stgrn = pymod.compute_static_grn(xarr, yarr, Length=20, converg_method='DCM') -stgrn = pymod.compute_static_grn(xarr, yarr, Length=20, converg_method='PTAM') -stgrn = pymod.compute_static_grn(xarr, yarr, Length=20, converg_method='none') +stgrn = pymod.compute_static_grn(norths, easts, Length=20, converg_method='DCM') +stgrn = pymod.compute_static_grn(norths, easts, Length=20, converg_method='PTAM') +stgrn = pymod.compute_static_grn(norths, easts, Length=20, converg_method='none') -stgrn = pymod.compute_static_grn(xarr, yarr, k0=4, keps=1e-3) -stgrn = pymod.compute_static_grn(xarr, yarr, k0=4, statsfile="stgrt_stats") +stgrn = pymod.compute_static_grn(norths, easts, k0=4, keps=1e-3) +stgrn = pymod.compute_static_grn(norths, easts, k0=4, statsfile="stgrt_stats") # boundary condition pymod = pygrt.PyModel1D(modarr, depsrc, deprcv, topbound='free', botbound='free') -stgrn = pymod.compute_static_grn(xarr, yarr) +stgrn = pymod.compute_static_grn(norths, easts) pymod = pygrt.PyModel1D(modarr, depsrc, deprcv, topbound='halfspace', botbound='free') -stgrn = pymod.compute_static_grn(xarr, yarr) +stgrn = pymod.compute_static_grn(norths, easts) pymod = pygrt.PyModel1D(modarr, depsrc, deprcv, topbound='rigid', botbound='rigid') -stgrn = pymod.compute_static_grn(xarr, yarr) \ No newline at end of file +stgrn = pymod.compute_static_grn(norths, easts) \ No newline at end of file diff --git a/test/static_syn/test_static_syn.py b/test/static_syn/test_static_syn.py index e5d6aef6..4c0e78e5 100644 --- a/test/static_syn/test_static_syn.py +++ b/test/static_syn/test_static_syn.py @@ -3,15 +3,15 @@ depsrc=2 deprcv=0 -xarr = np.arange(-3., 3.1, 0.2) -yarr = np.arange(-2., 2.1, 0.2) +norths = np.arange(-3., 3.1, 0.2) +easts = np.arange(-2., 2.1, 0.2) modname="../milrow" modarr = np.loadtxt(modname) pymod = pygrt.PyModel1D(modarr, depsrc, deprcv) -stgrn = pymod.compute_static_grn(xarr, yarr, calc_upar=True) +stgrn = pymod.compute_static_grn(norths, easts, calc_upar=True) stsyn = pygrt.utils.gen_syn_from_gf_EX(stgrn, 1e20, 22) stsyn = pygrt.utils.gen_syn_from_gf_SF(stgrn, 1e16, fN=-1, fE=2, fZ=-4, az=22) diff --git a/test/static_tensors/test_static_tensors.py b/test/static_tensors/test_static_tensors.py index d12856ba..78fd333c 100644 --- a/test/static_tensors/test_static_tensors.py +++ b/test/static_tensors/test_static_tensors.py @@ -3,15 +3,15 @@ depsrc=2 deprcv=0 -xarr = np.arange(-3., 3.1, 0.2) -yarr = np.arange(-2., 2.1, 0.2) +norths = np.arange(-3., 3.1, 0.2) +easts = np.arange(-2., 2.1, 0.2) modname="../milrow" modarr = np.loadtxt(modname) pymod = pygrt.PyModel1D(modarr, depsrc, deprcv) -stgrn = pymod.compute_static_grn(xarr, yarr, calc_upar=True) +stgrn = pymod.compute_static_grn(norths, easts, calc_upar=True) stsyn = pygrt.utils.gen_syn_from_gf_EX(stgrn, 1e20, 22, calc_upar=True)