From d5d7521eb27565ae48a9e1cee51edd0d2e986490 Mon Sep 17 00:00:00 2001 From: Dengda98 Date: Wed, 31 Dec 2025 14:48:54 +0800 Subject: [PATCH 01/12] FEAT: support tensile source --- .../include/grt/common/checkerror.h | 1 + .../include/grt/common/radiation.h | 32 ++- pygrt/C_extension/include/grt/dynamic/syn.h | 5 +- pygrt/C_extension/src/common/radiation.c | 47 ++++- pygrt/C_extension/src/dynamic/grt_syn.c | 91 +++++---- pygrt/C_extension/src/dynamic/syn.c | 11 +- pygrt/C_extension/src/static/grt_static_syn.c | 75 ++++--- pygrt/c_interfaces.py | 9 + pygrt/c_structures.py | 3 +- pygrt/utils.py | 191 +++++++----------- 10 files changed, 266 insertions(+), 199 deletions(-) diff --git a/pygrt/C_extension/include/grt/common/checkerror.h b/pygrt/C_extension/include/grt/common/checkerror.h index db7681f8..30574fed 100644 --- a/pygrt/C_extension/include/grt/common/checkerror.h +++ b/pygrt/C_extension/include/grt/common/checkerror.h @@ -9,6 +9,7 @@ #pragma once +#include #include #include #include diff --git a/pygrt/C_extension/include/grt/common/radiation.h b/pygrt/C_extension/include/grt/common/radiation.h index b46a6a3e..947360d5 100644 --- a/pygrt/C_extension/include/grt/common/radiation.h +++ b/pygrt/C_extension/include/grt/common/radiation.h @@ -13,20 +13,38 @@ #include "grt/common/const.h" -#define GRT_SYN_COMPUTE_EX 0 ///< 计算爆炸源 -#define GRT_SYN_COMPUTE_SF 1 ///< 计算单力源 -#define GRT_SYN_COMPUTE_DC 2 ///< 计算剪切源 -#define GRT_SYN_COMPUTE_MT 3 ///< 计算矩张量源 +#define __FOR_ALL_SYN_TYPE \ + X(GRT_SYN_EX, "Explosion") \ + X(GRT_SYN_SF, "Single Force") \ + X(GRT_SYN_DC, "Shear") \ + X(GRT_SYN_TS, "Tensile") \ + X(GRT_SYN_MT, "Moment Tensor") \ +// 要合成的震源类型 +typedef enum { + #define X(f, name) f, + __FOR_ALL_SYN_TYPE + #undef X +} GRT_SYN_TYPE; + +// 对应的震源类型全称 +static const char *srcTypeFullName[] = { + #define X(f, name) name, + __FOR_ALL_SYN_TYPE + #undef X +}; +#undef __FOR_ALL_SYN_TYPE + /** * 设置每个震源的方向因子 * * @param[out] srcRadi 方向因子，[3]表示ZRT三分量，[6]表示6个震源(EX,VF,HF,DD,DS,SS) - * @param[in] computeType 要计算的震源类型，使用宏定义 + * @param[in] computeType 要计算的震源类型 * @param[in] par_theta 方向因子中是否对theta(az)求导 * @param[in] M0 放大系数，对于剪切源、爆炸源、张量震源，M0是标量地震矩；对于单力源，M0是放大系数 * @param[in] coef 放大系数，用于位移空间导数的计算 + * @param[in] VpVs_ratio 震源层的 Vp/Vs 比值，用于张位错源 * @param[in] azrad 弧度制的方位角 * @param[in] mchn 震源机制参数， * 对于单力源，mchn={fn, fe, fz}， @@ -34,6 +52,6 @@ * 对于张量源，mchn={Mxx, Mxy, Mxz, Myy, Myz, Mzz} */ void grt_set_source_radiation( - realChnlGrid srcRadi, const int computeType, const bool par_theta, - const real_t M0, const real_t coef, const real_t azrad, const real_t mchn[GRT_MECHANISM_NUM] + realChnlGrid srcRadi, const GRT_SYN_TYPE computeType, const bool par_theta, + const real_t M0, const real_t coef, const real_t VpVs_ratio, const real_t azrad, const real_t mchn[GRT_MECHANISM_NUM] ); \ No newline at end of file diff --git a/pygrt/C_extension/include/grt/dynamic/syn.h b/pygrt/C_extension/include/grt/dynamic/syn.h index 54fb46d9..340fde79 100644 --- a/pygrt/C_extension/include/grt/dynamic/syn.h +++ b/pygrt/C_extension/include/grt/dynamic/syn.h @@ -11,14 +11,15 @@ #include "grt/common/const.h" #include "grt/common/sacio2.h" +#include "grt/common/radiation.h" /** * 根据已有系数（方向因子）合成理论地震图 * * @param[in] srcRadi 不同震源不同分量的系数 - * @param[in] computeType 要计算的震源类型，使用宏定义 + * @param[in] computeType 要计算的震源类型 * @param[in] dirpath 格林函数所在目录路径 * @param[in] prefix 格林函数文件名前缀 * @param[out] synsac 三分量 SACTRACE */ -void grt_syn(const realChnlGrid srcRadi, const int computeType, const char *dirpath, const char *prefix, SACTRACE *synsac[GRT_CHANNEL_NUM]); \ No newline at end of file +void grt_syn(const realChnlGrid srcRadi, const GRT_SYN_TYPE computeType, const char *dirpath, const char *prefix, SACTRACE *synsac[GRT_CHANNEL_NUM]); \ No newline at end of file diff --git a/pygrt/C_extension/src/common/radiation.c b/pygrt/C_extension/src/common/radiation.c index 538f953a..b11788c5 100644 --- a/pygrt/C_extension/src/common/radiation.c +++ b/pygrt/C_extension/src/common/radiation.c @@ -11,13 +11,14 @@ #include "grt/common/radiation.h" #include "grt/common/const.h" +#include "grt/common/checkerror.h" void grt_set_source_radiation( - realChnlGrid srcRadi, const int computeType, const bool par_theta, - const real_t M0, const real_t coef, const real_t azrad, const real_t mchn[GRT_MECHANISM_NUM] + realChnlGrid srcRadi, const GRT_SYN_TYPE computeType, const bool par_theta, + const real_t M0, const real_t coef, const real_t VpVs_ratio, const real_t azrad, const real_t mchn[GRT_MECHANISM_NUM] ){ real_t mult; - if(computeType == GRT_SYN_COMPUTE_SF){ + if(computeType == GRT_SYN_SF){ mult = 1e-15*M0*coef; } else { mult = 1e-20*M0*coef; @@ -27,11 +28,11 @@ void grt_set_source_radiation( saz = sin(azrad); caz = cos(azrad); - if(computeType == GRT_SYN_COMPUTE_EX){ + if(computeType == GRT_SYN_EX){ srcRadi[0][0] = srcRadi[0][1] = (par_theta)? 0.0 : mult; // Z/R srcRadi[0][2] = 0.0; // T } - else if(computeType == GRT_SYN_COMPUTE_SF){ + else if(computeType == GRT_SYN_SF){ real_t A0, A1, A4; real_t fn, fe, fz; fn=mchn[0]; fe=mchn[1]; fz=mchn[2]; @@ -45,7 +46,7 @@ void grt_set_source_radiation( srcRadi[1][2] = 0.0; // VF, T srcRadi[2][2] = (par_theta)? -A1 : A4; // HF, T } - else if(computeType == GRT_SYN_COMPUTE_DC){ + else if(computeType == GRT_SYN_DC){ real_t strike, dip, rake; strike=mchn[0]; dip=mchn[1]; rake=mchn[2]; // 公式(4.8.35) @@ -74,7 +75,36 @@ void grt_set_source_radiation( srcRadi[4][2] = (par_theta)? -A1 : A4; // DS, T srcRadi[5][2] = (par_theta)? -2.0*A2 : A5; // DS, T } - else if(computeType == GRT_SYN_COMPUTE_MT){ + else if(computeType == GRT_SYN_TS){ + real_t strike, dip; + strike=mchn[0]; dip=mchn[1]; + // 公式(4.8.35) + real_t stkrad = strike*DEG1; + real_t diprad = dip*DEG1; + real_t therad = azrad - stkrad; + real_t sdip, cdip, sdip2, sthe, cthe, sthe2, cthe2; + sdip = sin(diprad); cdip = cos(diprad); + sdip2 = 2.0*sdip*cdip; + sthe = sin(therad); cthe = cos(therad); + sthe2 = 2.0*sthe*cthe; cthe2 = 2.0*cthe*cthe - 1.0; + + real_t A0, A1, A2, A4, A5; + A0 = mult * (2.0/3.0 - sdip*sdip); + A1 = mult * sdip2*sthe; + A2 = mult * ( - sdip*sdip*cthe2); + A4 = mult * sdip2*cthe; + A5 = mult * sdip*sdip*sthe2; + + srcRadi[0][0] = srcRadi[0][1] = (par_theta)? 0.0 : mult * (GRT_SQUARE(VpVs_ratio) - 4.0/3.0); // EX, Z/R + srcRadi[3][0] = srcRadi[3][1] = (par_theta)? 0.0 : A0; // DD, Z/R + srcRadi[4][0] = srcRadi[4][1] = (par_theta)? A4 : A1; // DS, Z/R + srcRadi[5][0] = srcRadi[5][1] = (par_theta)? 2.0*A5 : A2; // SS, Z/R + srcRadi[0][2] = 0.0; // EX, T + srcRadi[3][2] = 0.0; // DD, T + srcRadi[4][2] = (par_theta)? -A1 : A4; // DS, T + srcRadi[5][2] = (par_theta)? -2.0*A2 : A5; // DS, T + } + else if(computeType == GRT_SYN_MT){ // 公式(4.9.7)但修改了各向同性的量 real_t M11, M12, M13, M22, M23, M33; M11 = mchn[0]; M12 = mchn[1]; M13 = mchn[2]; @@ -105,4 +135,7 @@ void grt_set_source_radiation( srcRadi[4][2] = (par_theta)? -A1 : A4; // DS, T srcRadi[5][2] = (par_theta)? -2.0*A2 : A5; // DS, T } + else { + GRTRaiseError("Unsupported source type."); + } } \ No newline at end of file diff --git a/pygrt/C_extension/src/dynamic/grt_syn.c b/pygrt/C_extension/src/dynamic/grt_syn.c index 8768b493..9523f8f2 100644 --- a/pygrt/C_extension/src/dynamic/grt_syn.c +++ b/pygrt/C_extension/src/dynamic/grt_syn.c @@ -21,9 +21,6 @@ // 防止被替换为虚数单位 #undef I -// 和宏命令对应的震源类型全称 -static const char *sourceTypeFullName[] = {"Explosion", "Single Force", "Shear", "Moment Tensor"}; - /** 该子模块的参数控制结构体 */ typedef struct { /** 格林函数路径 */ @@ -97,11 +94,14 @@ typedef struct { // 震中距 real_t dist; + // 震源层 Vp/Vs + real_t VpVs_ratio; + // 方向因子数组 realChnlGrid srcRadi; // 最终要计算的震源类型 - int computeType; + GRT_SYN_TYPE computeType; char s_computeType[3]; } GRT_MODULE_CTRL; @@ -137,7 +137,7 @@ printf("\n" "Usage:\n" "----------------------------------------------------------------\n" " grt syn -G -A -S[u] -O \n" -" [-M//]\n" +" [-M/[/]]\n" " [-T/////]\n" " [-F//] \n" " [-D/] [-I] [-J]\n" @@ -162,9 +162,10 @@ printf("\n" " For source type, you can only set at most one of\n" " '-M', '-T' and '-F'. If none, an Explosion is used.\n" "\n" -" -M//\n" -" Three angles to define a fault. \n" +" -M/[/]\n" +" Three angles to define a shear fault. \n" " The angles are in degree.\n" +" If not set, then define a tensile fault.\n" "\n" " -T/////\n" " Six elements of Moment Tensor. \n" @@ -234,6 +235,9 @@ printf("\n" " or Shear\n" " grt syn -Gres/milrow_2_0_10 -Osyn_dc -A30 -S1e24 -M100/20/80\n" "\n" +" or Tension\n" +" grt syn -Gres/milrow_2_0_10 -Osyn_dc -A30 -S1e24 -M100/20\n" +"\n" " or Single Force\n" " grt syn -Gres/milrow_2_0_10 -Osyn_sf -A30 -S1e24 -F0.5/-1.2/3.3\n" "\n" @@ -247,7 +251,7 @@ printf("\n" /** 从命令行中读取选项，处理后记录到全局变量中 */ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ // 先为个别参数设置非0初始值 - Ctrl->computeType = GRT_SYN_COMPUTE_EX; + Ctrl->computeType = GRT_SYN_EX; sprintf(Ctrl->s_computeType, "%s", "EX"); int opt; @@ -291,25 +295,33 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ }; break; - // 剪切震源 + // 剪切震源， 张位错源 case 'M': Ctrl->M.active = true; - Ctrl->computeType = GRT_SYN_COMPUTE_DC; { - real_t strike, dip, rake; - sprintf(Ctrl->s_computeType, "%s", "DC"); - if(3 != sscanf(optarg, "%lf/%lf/%lf", &strike, &dip, &rake)){ + real_t strike=0.0, dip=0.0, rake=0.0; + int nscan = sscanf(optarg, "%lf/%lf/%lf", &strike, &dip, &rake); + if(nscan >= 2){ + Ctrl->computeType = GRT_SYN_TS; + sprintf(Ctrl->s_computeType, "%s", "TS"); + if(strike < 0.0 || strike > 360.0){ + GRTBadOptionError(M, "Strike must be in [0, 360]."); + } + if(dip < 0.0 || dip > 90.0){ + GRTBadOptionError(M, "Dip must be in [0, 90]."); + } + if(nscan == 3){ + Ctrl->computeType = GRT_SYN_DC; + sprintf(Ctrl->s_computeType, "%s", "DC"); + if(rake < -180.0 || rake > 180.0){ + GRTBadOptionError(M, "Rake must be in [-180, 180]."); + } + } + } else { GRTBadOptionError(M, ""); }; - if(strike < 0.0 || strike > 360.0){ - GRTBadOptionError(M, "Strike must be in [0, 360]."); - } - if(dip < 0.0 || dip > 90.0){ - GRTBadOptionError(M, "Dip must be in [0, 90]."); - } - if(rake < -180.0 || rake > 180.0){ - GRTBadOptionError(M, "Rake must be in [-180, 180]."); - } + + Ctrl->mchn[0] = strike; Ctrl->mchn[1] = dip; Ctrl->mchn[2] = rake; @@ -319,7 +331,7 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ // 单力源 case 'F': Ctrl->F.active = true; - Ctrl->computeType = GRT_SYN_COMPUTE_SF; + Ctrl->computeType = GRT_SYN_SF; { real_t fn, fe, fz; sprintf(Ctrl->s_computeType, "%s", "SF"); @@ -335,7 +347,7 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ // 张量震源 case 'T': Ctrl->T.active = true; - Ctrl->computeType = GRT_SYN_COMPUTE_MT; + Ctrl->computeType = GRT_SYN_MT; { real_t Mxx, Mxy, Mxz, Myy, Myz, Mzz; sprintf(Ctrl->s_computeType, "%s", "MT"); @@ -442,19 +454,22 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ Ctrl->dist = sac->hd.dist; + float va, vb, rho; + va = sac->hd.user6; + vb = sac->hd.user7; + rho = sac->hd.user8; + if(va <= 0.0 || vb < 0.0 || rho <= 0.0){ + GRTRaiseError("Bad src_va, src_vb or src_rho in \"%s\" header.\n", entry->d_name); + } + if(vb == 0.0){ + GRTRaiseError("Zero src_vb in \"%s\" header. " + "Maybe you try to use -Su but the source is in the liquid. " + "Use -S instead.\n" , entry->d_name); + } + + Ctrl->VpVs_ratio = ((real_t)va)/vb; + if (Ctrl->S.mult_src_mu) { - float va, vb, rho; - va = sac->hd.user6; - vb = sac->hd.user7; - rho = sac->hd.user8; - if(va <= 0.0 || vb < 0.0 || rho <= 0.0){ - GRTRaiseError("Bad src_va, src_vb or src_rho in \"%s\" header.\n", entry->d_name); - } - if(vb == 0.0){ - GRTRaiseError("Zero src_vb in \"%s\" header. " - "Maybe you try to use -Su but the source is in the liquid. " - "Use -S instead.\n" , entry->d_name); - } Ctrl->S.src_mu = vb*vb*rho*1e10; Ctrl->S.M0 *= Ctrl->S.src_mu; } @@ -572,7 +587,7 @@ int syn_main(int argc, char **argv){ } // 重新计算方向因子 - grt_set_source_radiation(Ctrl->srcRadi, Ctrl->computeType, (ityp==3), Ctrl->S.M0, upar_scale, Ctrl->A.azrad, Ctrl->mchn); + grt_set_source_radiation(Ctrl->srcRadi, Ctrl->computeType, (ityp==3), Ctrl->S.M0, upar_scale, Ctrl->VpVs_ratio, Ctrl->A.azrad, Ctrl->mchn); // 合成地震图 if (ityp==0 || ityp==3) { @@ -662,7 +677,7 @@ int syn_main(int argc, char **argv){ if(! Ctrl->s.active) { GRTRaiseInfo("Under \"%s\"", Ctrl->O.s_output_dir); - GRTRaiseInfo("Synthetic Seismograms of %-13s source done.", sourceTypeFullName[Ctrl->computeType]); + GRTRaiseInfo("Synthetic Seismograms of %-13s source done.", srcTypeFullName[Ctrl->computeType]); if(tfsac != NULL) GRTRaiseInfo("Time Function saved."); } diff --git a/pygrt/C_extension/src/dynamic/syn.c b/pygrt/C_extension/src/dynamic/syn.c index ad7af15d..af0aabde 100644 --- a/pygrt/C_extension/src/dynamic/syn.c +++ b/pygrt/C_extension/src/dynamic/syn.c @@ -10,21 +10,20 @@ #include #include #include "grt/dynamic/syn.h" -#include "grt/common/radiation.h" -void grt_syn(const realChnlGrid srcRadi, const int computeType, const char *dirpath, const char *prefix, SACTRACE *synsac[GRT_CHANNEL_NUM]) +void grt_syn(const realChnlGrid srcRadi, const GRT_SYN_TYPE computeType, const char *dirpath, const char *prefix, SACTRACE *synsac[GRT_CHANNEL_NUM]) { GRT_LOOP_ChnlGrid(im, c) { int modr = GRT_SRC_M_ORDERS[im]; if(modr == 0 && GRT_QWV_CODES[c] == 'v') continue; - if (computeType == GRT_SYN_COMPUTE_EX) { + if (computeType == GRT_SYN_EX) { if (im != GRT_SRC_M_EX_INDEX) continue; - } else if (computeType == GRT_SYN_COMPUTE_SF) { + } else if (computeType == GRT_SYN_SF) { if (im != GRT_SRC_M_VF_INDEX && im != GRT_SRC_M_HF_INDEX) continue; - } else if (computeType == GRT_SYN_COMPUTE_DC) { + } else if (computeType == GRT_SYN_DC) { if (im < GRT_SRC_M_DD_INDEX) continue; - } else if (computeType == GRT_SYN_COMPUTE_MT) { + } else if (computeType == GRT_SYN_TS || computeType == GRT_SYN_MT) { if (im < GRT_SRC_M_DD_INDEX && im != GRT_SRC_M_EX_INDEX) continue; } else { GRTRaiseError("Not Supported."); diff --git a/pygrt/C_extension/src/static/grt_static_syn.c b/pygrt/C_extension/src/static/grt_static_syn.c index 0c43e5c6..5b39b466 100644 --- a/pygrt/C_extension/src/static/grt_static_syn.c +++ b/pygrt/C_extension/src/static/grt_static_syn.c @@ -51,6 +51,10 @@ typedef struct { bool active; char *s_outgrid; } O; + /** 静默输出 */ + struct { + bool active; + } s; /** 是否计算空间导数 */ struct { bool active; @@ -63,7 +67,7 @@ typedef struct { realChnlGrid srcRadi; // 最终要计算的震源类型 - int computeType; + GRT_SYN_TYPE computeType; char s_computeType[3]; } GRT_MODULE_CTRL; @@ -94,10 +98,10 @@ printf("\n" "Usage:\n" "----------------------------------------------------------------\n" " grt static syn -G -S[u] \n" -" [-M//]\n" +" [-M/[/]]\n" " [-T/////]\n" " [-F//] \n" -" [-N] [-e]\n" +" [-N] [-e] [-s]\n" " -O \n" "\n" "\n\n" @@ -119,9 +123,10 @@ printf("\n" " For source type, you can only set at most one of\n" " '-M', '-T' and '-F'. If none, an Explosion is used.\n" "\n" -" -M//\n" -" Three angles to define a fault. \n" +" -M/[/]\n" +" Three angles to define a shear fault. \n" " The angles are in degree.\n" +" If not set, then define a tensile fault.\n" "\n" " -T/////\n" " Six elements of Moment Tensor. \n" @@ -138,6 +143,8 @@ printf("\n" " of displacement u. In filenames, prefix \"r\" means \n" " ui_r and \"z\" means ui_z. \n" "\n" +" -s Silence all outputs.\n" +"\n" " -h Display this help message.\n" "\n\n" "Examples:\n" @@ -151,6 +158,9 @@ printf("\n" " or Shear\n" " grt static syn -Gstgrn.nc -Su1e16 -M100/20/80 -Ostsyn_dc.nc\n" "\n" +" or Tension\n" +" grt static syn -Gstgrn.nc -Su1e16 -M100/20 -Ostsyn_dc.nc\n" +"\n" " or Single Force\n" " grt static syn -Gstgrn.nc -S1e20 -F0.5/-1.2/3.3 -Ostsyn_sf.nc\n" "\n" @@ -165,11 +175,11 @@ printf("\n" /** 从命令行中读取选项，处理后记录到全局变量中 */ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ // 先为个别参数设置非0初始值 - Ctrl->computeType = GRT_SYN_COMPUTE_EX; + Ctrl->computeType = GRT_SYN_EX; sprintf(Ctrl->s_computeType, "%s", "EX"); int opt; - while ((opt = getopt(argc, argv, ":G:O:S:M:F:T:Neh")) != -1) { + while ((opt = getopt(argc, argv, ":G:O:S:M:F:T:Nesh")) != -1) { switch (opt) { // 输入 nc 文件名 case 'G': @@ -199,25 +209,33 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ }; break; - // 剪切震源 + // 剪切震源， 张位错源 case 'M': Ctrl->M.active = true; - Ctrl->computeType = GRT_SYN_COMPUTE_DC; { - real_t strike, dip, rake; - sprintf(Ctrl->s_computeType, "%s", "DC"); - if(3 != sscanf(optarg, "%lf/%lf/%lf", &strike, &dip, &rake)){ + real_t strike=0.0, dip=0.0, rake=0.0; + int nscan = sscanf(optarg, "%lf/%lf/%lf", &strike, &dip, &rake); + if(nscan >= 2){ + Ctrl->computeType = GRT_SYN_TS; + sprintf(Ctrl->s_computeType, "%s", "TS"); + if(strike < 0.0 || strike > 360.0){ + GRTBadOptionError(M, "Strike must be in [0, 360]."); + } + if(dip < 0.0 || dip > 90.0){ + GRTBadOptionError(M, "Dip must be in [0, 90]."); + } + if(nscan == 3){ + Ctrl->computeType = GRT_SYN_DC; + sprintf(Ctrl->s_computeType, "%s", "DC"); + if(rake < -180.0 || rake > 180.0){ + GRTBadOptionError(M, "Rake must be in [-180, 180]."); + } + } + } else { GRTBadOptionError(M, ""); }; - if(strike < 0.0 || strike > 360.0){ - GRTBadOptionError(M, "Strike must be in [0, 360]."); - } - if(dip < 0.0 || dip > 90.0){ - GRTBadOptionError(M, "Dip must be in [0, 90]."); - } - if(rake < -180.0 || rake > 180.0){ - GRTBadOptionError(M, "Rake must be in [-180, 180]."); - } + + Ctrl->mchn[0] = strike; Ctrl->mchn[1] = dip; Ctrl->mchn[2] = rake; @@ -227,7 +245,7 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ // 单力源 case 'F': Ctrl->F.active = true; - Ctrl->computeType = GRT_SYN_COMPUTE_SF; + Ctrl->computeType = GRT_SYN_SF; { real_t fn, fe, fz; sprintf(Ctrl->s_computeType, "%s", "SF"); @@ -243,7 +261,7 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ // 张量震源 case 'T': Ctrl->T.active = true; - Ctrl->computeType = GRT_SYN_COMPUTE_MT; + Ctrl->computeType = GRT_SYN_MT; { real_t Mxx, Mxy, Mxz, Myy, Myz, Mzz; sprintf(Ctrl->s_computeType, "%s", "MT"); @@ -269,6 +287,11 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ Ctrl->N.active = true; break; + // 不打印在终端 + case 's': + Ctrl->s.active = true; + break; + GRT_Common_Options_in_Switch((char)(optopt)); } } @@ -500,7 +523,7 @@ int static_syn_main(int argc, char **argv){ memset(tmpsyn, 0, sizeof(real_t)*GRT_CHANNEL_NUM); // 计算震源辐射因子 - grt_set_source_radiation(Ctrl->srcRadi, Ctrl->computeType, (ityp > 0) && GRT_ZRT_CODES[ityp-1]=='T', Ctrl->S.M0, upar_scale, azrad, Ctrl->mchn); + grt_set_source_radiation(Ctrl->srcRadi, Ctrl->computeType, (ityp > 0) && GRT_ZRT_CODES[ityp-1]=='T', Ctrl->S.M0, upar_scale, src_va/src_vb, azrad, Ctrl->mchn); // 合成 GRT_LOOP_ChnlGrid(im, c){ @@ -557,6 +580,10 @@ int static_syn_main(int argc, char **argv){ NC_CHECK(nc_close(in_ncid)); NC_CHECK(nc_close(out_ncid)); + if(! Ctrl->s.active) { + GRTRaiseInfo("Synthetic static displacements of %s source saved in \"%s\".", srcTypeFullName[Ctrl->computeType], Ctrl->O.s_outgrid); + } + // 释放内存 GRT_LOOP_ChnlGrid(im, c){ GRT_SAFE_FREE_PTR(u[im][c]); diff --git a/pygrt/c_interfaces.py b/pygrt/c_interfaces.py index a546c71b..96ae8455 100755 --- a/pygrt/c_interfaces.py +++ b/pygrt/c_interfaces.py @@ -142,6 +142,15 @@ def set_num_threads(n): C_grt_py_attenuation_law.argtypes = [REAL, DPOINTER, DPOINTER] # double, double[2], double[2] +# ------------------------------------------------------------------- +# C 函数定义的方向因子 +# ------------------------------------------------------------------- +C_grt_set_source_radiation = libgrt.grt_set_source_radiation +C_grt_set_source_radiation.restype = None +C_grt_set_source_radiation.argtypes = [ + (REAL*CHANNEL_NUM)*SRC_M_NUM, c_int, c_bool, + REAL, REAL, REAL, REAL, REAL*MECHANISM_NUM +] diff --git a/pygrt/c_structures.py b/pygrt/c_structures.py index a94455ca..ce17a517 100755 --- a/pygrt/c_structures.py +++ b/pygrt/c_structures.py @@ -22,6 +22,7 @@ "ZRTchs", "ZNEchs", "qwvchs", + "MECHANISM_NUM", "NPCT_REAL_TYPE", "NPCT_CMPLX_TYPE", @@ -44,7 +45,7 @@ ZRTchs = ['Z', 'R', 'T'] ZNEchs = ['Z', 'N', 'E'] qwvchs = ['q', 'w', 'v'] - +MECHANISM_NUM = 6 NPCT_REAL_TYPE = 'f8' NPCT_CMPLX_TYPE = 'c16' diff --git a/pygrt/utils.py b/pygrt/utils.py index 79b7276b..648ef409 100755 --- a/pygrt/utils.py +++ b/pygrt/utils.py @@ -4,7 +4,7 @@ :date: 2024-07-24 该文件包含一些数据处理操作上的补充: - 1、剪切源、单力源、爆炸源、矩张量源 通过格林函数合成理论地震图的函数\n + 1、剪切源、张位错源、单力源、爆炸源、矩张量源 通过格林函数合成理论地震图的函数\n 2、Stream类型的时域卷积、微分、积分 (基于numpy和scipy) \n 3、Stream类型写到本地sac文件，自定义名称 \n 4、读取波数积分和峰谷平均法过程文件 \n @@ -35,6 +35,17 @@ from .c_interfaces import * +from enum import Enum, unique + +@unique +class GRT_SYN_TYPE(Enum): + GRT_SYN_EX = 0 + GRT_SYN_SF = 1 + GRT_SYN_DC = 2 + GRT_SYN_TS = 3 + GRT_SYN_MT = 4 + + __all__ = [ "gen_syn_from_gf_DC", @@ -67,15 +78,13 @@ # #================================================================================================================= -def _gen_syn_from_gf(st:Stream, calc_upar:bool, compute_type:str, M0:float, az:float, ZNE=False, **kwargs): +def _gen_syn_from_gf(st:Stream, calc_upar:bool, compute_type:GRT_SYN_TYPE, M0:float, az:float, ZNE=False, **kwargs): r""" 一个发起函数，根据不同震源参数，从格林函数中合成理论地震图 :param st: 计算好的时域格林函数, :class:`obspy.Stream` 类型 :param calc_upar: 是否计算位移u的空间导数 - :param compute_type: 计算类型，应为以下之一: - 'COMPUTE_EX'(爆炸源), 'COMPUTE_SF'(单力源), - 'COMPUTE_DC'(剪切源), 'COMPUTE_MT'(矩张量源) + :param compute_type: 计算震源类型 :param M0: 标量地震矩, 单位dyne*cm :param az: 方位角(度) :param ZNE: 是否以ZNE分量输出? @@ -87,6 +96,11 @@ def _gen_syn_from_gf(st:Stream, calc_upar:bool, compute_type:str, M0:float, az:f srcName = ["EX", "VF", "HF", "DD", "DS", "SS"] allchs = [tr.stats.channel for tr in st] + # 为张位错计算 Vp/Vs + src_va = st[0].stats.sac['user6'] + src_vb = st[0].stats.sac['user6'] + VpVs_ratio = src_va / src_vb + baz = 180 + az if baz > 360: baz -= 360 @@ -105,7 +119,7 @@ def _gen_syn_from_gf(st:Stream, calc_upar:bool, compute_type:str, M0:float, az:f if ityp == 3: upar_scale /= dist - srcRadi = _set_source_radi(ityp==3, upar_scale, compute_type, M0, azrad, **kwargs) + srcRadi = _set_source_radi(ityp==3, upar_scale, compute_type, M0, azrad, VpVs_ratio=VpVs_ratio, **kwargs) inpref = sacin_prefixes[ityp] outpref = sacout_prefixes[ityp] @@ -138,15 +152,13 @@ def _gen_syn_from_gf(st:Stream, calc_upar:bool, compute_type:str, M0:float, az:f return stall -def _gen_syn_from_static_gf(grn:dict, calc_upar:bool, compute_type:str, M0:float, ZNE=False, **kwargs): +def _gen_syn_from_static_gf(grn:dict, calc_upar:bool, compute_type:GRT_SYN_TYPE, M0:float, ZNE=False, **kwargs): r""" 一个发起函数，根据不同震源参数，从静态格林函数中合成理论静态场 :param grn: 计算好的静态格林函数, 字典类型 :param calc_upar: 是否计算位移u的空间导数 - :param compute_type: 计算类型，应为以下之一: - 'COMPUTE_EX'(爆炸源), 'COMPUTE_SF'(单力源), - 'COMPUTE_DC'(剪切源), 'COMPUTE_MT'(矩张量源) + :param compute_type: 计算震源类型 :param M0: 标量地震矩, 单位dyne*cm :param ZNE: 是否以ZNE分量输出? @@ -156,6 +168,11 @@ def _gen_syn_from_static_gf(grn:dict, calc_upar:bool, compute_type:str, M0:float srcName = ["EX", "VF", "HF", "DD", "DS", "SS"] allchs = list(grn.keys()) + # 为张位错计算 Vp/Vs + src_va = grn['_src_va'] + src_vb = grn['_src_vb'] + VpVs_ratio = src_va / src_vb + calcUTypes = 4 if calc_upar else 1 xarr:np.ndarray = grn['_xarr'] @@ -189,7 +206,7 @@ def _gen_syn_from_static_gf(grn:dict, calc_upar:bool, compute_type:str, M0:float if ityp == 3: upar_scale /= dist - srcRadi = _set_source_radi(ityp==3, upar_scale, compute_type, M0, azrad, **kwargs) + srcRadi = _set_source_radi(ityp==3, upar_scale, compute_type, M0, azrad, VpVs_ratio=VpVs_ratio, **kwargs) inpref = sacin_prefixes[ityp] @@ -319,18 +336,16 @@ def _data_zrt2zne(stall:Stream): def _set_source_radi( - par_theta:bool, coef:float, compute_type:str, M0:float, azrad:float, + par_theta:bool, coef:float, compute_type:GRT_SYN_TYPE, M0:float, azrad:float, fZ=None, fN=None, fE=None, strike=None, dip=None, rake=None, MT=None, **kwargs): r""" - 设置不同震源的方向因子矩阵 + 使用 C 函数计算不同震源的方向因子矩阵 :param par_theta: 是否求对theta的偏导 :param coef: 比例系数 - :param compute_type: 计算类型，应为以下之一: - 'COMPUTE_EX'(爆炸源), 'COMPUTE_SF'(单力源), - 'COMPUTE_DC'(剪切源), 'COMPUTE_MT'(矩张量源) + :param compute_type: 计算震源类型 :param M0: 地震矩 :param azrad: 方位角(弧度) @@ -339,100 +354,24 @@ def _set_source_radi( - 剪切源需要: strike, dip, rake - 矩张量源需要: MT=(Mxx, Mxy, Mxz, Myy, Myz, Mzz) """ - - caz = np.cos(azrad) - saz = np.sin(azrad) - src_coef = np.zeros((SRC_M_NUM, CHANNEL_NUM), dtype='f8') - - # 计算乘法因子 - if compute_type == 'COMPUTE_SF': - mult = 1e-15 * M0 * coef + src_coef = np.zeros((SRC_M_NUM, CHANNEL_NUM), dtype=NPCT_REAL_TYPE) + mchn = np.zeros((MECHANISM_NUM,), dtype=NPCT_REAL_TYPE) + if compute_type == GRT_SYN_TYPE.GRT_SYN_SF: + mchn[:3] = [fN, fE, fZ] + elif compute_type == GRT_SYN_TYPE.GRT_SYN_DC: + mchn[:3] = [strike, dip, rake] + elif compute_type == GRT_SYN_TYPE.GRT_SYN_TS: + mchn[:2] = [strike, dip] + elif compute_type == GRT_SYN_TYPE.GRT_SYN_MT: + mchn[:] = MT[:] else: - mult = 1e-20 * M0 * coef - - # 根据不同计算类型处理 - if compute_type == 'COMPUTE_EX': - # 爆炸源情况 - src_coef[0, 0] = src_coef[0, 1] = 0.0 if par_theta else mult # Z/R分量 - src_coef[0, 2] = 0.0 # T分量 - - elif compute_type == 'COMPUTE_SF': - # 单力源情况 - # 计算各向异性系数 - A0 = fZ * mult - A1 = (fN * caz + fE * saz) * mult - A4 = (-fN * saz + fE * caz) * mult - - # 设置震源系数矩阵 (公式4.6.20) - src_coef[1, 0] = src_coef[1, 1] = 0.0 if par_theta else A0 # VF, Z/R - src_coef[2, 0] = src_coef[2, 1] = A4 if par_theta else A1 # HF, Z/R - src_coef[1, 2] = 0.0 # VF, T - src_coef[2, 2] = -A1 if par_theta else A4 # HF, T + raise ValueError("Unsupported source type.") - elif compute_type == 'COMPUTE_DC': - # 剪切源情况 (公式4.8.35) - # 计算各种角度值(转为弧度) - stkrad = np.deg2rad(strike) # 走向角 - diprad = np.deg2rad(dip) # 倾角 - rakrad = np.deg2rad(rake) # 滑动角 - therad = azrad - stkrad # 方位角与走向角差 - - # 计算各种三角函数值 - srak = np.sin(rakrad); crak = np.cos(rakrad) - sdip = np.sin(diprad); cdip = np.cos(diprad) - sdip2 = 2.0 * sdip * cdip; cdip2 = 2.0 * cdip**2 - 1.0 - sthe = np.sin(therad); cthe = np.cos(therad) - sthe2 = 2.0 * sthe * cthe; cthe2 = 2.0 * cthe**2 - 1.0 - - # 计算各向异性系数 - A0 = mult * (0.5 * sdip2 * srak) - A1 = mult * (cdip * crak * cthe - cdip2 * srak * sthe) - A2 = mult * (0.5 * sdip2 * srak * cthe2 + sdip * crak * sthe2) - A4 = mult * (-cdip2 * srak * cthe - cdip * crak * sthe) - A5 = mult * (sdip * crak * cthe2 - 0.5 * sdip2 * srak * sthe2) - - # 设置震源系数矩阵 - src_coef[3, 0] = src_coef[3, 1] = 0.0 if par_theta else A0 # DD, Z/R - src_coef[4, 0] = src_coef[4, 1] = A4 if par_theta else A1 # DS, Z/R - src_coef[5, 0] = src_coef[5, 1] = 2.0 * A5 if par_theta else A2 # SS, Z/R - src_coef[3, 2] = 0.0 # DD, T - src_coef[4, 2] = -A1 if par_theta else A4 # DS, T - src_coef[5, 2] = -2.0 * A2 if par_theta else A5 # DS, T + C_grt_set_source_radiation( + src_coef, compute_type.value, par_theta, + M0, coef, kwargs['VpVs_ratio'], azrad, mchn) - elif compute_type == 'COMPUTE_MT': - # 矩张量源情况 (公式4.9.7，修改了各向同性项) - # 初始化矩张量分量 - M11, M12, M13, M22, M23, M33 = MT - - # 计算各向同性部分并减去 - Mexp = (M11 + M22 + M33) / 3.0 - M11 -= Mexp - M22 -= Mexp - M33 -= Mexp - - # 计算方位角的2倍角三角函数 - caz2 = np.cos(2 * azrad) - saz2 = np.sin(2 * azrad) - - # 计算各向异性系数 - A0 = mult * ((2.0 * M33 - M11 - M22) / 6.0) - A1 = mult * (- (M13 * caz + M23 * saz)) - A2 = mult * (0.5 * (M11 - M22) * caz2 + M12 * saz2) - A4 = mult * (M13 * saz - M23 * caz) - A5 = mult * (-0.5 * (M11 - M22) * saz2 + M12 * caz2) - - # 设置震源系数矩阵 - src_coef[0, 0] = src_coef[0, 1] = 0.0 if par_theta else mult * Mexp # EX, Z/R - src_coef[3, 0] = src_coef[3, 1] = 0.0 if par_theta else A0 # DD, Z/R - src_coef[4, 0] = src_coef[4, 1] = A4 if par_theta else A1 # DS, Z/R - src_coef[5, 0] = src_coef[5, 1] = 2.0 * A5 if par_theta else A2 # SS, Z/R - src_coef[0, 2] = 0.0 # EX, T - src_coef[3, 2] = 0.0 # DD, T - src_coef[4, 2] = -A1 if par_theta else A4 # DS, T - src_coef[5, 2] = -2.0 * A2 if par_theta else A5 # DS, T - - return src_coef @@ -455,9 +394,33 @@ def gen_syn_from_gf_DC(st:Union[Stream,dict], M0:float, strike:float, dip:float, if isinstance(st, Stream): if az > 360 or az < -360: raise ValueError(f"WRONG azimuth ({az})") - return _gen_syn_from_gf(st, calc_upar, "COMPUTE_DC", M0, az, ZNE, strike=strike, dip=dip, rake=rake) + return _gen_syn_from_gf(st, calc_upar, GRT_SYN_TYPE.GRT_SYN_DC, M0, az, ZNE, strike=strike, dip=dip, rake=rake) + elif isinstance(st, dict): + return _gen_syn_from_static_gf(st, calc_upar, GRT_SYN_TYPE.GRT_SYN_DC, M0, ZNE, strike=strike, dip=dip, rake=rake) + else: + raise NotImplementedError + +def gen_syn_from_gf_TS(st:Union[Stream,dict], M0:float, strike:float, dip:float, az:float=-999, ZNE=False, calc_upar:bool=False): + ''' + Tension source, the unit of angles is all degrees(°) + + :param st: Green's functions in a :class:`obspy.Stream` (dynamic-case) or a dict (static-case) + :param M0: scalar seismic moment (dyne*cm) + :param strike: 0 <= strike <= 360 (north=0, clockwise as positive) + :param dip: 0 <= dip <= 90 + :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. + + :return: + - **stream** - :class:`obspy.Stream` + ''' + if isinstance(st, Stream): + if az > 360 or az < -360: + raise ValueError(f"WRONG azimuth ({az})") + return _gen_syn_from_gf(st, calc_upar, GRT_SYN_TYPE.GRT_SYN_TS, M0, az, ZNE, strike=strike, dip=dip) elif isinstance(st, dict): - return _gen_syn_from_static_gf(st, calc_upar, "COMPUTE_DC", M0, ZNE, strike=strike, dip=dip, rake=rake) + return _gen_syn_from_static_gf(st, calc_upar, GRT_SYN_TYPE.GRT_SYN_TS, M0, ZNE, strike=strike, dip=dip) else: raise NotImplementedError @@ -481,9 +444,9 @@ def gen_syn_from_gf_SF(st:Union[Stream,dict], S:float, fN:float, fE:float, fZ:fl if isinstance(st, Stream): if az > 360 or az < -360: raise ValueError(f"WRONG azimuth ({az})") - return _gen_syn_from_gf(st, calc_upar, "COMPUTE_SF", S, az, ZNE, fN=fN, fE=fE, fZ=fZ) + return _gen_syn_from_gf(st, calc_upar, GRT_SYN_TYPE.GRT_SYN_SF, S, az, ZNE, fN=fN, fE=fE, fZ=fZ) elif isinstance(st, dict): - return _gen_syn_from_static_gf(st, calc_upar, "COMPUTE_SF", S, ZNE, fN=fN, fE=fE, fZ=fZ) + return _gen_syn_from_static_gf(st, calc_upar, GRT_SYN_TYPE.GRT_SYN_SF, S, ZNE, fN=fN, fE=fE, fZ=fZ) else: raise NotImplementedError @@ -504,9 +467,9 @@ def gen_syn_from_gf_EX(st:Union[Stream,dict], M0:float, az:float=-999, ZNE=False if isinstance(st, Stream): if az > 360 or az < -360: raise ValueError(f"WRONG azimuth ({az})") - return _gen_syn_from_gf(st, calc_upar, "COMPUTE_EX", M0, az, ZNE) + return _gen_syn_from_gf(st, calc_upar, GRT_SYN_TYPE.GRT_SYN_EX, M0, az, ZNE) elif isinstance(st, dict): - return _gen_syn_from_static_gf(st, calc_upar, "COMPUTE_EX", M0, ZNE) + return _gen_syn_from_static_gf(st, calc_upar, GRT_SYN_TYPE.GRT_SYN_EX, M0, ZNE) else: raise NotImplementedError @@ -528,9 +491,9 @@ def gen_syn_from_gf_MT(st:Union[Stream,dict], M0:float, MT:ArrayLike, az:float=- if isinstance(st, Stream): if az > 360 or az < -360: raise ValueError(f"WRONG azimuth ({az})") - return _gen_syn_from_gf(st, calc_upar, "COMPUTE_MT", M0, az, ZNE, MT=MT) + return _gen_syn_from_gf(st, calc_upar, GRT_SYN_TYPE.GRT_SYN_MT, M0, az, ZNE, MT=MT) elif isinstance(st, dict): - return _gen_syn_from_static_gf(st, calc_upar, "COMPUTE_MT", M0, ZNE, MT=MT) + return _gen_syn_from_static_gf(st, calc_upar, GRT_SYN_TYPE.GRT_SYN_MT, M0, ZNE, MT=MT) else: raise NotImplementedError From 2a50aded208374e01f9716de89991ee5b345cc84 Mon Sep 17 00:00:00 2001 From: Dengda98 Date: Wed, 31 Dec 2025 15:04:17 +0800 Subject: [PATCH 02/12] fix: update _set_source_radi --- pygrt/utils.py | 10 +++++++--- 1 file changed, 7 insertions(+), 3 deletions(-) diff --git a/pygrt/utils.py b/pygrt/utils.py index 648ef409..b81fe9ec 100755 --- a/pygrt/utils.py +++ b/pygrt/utils.py @@ -355,9 +355,13 @@ def _set_source_radi( - 矩张量源需要: MT=(Mxx, Mxy, Mxz, Myy, Myz, Mzz) """ + VpVs_ratio = kwargs['VpVs_ratio'] if 'VpVs_ratio' in kwargs else 0.0 + src_coef = np.zeros((SRC_M_NUM, CHANNEL_NUM), dtype=NPCT_REAL_TYPE) mchn = np.zeros((MECHANISM_NUM,), dtype=NPCT_REAL_TYPE) - if compute_type == GRT_SYN_TYPE.GRT_SYN_SF: + if compute_type == GRT_SYN_TYPE.GRT_SYN_EX: + pass + elif compute_type == GRT_SYN_TYPE.GRT_SYN_SF: mchn[:3] = [fN, fE, fZ] elif compute_type == GRT_SYN_TYPE.GRT_SYN_DC: mchn[:3] = [strike, dip, rake] @@ -369,8 +373,8 @@ def _set_source_radi( raise ValueError("Unsupported source type.") C_grt_set_source_radiation( - src_coef, compute_type.value, par_theta, - M0, coef, kwargs['VpVs_ratio'], azrad, mchn) + npct.as_ctypes(src_coef), compute_type.value, par_theta, + M0, coef, VpVs_ratio, azrad, npct.as_ctypes(mchn)) return src_coef From 28ae7acd984d4fb9177bbcb01e1b7f486c5985f4 Mon Sep 17 00:00:00 2001 From: Dengda98 Date: Wed, 31 Dec 2025 19:38:22 +0800 Subject: [PATCH 03/12] fix --- pygrt/utils.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pygrt/utils.py b/pygrt/utils.py index b81fe9ec..ad2afbe2 100755 --- a/pygrt/utils.py +++ b/pygrt/utils.py @@ -98,7 +98,7 @@ def _gen_syn_from_gf(st:Stream, calc_upar:bool, compute_type:GRT_SYN_TYPE, M0:fl # 为张位错计算 Vp/Vs src_va = st[0].stats.sac['user6'] - src_vb = st[0].stats.sac['user6'] + src_vb = st[0].stats.sac['user7'] VpVs_ratio = src_va / src_vb baz = 180 + az From 54a16200db37043ac36329e5baff611bde42456a Mon Sep 17 00:00:00 2001 From: Dengda98 Date: Wed, 31 Dec 2025 20:12:18 +0800 Subject: [PATCH 04/12] fix --- pygrt/C_extension/src/dynamic/grt_syn.c | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/pygrt/C_extension/src/dynamic/grt_syn.c b/pygrt/C_extension/src/dynamic/grt_syn.c index 9523f8f2..02ea4894 100644 --- a/pygrt/C_extension/src/dynamic/grt_syn.c +++ b/pygrt/C_extension/src/dynamic/grt_syn.c @@ -461,15 +461,15 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ if(va <= 0.0 || vb < 0.0 || rho <= 0.0){ GRTRaiseError("Bad src_va, src_vb or src_rho in \"%s\" header.\n", entry->d_name); } - if(vb == 0.0){ - GRTRaiseError("Zero src_vb in \"%s\" header. " - "Maybe you try to use -Su but the source is in the liquid. " - "Use -S instead.\n" , entry->d_name); - } - Ctrl->VpVs_ratio = ((real_t)va)/vb; + Ctrl->VpVs_ratio = (vb == 0.0)? 0.0 : ((real_t)va)/vb; if (Ctrl->S.mult_src_mu) { + if(vb == 0.0){ + GRTRaiseError("Zero src_vb in \"%s\" header. " + "Maybe you try to use -Su but the source is in the liquid. " + "Use -S instead.\n" , entry->d_name); + } Ctrl->S.src_mu = vb*vb*rho*1e10; Ctrl->S.M0 *= Ctrl->S.src_mu; } From bfe7b83daa78d01ae754607cfade7f44d1575a3c Mon Sep 17 00:00:00 2001 From: Dengda98 Date: Wed, 31 Dec 2025 20:20:59 +0800 Subject: [PATCH 05/12] DOC: update Gallery ex13 --- docs/source/Gallery/ex13/plot.py | 33 ++++++++++++++++++-------------- docs/source/Gallery/ex13/run.sh | 7 +++++++ docs/source/Gallery/ex14/run.sh | 0 3 files changed, 26 insertions(+), 14 deletions(-) mode change 100644 => 100755 docs/source/Gallery/ex14/run.sh diff --git a/docs/source/Gallery/ex13/plot.py b/docs/source/Gallery/ex13/plot.py index d90b5a77..66b37d66 100644 --- a/docs/source/Gallery/ex13/plot.py +++ b/docs/source/Gallery/ex13/plot.py @@ -3,7 +3,7 @@ from scipy.io import netcdf_file def plot_disp(fLst): - fig, axs = plt.subplots(2, 3, figsize=(12, 6)) + fig, axs = plt.subplots(3, 3, figsize=(12, 8)) for i, f in enumerate(fLst): ax3 = axs[i] chs = ['E', 'N', 'Z'] @@ -18,11 +18,12 @@ def plot_disp(fLst): axs[0, 0].set_ylabel("Strike-slip", fontsize=12) axs[1, 0].set_ylabel("Dip-slip", fontsize=12) + axs[2, 0].set_ylabel("Tensile", fontsize=12) return fig, axs def plot_stress(fLst): - fig, axs = plt.subplots(2, 2, figsize=(9,7), gridspec_kw=dict(hspace=0.2)) - for i, f in enumerate([f1, f2]): + fig, axs = plt.subplots(2, 3, figsize=(12, 7), gridspec_kw=dict(hspace=0.2)) + for i, f in enumerate(fLst): ax = axs[0, i] Y = f.variables['north'][:] ax.plot(Y, f.variables['stress_ZZ'][:,0], 'bo-', ms=3, label='ZZ') @@ -38,33 +39,37 @@ def plot_stress(fLst): axs[0, 0].set_title("Strike-slip", fontsize=12) axs[0, 1].set_title("Dip-slip", fontsize=12) + axs[0, 2].set_title("Tensile", fontsize=12) axs[0, 0].set_ylabel("Normal", fontsize=12) axs[1, 0].set_ylabel("Shear", fontsize=12) - axs[0, 1].legend(bbox_to_anchor=(1.3, 1)) - axs[1, 1].legend(bbox_to_anchor=(1.3, 1)) + axs[0, 2].legend(bbox_to_anchor=(1.3, 1)) + axs[1, 2].legend(bbox_to_anchor=(1.3, 1)) return fig, axs # ====================================================================== with netcdf_file("stsyn_1_stkslip.nc", mmap=False) as f1, \ - netcdf_file("stsyn_1_dipslip.nc", mmap=False) as f2: - fig, axs = plot_disp([f1, f2]) + netcdf_file("stsyn_1_dipslip.nc", mmap=False) as f2, \ + netcdf_file("stsyn_1_tensile.nc", mmap=False) as f3: + fig, axs = plot_disp([f1, f2, f3]) fig.savefig("disp1.svg", bbox_inches='tight') - fig, axs = plot_stress([f1, f2]) + fig, axs = plot_stress([f1, f2, f3]) fig.savefig("stress1.svg", bbox_inches='tight') with netcdf_file("stsyn_2_stkslip.nc", mmap=False) as f1, \ - netcdf_file("stsyn_2_dipslip.nc", mmap=False) as f2: - fig, axs = plot_disp([f1, f2]) + netcdf_file("stsyn_2_dipslip.nc", mmap=False) as f2, \ + netcdf_file("stsyn_2_tensile.nc", mmap=False) as f3: + fig, axs = plot_disp([f1, f2, f3]) fig.savefig("disp2.svg", bbox_inches='tight') - fig, axs = plot_stress([f1, f2]) + fig, axs = plot_stress([f1, f2, f3]) fig.savefig("stress2.svg", bbox_inches='tight') with netcdf_file("stsyn_3_stkslip.nc", mmap=False) as f1, \ - netcdf_file("stsyn_3_dipslip.nc", mmap=False) as f2: - fig, axs = plot_disp([f1, f2]) + netcdf_file("stsyn_3_dipslip.nc", mmap=False) as f2, \ + netcdf_file("stsyn_3_tensile.nc", mmap=False) as f3: + fig, axs = plot_disp([f1, f2, f3]) fig.savefig("disp3.svg", bbox_inches='tight') - fig, axs = plot_stress([f1, f2]) + fig, axs = plot_stress([f1, f2, f3]) fig.savefig("stress3.svg", bbox_inches='tight') \ No newline at end of file diff --git a/docs/source/Gallery/ex13/run.sh b/docs/source/Gallery/ex13/run.sh index 2e5c3d2c..b6624d18 100755 --- a/docs/source/Gallery/ex13/run.sh +++ b/docs/source/Gallery/ex13/run.sh @@ -13,6 +13,9 @@ grt static stress stsyn_1_stkslip.nc grt static syn -M90/70/90 -Su1e6 -e -N -Gstgrn1.nc -Ostsyn_1_dipslip.nc grt static stress stsyn_1_dipslip.nc +grt static syn -M90/70 -Su1e6 -e -N -Gstgrn1.nc -Ostsyn_1_tensile.nc +grt static stress stsyn_1_tensile.nc + grt static greenfn -M$mod -D10/14 -Y3/3/1 -X-10/10/0.5 -e -Ostgrn2.nc @@ -22,6 +25,8 @@ grt static stress stsyn_2_stkslip.nc grt static syn -M90/70/90 -Su1e6 -e -N -Gstgrn2.nc -Ostsyn_2_dipslip.nc grt static stress stsyn_2_dipslip.nc +grt static syn -M90/70 -Su1e6 -e -N -Gstgrn2.nc -Ostsyn_2_tensile.nc +grt static stress stsyn_2_tensile.nc mod="mod" grt static greenfn -M$mod -D5/0 -Y3/3/1 -X-10/10/0.5 -e -Ostgrn3.nc @@ -32,6 +37,8 @@ grt static stress stsyn_3_stkslip.nc grt static syn -M90/70/90 -Su1e6 -e -N -Gstgrn3.nc -Ostsyn_3_dipslip.nc grt static stress stsyn_3_dipslip.nc +grt static syn -M90/70 -Su1e6 -e -N -Gstgrn3.nc -Ostsyn_3_tensile.nc +grt static stress stsyn_3_tensile.nc python plot.py diff --git a/docs/source/Gallery/ex14/run.sh b/docs/source/Gallery/ex14/run.sh old mode 100644 new mode 100755 From 1681bfa385a9019098353c55594de5e117675178 Mon Sep 17 00:00:00 2001 From: Dengda98 Date: Wed, 31 Dec 2025 20:21:51 +0800 Subject: [PATCH 06/12] DOC: update tutorial (dynamic) --- docs/source/Tutorial/dynamic/run/run.py | 17 +++++++++++++++++ docs/source/Tutorial/dynamic/run/run.sh | 7 +++++++ docs/source/Tutorial/dynamic/syn.rst | 25 +++++++++++++++++++++++++ 3 files changed, 49 insertions(+) diff --git a/docs/source/Tutorial/dynamic/run/run.py b/docs/source/Tutorial/dynamic/run/run.py index dd880c9f..6ad3d1ea 100755 --- a/docs/source/Tutorial/dynamic/run/run.py +++ b/docs/source/Tutorial/dynamic/run/run.py @@ -171,6 +171,23 @@ def plot_int_dif(stsyn:Stream, stsyn_int:Stream, stsyn_dif:Stream, comp:str, out # ----------------------------------------------------------------------------------- +# ----------------------------------------------------------------------------------- +# BEGIN SYN TS +# 接之前的代码 +idx = 2 +stgrn = stgrnLst[idx] # 选择格林函数 + +stsyn = pygrt.utils.gen_syn_from_gf_TS(stgrn, M0=1e24, strike=33, dip=50, az=30) +print(stsyn) +# 3 Trace(s) in Stream: +# .SYN..Z | 1970-01-01T00:00:00.000000Z - 1970-01-01T00:00:09.980000Z | 50.0 Hz, 500 samples +# .SYN..R | 1970-01-01T00:00:00.000000Z - 1970-01-01T00:00:09.980000Z | 50.0 Hz, 500 samples +# .SYN..T | 1970-01-01T00:00:00.000000Z - 1970-01-01T00:00:09.980000Z | 50.0 Hz, 500 samples +plot_syn(stsyn, "syn_ts.svg") +# END SYN TS +# ----------------------------------------------------------------------------------- + + # ----------------------------------------------------------------------------------- # BEGIN SYN MT idx = 2 diff --git a/docs/source/Tutorial/dynamic/run/run.sh b/docs/source/Tutorial/dynamic/run/run.sh index e4bffc46..d9418148 100755 --- a/docs/source/Tutorial/dynamic/run/run.sh +++ b/docs/source/Tutorial/dynamic/run/run.sh @@ -48,6 +48,13 @@ grt syn -GGRN/milrow_2_0_10 -S1e24 -A30 -M33/50/120 -Osyn_dc # END SYN DC # ----------------------------------------------------------------------------------- +# ----------------------------------------------------------------------------------- +# BEGIN SYN TS +# 合成结果在 syn_ts/ 目录下，以SAC格式保存。 +grt syn -GGRN/milrow_2_0_10 -S1e24 -A30 -M33/50 -Osyn_ts +# END SYN TS +# ----------------------------------------------------------------------------------- + # ----------------------------------------------------------------------------------- # BEGIN SYN MT # 合成结果在 syn_mt/ 目录下，以SAC格式保存。 diff --git a/docs/source/Tutorial/dynamic/syn.rst b/docs/source/Tutorial/dynamic/syn.rst index a8ad7fbe..0b261c40 100644 --- a/docs/source/Tutorial/dynamic/syn.rst +++ b/docs/source/Tutorial/dynamic/syn.rst @@ -126,6 +126,31 @@ Python中合成动态位移的主函数为 :func:`gen_syn_from_gf_*() Date: Wed, 31 Dec 2025 20:22:29 +0800 Subject: [PATCH 07/12] DOC: update tutorial (static) --- docs/source/Tutorial/static/run/run.py | 18 +++++++ docs/source/Tutorial/static/run/run.sh | 35 ++++++++++++- docs/source/Tutorial/static/run/tension2mt.py | 34 ++++++++++++ docs/source/Tutorial/static/static_syn.rst | 52 +++++++++++++++++++ 4 files changed, 138 insertions(+), 1 deletion(-) create mode 100644 docs/source/Tutorial/static/run/tension2mt.py diff --git a/docs/source/Tutorial/static/run/run.py b/docs/source/Tutorial/static/run/run.py index 215c6b28..09558d8b 100755 --- a/docs/source/Tutorial/static/run/run.py +++ b/docs/source/Tutorial/static/run/run.py @@ -77,6 +77,24 @@ def plot_static(static_syn:dict, out:Union[str,None]=None): # END SYN DC2 # --------------------------------------------------------------------------------- +# --------------------------------------------------------------------------------- +# 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']) +plot_static(static_syn, "syn_ts.svg") +# END SYN TS +# --------------------------------------------------------------------------------- + + +# --------------------------------------------------------------------------------- +# 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']) +plot_static(static_syn, "syn_ts2.svg") +# END SYN TS2 +# --------------------------------------------------------------------------------- # --------------------------------------------------------------------------------- diff --git a/docs/source/Tutorial/static/run/run.sh b/docs/source/Tutorial/static/run/run.sh index 62cb679f..431bd839 100755 --- a/docs/source/Tutorial/static/run/run.sh +++ b/docs/source/Tutorial/static/run/run.sh @@ -104,6 +104,39 @@ EOF gmt colorbar -Bx+l"Z (cm)" gmt end +# --------------------------------------------------------------------------------- +# BEGIN SYN TS +# 从网格文件中读取格林函数，再将合成结果写入新网格 +grt static syn -S1e24 -M33/50 -N -Gstgrn.nc -Ostsyn_ts.nc +# END SYN TS +# --------------------------------------------------------------------------------- + +# 计算张位错的矩张量表示，仅绘制其中的 DC+CLVD 分量 +gmt begin syn_ts pdf + gmtplot_static stsyn_ts.nc -Si0.03c + gmt meca -Sz0.5c < Date: Wed, 31 Dec 2025 20:22:57 +0800 Subject: [PATCH 08/12] TEST: update syn --- test/static_syn/test_static_syn.py | 1 + test/static_syn/test_static_syn.sh | 1 + test/syn/test_syn.py | 1 + test/syn/test_syn.sh | 1 + 4 files changed, 4 insertions(+) diff --git a/test/static_syn/test_static_syn.py b/test/static_syn/test_static_syn.py index d15af83e..e5d6aef6 100644 --- a/test/static_syn/test_static_syn.py +++ b/test/static_syn/test_static_syn.py @@ -16,6 +16,7 @@ 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) stsyn = pygrt.utils.gen_syn_from_gf_DC(stgrn, 1e20, strike=33, dip=44, rake=55, az=22) +stsyn = pygrt.utils.gen_syn_from_gf_TS(stgrn, 1e20, strike=33, dip=44, az=22) stsyn = pygrt.utils.gen_syn_from_gf_MT(stgrn, 1e20, [1, -2, -5, 0.5, 3, 1.2], az=22) diff --git a/test/static_syn/test_static_syn.sh b/test/static_syn/test_static_syn.sh index b5c9f438..5afed0ca 100755 --- a/test/static_syn/test_static_syn.sh +++ b/test/static_syn/test_static_syn.sh @@ -11,6 +11,7 @@ grt static syn -S1e20 -Gstgrn.nc -Ostsyn.nc grt static syn -S1e20 -F2/-1/4 -Gstgrn.nc -Ostsyn.nc grt static syn -S1e20 -M77/88/111 -Gstgrn.nc -Ostsyn.nc grt static syn -Su1e6 -M77/88/111 -Gstgrn.nc -Ostsyn.nc +grt static syn -Su1e6 -M77/88 -Gstgrn.nc -Ostsyn.nc grt static syn -S1e20 -T1/-2/-5/0.5/3/1.2 -Gstgrn.nc -Ostsyn.nc grt static syn -S1e20 -F2/-1/4 -e -Gstgrn.nc -Ostsyn.nc diff --git a/test/syn/test_syn.py b/test/syn/test_syn.py index 9bf9ea27..880bef45 100644 --- a/test/syn/test_syn.py +++ b/test/syn/test_syn.py @@ -18,6 +18,7 @@ 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) stsyn = pygrt.utils.gen_syn_from_gf_DC(stgrn, 1e20, strike=33, dip=44, rake=55, az=22) +stsyn = pygrt.utils.gen_syn_from_gf_TS(stgrn, 1e20, strike=33, dip=44, az=22) stsyn = pygrt.utils.gen_syn_from_gf_MT(stgrn, 1e20, [1, -2, -5, 0.5, 3, 1.2], az=22) diff --git a/test/syn/test_syn.sh b/test/syn/test_syn.sh index 60960841..1e8b7904 100755 --- a/test/syn/test_syn.sh +++ b/test/syn/test_syn.sh @@ -10,6 +10,7 @@ grt syn -GGRN/milrow_2_0_10 -A22 -S1e20 -Osyn grt syn -GGRN/milrow_2_0_10 -A22 -S1e16 -F-1/2/-4 -Osyn grt syn -GGRN/milrow_2_0_10 -A22 -S1e20 -M33/44/55 -Osyn grt syn -GGRN/milrow_2_0_10 -A22 -Su1e10 -M33/44/55 -Osyn +grt syn -GGRN/milrow_2_0_10 -A22 -Su1e10 -M33/44 -Osyn grt syn -GGRN/milrow_2_0_10 -A22 -S1e20 -T1/-2/-5/0.5/3/1.2 -Osyn grt syn -GGRN/milrow_2_0_10 -A22 -S1e20 -Dp/0.6 -Osyn From 9b80f25dc24fb250c743671ba4c344167af468ee Mon Sep 17 00:00:00 2001 From: Dengda98 Date: Wed, 31 Dec 2025 20:24:04 +0800 Subject: [PATCH 09/12] TEST: update _correct* --- test/_correct_full_wave/_Ref/syn_ts/R.sac | Bin 0 -> 3032 bytes test/_correct_full_wave/_Ref/syn_ts/T.sac | Bin 0 -> 3032 bytes test/_correct_full_wave/_Ref/syn_ts/Z.sac | Bin 0 -> 3032 bytes test/_correct_full_wave/_Ref/syn_ts/rR.sac | Bin 0 -> 3032 bytes test/_correct_full_wave/_Ref/syn_ts/rT.sac | Bin 0 -> 3032 bytes test/_correct_full_wave/_Ref/syn_ts/rZ.sac | Bin 0 -> 3032 bytes test/_correct_full_wave/_Ref/syn_ts/sig.sac | Bin 0 -> 756 bytes test/_correct_full_wave/_Ref/syn_ts/tR.sac | Bin 0 -> 3032 bytes test/_correct_full_wave/_Ref/syn_ts/tT.sac | Bin 0 -> 3032 bytes test/_correct_full_wave/_Ref/syn_ts/tZ.sac | Bin 0 -> 3032 bytes test/_correct_full_wave/_Ref/syn_ts/zR.sac | Bin 0 -> 3032 bytes test/_correct_full_wave/_Ref/syn_ts/zT.sac | Bin 0 -> 3032 bytes test/_correct_full_wave/_Ref/syn_ts/zZ.sac | Bin 0 -> 3032 bytes test/_correct_full_wave/correctness.sh | 2 ++ test/_correct_static_disp/_Ref/stsyn_ts.nc | Bin 0 -> 68880 bytes test/_correct_static_disp/correctness.sh | 2 ++ 16 files changed, 4 insertions(+) create mode 100644 test/_correct_full_wave/_Ref/syn_ts/R.sac create mode 100644 test/_correct_full_wave/_Ref/syn_ts/T.sac create mode 100644 test/_correct_full_wave/_Ref/syn_ts/Z.sac create mode 100644 test/_correct_full_wave/_Ref/syn_ts/rR.sac create mode 100644 test/_correct_full_wave/_Ref/syn_ts/rT.sac create mode 100644 test/_correct_full_wave/_Ref/syn_ts/rZ.sac create mode 100644 test/_correct_full_wave/_Ref/syn_ts/sig.sac create mode 100644 test/_correct_full_wave/_Ref/syn_ts/tR.sac create mode 100644 test/_correct_full_wave/_Ref/syn_ts/tT.sac create mode 100644 test/_correct_full_wave/_Ref/syn_ts/tZ.sac create mode 100644 test/_correct_full_wave/_Ref/syn_ts/zR.sac create mode 100644 test/_correct_full_wave/_Ref/syn_ts/zT.sac create mode 100644 test/_correct_full_wave/_Ref/syn_ts/zZ.sac create mode 100644 test/_correct_static_disp/_Ref/stsyn_ts.nc diff --git a/test/_correct_full_wave/_Ref/syn_ts/R.sac b/test/_correct_full_wave/_Ref/syn_ts/R.sac new file mode 100644 index 0000000000000000000000000000000000000000..023d75d0b4d3798671ed81c4f5d2db9f2b86297f GIT binary patch literal 3032 zcmb_YdpK2D8$X#`hNC8z5l1=YQ0|7k*511#xurp-LN0}r5M@e6GEI$Z8YU$YIYMHF zBGuVz?Hy*25fP`-Al*k)nyIM_6`H=8@qB&F_s91=&#dQp*YEwk@4J5Mx0Z5ui4Cvc zuK7QL$Kz>KEwXLJ7 z3yHAQhf_Q7UIQO7y<> zN}BiN8sRgf%p@vw5mhR4n6B0B#3K_SgorN3C=Mql{Gkv4f@O<(I!rJL^a=t_fns_f2D z>Sup9sx2&^oTlSJ%C2t5dD(TC*s_4#xN#5rpfZm&9lpkHyRXURJP>gTYC_JfU6ZR# zmgkv&N>DqZA^m_kfS|cNZ&c5YE zuR9$=pZqa~9{M(eo+&D$uh#d{o?|667q*vv>T#MXYt5$0eJ7}(?;%~mGojaRvZO+! zBqiCADh;o1WCN0(VZlrZYmynq)|}kOibm7e+aG(gvjb`BQdS<#x7p1g#6K=O)H|v<4z^b(evbFs^Y*&jW_u1zH?$x0w zoW-=H7$2Vll~HQ^>`2J->i;^$;sgq>$C-cC77=7)nz|iLRIN=~UAs>Ps$xQn;8% z?dCa9)7Kc2BjeBUqM}r`eqfY+^|>N9@0K2S!@r9SJ*dR>KC|Tfh4VT2JF_`Ai$z=x zt-~eu8gto6I$ZMmx?Hl7j4gUp%N94Jv*cv~Ik@U;GOunvr9EpiT-hT>hor5frK1XrPnI5dVF9V0Sc3}N&PjLU-8j>&0B!_egNt+>-tTnQv z>d&iD_r~v$9B6RX8?g4tDBWh3z+KV?X0+oUw2}x*{HB=A7}teTnn& zr^Q#8qPd%~%6=nK<0X&!Wp5FANXY?5oHoFMUyGPR+hkM})QtLL7U3^G$$}@Q?m$Kt zlTeMRiB#q5L?jN<#oHvgXzsNt}EmyL`&$b5P*zrhgbvYV8M{lhEy*3#OUURo)?^&EA04ge*Y4*-ChV=y%?xnErG*&-Eirt0?5hUvN@u-ZtA zSnJI&x`p!S^1uR=(&~z$r?Ti|S~?=z>yg>QEYv_ZBPh#5uIn1ntNv!>+;Sax)_#xT z4z!{r{fE)+du=GfJO-7Wtw)=G%RnaOrHHTb0UDg%16J>8I>V{*r8mP2n7!kvuqjOe zjx6~K+**_cmM%lWF_Bt*_`GNAN>04(s)1;Z<-CBWU{xgu%jMeDK zY(9F{7mrdWXq0kw3%VaN8@bFVM1v(ZD2q6U{F38Pujn{ZS)0ym>{NyWE=Az%{X(GC za0-;B?*kh+btrIu1Y{Bi_=QV9C>FZIKHc}=W=$tJ8dL?eES=#Jxh_D2`N2mwW#C|+ zH*~Y=01;_Eu(;p?cs|1sW@waxwTUipG4KLCa50Sk(H!gwhVZWNHK8$X3r$}|i#E+T5>vC`do{GLP>;O+U@X^?;8(?C(GID76jyiS~O3K)V(@)R?se`my#XvU?UR+mnZiFB?O< SC8a2Bl`fnn%S2+}3V#RSSc)|O literal 0 HcmV?d00001 diff --git a/test/_correct_full_wave/_Ref/syn_ts/T.sac b/test/_correct_full_wave/_Ref/syn_ts/T.sac new file mode 100644 index 0000000000000000000000000000000000000000..b342e695c1b2ce10ef9ddfbb42d560a19950fdd8 GIT binary patch literal 3032 zcmb_XXH-?B(O~`Hx;xHoE^}Q5`QM82si*#;CldX>YyS5Oob2uMTy|zlUCy4F z=<>U#ch}i73Kx5Odza8q%RfAc#=P@;_+5PF_7|-W|Mx%enudXa{}z&Wy9euj;(zCZ zwSS??pEFwT&A#9NfAc}T_uoyjwXKjyvv#1H#_H}?^JDSF@E zhx)taBDPP}mD%d9!QO~AW3sXh*-=;Q8JhqTcF8o%Os{ulk0@=J#yc}vpW^9EWS%2? zJ1&qh+Hb=qSQ(SHGnweG3y(n3#K+aW4|yq1PcXqaxz;FM44(4tZEFdh(9LJ-s*yY z%#kE*o&l)bH-`Lp+z`0>S`w>L6EJtiB$8!e33!J|q$+I!*fZIV#7(vVKOMwGJqrNi zbxy=Qbt)*j?@Ej{+`uVqcT)f5ETC|jOX^L%!4VrTQk@qCb_eFO^syK|TecfC7$(6u zzc84ne;sz%UV^USmZ)#7B^uYT6j`qcK^rD)L92#tM#~*akhr)ItrS7xAzG*?7PrA|FEmzUe)x9Y8_7zmWtQ#$Azkn=iZlbcMr;vVCBl5m}6y5Bv zKwIka(A$wmk+fzTDyzvsX~C<|^$oEIdm^N$_eOUcMkBVt4V~`NLGn-sq}?q+AsJJU z{1etbMik29-^`fnN{5;NEcI zt|k^-Ej@+n^~;3|s%r$zu}MN!Y^rcamM#QEWe7cI)(P$3WeHa6ZG|$!`NHFI!Gimj z+JgSXUwDh*cD$YXT~=-2Q(Vxp23L4#NuS0~#4~r);^E*GRv+gmZI50gZ5Sm+T3J<&uS!yJ@f2@tS~&~1=$4`}nQmxge+K7QxR_6> zZbGjgF5%p{ZTyf2O1gfv0iAq#G{)x;J`|pdL#N4cPO_#HQftm`0YXj|(`jLoKJIG2 z4HdtP#dSRxCH-K5O*BWMl5+;Qp>G&%s#L)R8`pD%;T7!AVYbzD>`QzBm*0fBgI{ZdUj?bMIbgC$#~-gOtH2DV*No5^}I z=XJ)i;d32j8$X#NGy2|NW_Hk5Hv6!%tX>Pr7QNDyP2DDujT(IpYwXp;%O~H!wf$Y# zG;XqVc>k~XQRhSS=%NvxaAO5NCGW$5#lvt&ainyrsj+f$V! zX*|FgM9<`w_8p@pM&>js_CCkd#&CB!in*3vGw#$`1sCeA&b>)WxCa|F zY3-XPKGJM8|7lbrKb_k0m78brdIj2iS|7uMy60T(C==eK?Gbm-V+e1f%;k2O>+pV; zOF1v=My{*dlm5I-#_7KTL~uTb zTj{>NX&lLWN!1S&agp3{;&R~{HTCbYcWmpXs**4exTJ+z|I`Ll7RJ!%n6II*FN}^* z+YDt9sWeIT2l)Do0fqW0aI1IW=K^qExHzF$J`&n$y0whw2+ z4Dj216Mts4UAULM_x@iFfh>xot?A? zT4^H`dMFR>zhi?M?Dg=JK;sB+s-pNDM1(Wc zZBsXMRP0Q1x?WzMv zf6i<0tWJ$={BvaWssj4(Q$Q4XJR8I`DlfYu@DiTz(2dvU4C3=hqJb84E zMAb|r*;~zMY)cI50adAadJvm$)lQ}jt2~-GrGg}&W;^jqN){jA$!M1kqp>41$jvj4 z$X_74FQnprR}WgIZR zUC8JeC6bK}8lbg5oFo^yfn_*_*pzsI_@M_#kIHWRH4$QL;9#E&b;S0g);C8Y! eJstd%FcES|B*v|wWbeJ!-V}1Fq~y4iOL{VcshLS8D#Q#5@e56Ek#&+2h2&U`h`ySbL+pV#i=6@^3T3`8#o zMwdxbpLn7cuOx@ceMMi^-6Zag@gmcvdO{F-UgRKCNVFY4BhrxQAdDZDiIOCr6Q+u% zMFSm=i7hKnid;=}$&dH%7iCsDkjtNJ6$KX`CnFq!MR!iiQQj+ZEkActr+Vzl2*(jE z%DF0z5EyDvU*2^it~VP~C9WBSmaQ&TCTm0Zgv(L8`&C7r%CE@vs$QZURX?|yfQl>8=&ZM$=UyeS>Wrd8i2N6KQ^1d|!mA>SxAGe?hl zvv&_GPeDrgavVG0{54f{B9blX6H=-YA#7n(KGiV(IqOHHQsaJZaKg3+mbWO-*4rTM z{m_aYa`&J`9k%o)@746`@_F>I-6&KUZ-R#MfvD>nU9QnBl?%%2T_-mh@Aqc@Ie4aX1Z@385zI4sD?!fMVb_{B31 zoE7YYUndlEom-o@rXdMzBi@fI4UMtTeB!$HNF zw6&EbT<%+itZhO$cXv5FvE>u|S$_!Deq9MhZ8MrDvW@vQE!MoCFow4ZDd3gPR`GC$kdGS3;rkj>d16H@ze?c1e_^h~ zt6#0adT|BZJ>^VvMdSzS?Yp4jSPboQqM07mc}(N_C}yG4&y4Q!R%ZL!T&Az)4D-Y< zlc^0WUX_Z;qD;%Y+erZtA}`Cbs2`O zwfLAt7?vr@#A&=KHq!LMRWW9GY<&hUJ|`wUjSgXbqXS%Z%tM$2x%Q2hc z&cho(kuJ+v!)-nmnL%yOV#?xqL;ElZ$(mAf{ z9yk2FfV&0f4{hjSLbv1 zA2c)gbl-h^LYfV~E=8X2iudP7OLFLDgHG7s(@ys~`_Rs_<}*i95}AC(%Zz^gNoKu; zkeP{3FiveeQ*{k76Fr%%*e~tOzXPYbWqGC;PyNKkGZeH z2}5r59tC0M9EYXexM11Suee*w7V~yp$7r)l>$$MHRZzQT7k5BvCNvUD zz$d+RFg76{Kh80eetspye8W1N{IZ9u)3U@K`kQe>ku^?o+=IPpm9W#=5`2DVD%NQH z3=i+Q#3!;Gl!Avz44Y&gA?) znZcQ9g`uU9vyjJRC=1n=0v(NJz-x4Zo^uow4kbgrB@7ON1MskZ8GQCm1Q+eg0=miw zT5Ju6JB1~XRtkVyOfq5q^4)NG?>?Bc>maOwxlnr3b|~~e1zlPrpm@s>_-fuBs8VnU z`dy2HJ`JVt_P~BH#jc41up3TaDIt9 ztP3)L>V8E~YdjSqjXc=7@>`g%oDOds2!=BL(NH_W5%!#Op89iIFz~q*e0qBr*k9Fv zlQ|tA?bAW9yXgSX73Klye>PIf+nvB&b`e=2G^N_Sb&;}6GFyK&6Ww27iu_JoMBU&F zvVShab>_4oS}e;2$6rLcPRg96UKe_iDb2|i+(O!2v$!<=1o9kK;%bZfQ6$vibX#7c z*AY4#ul)esi&o_(c8s756&0?puorc4id@g5P89E=$c>rZM6pWf_Q*dqJV(Q>MH9*InCaDCvllcTB7jd<36Ko72L&67fxn3 zSGy5BicVTtR&rBET+7i|KPloLdQ&m_QJ!&ugR zK?}vW)w4_Y?xkpz8EELB6y^SEE;4g|Nf;ipKzGHoC^2g-y41O!mCxUes%FU{6;(Ub zGFm~6?lM5;KkuRJl@_BAiJ6q-QWK;QZ9qMmp^XBwjHwG6v8+X|0eF%V#qJ$Z0fw%* zEEg;XRI)4B_Alf>VD{&1XHHFZc3nJjoiDELR*OSJL6cP?|4@`Tp+PB?dLZq4@>GzE z4=PB^BQIB1?6_?{|9d&skEDc#4&Q!8F#P~S8s zkoA`@B~CM!v)q|0$^7F%Y;9{8S(m$>$){+I+-?~*_~ z4kWRY2UAG+!!-8%>kQ(myPN%LlS9xfhc#$sljNp6_VEwfenT5 z$HTBsQw>~JtE1CiO^|ww#UArjNg&n+ z7xE6!bLDlR zf5}YIoUfbc5cg&BC%zYWr(qM1i#kYdZ=DAd^t9wl`7n`LrVY2aGfIBFB7wF9D9UT} z^=aSZWbXBBER8k~=CQ?3`DXLoXrxsGb|((=^r%+WcsQAc_cWqPXBcZ}%j0*iFTlh_ z9`Y47R#58DK?lB95%s-J6s1+_$;m@8-<$H1JD+nF1@+ZY-EX&HYQC$eOkpTC7B-{0 z38TR}9=uf)L~FVjdaoVM=f)nCCtunp-)((Kgl`sz%9s)lIHoKT8C25@dm`H7r(W#$ z;}_J+xtC`)$MTjYH@WJ}5^7?yh3}OJMX=D0J1l5|MKe|9Z9Wau`tc^ZD}5b3Fe(nd z>2DWiD(g{`@RZG zehTA>yFaBF+Rw>vfxG1$DHhzXeG?+FS-$_PGb_?L=TlC za6@l44OBbC^I8e_QBW00*&=$Ftl^!eoXWSiv$f}1uq{8WaSj!77kVT7+b<@DHnGdeN*C^UqG(Wj>j=^1M^7<^NY##TH>ZzVBp zJDS!SJ^^3fkGaaM&Ghb^PF6oNlm<^V<85xKG=0cE z)amGj_)jy@{iY}G)b%D?&IZt+l4ve^kViYt`tkz_i@9>aejKWvMSr%jr?+pE;FGL4 zR7tHQn4k*HEZ?#@ehw*jenXshm~)j^+AwS5cz$*;eNuS^g1uZt}nvD(5MI+ zi-zMbUE|oG4f;@UXNc@*pc$78iv{EHwV*IrpND^x521NYtY}OIADn&&!}c-iSUZM# zb^VA5I}=!xz7nW>5(ett@eF;pW9Nx0Ab6FJrWXs~xI?McW!75aO1}0vEmV_U zC-(KTP*|8l)Z1S(_lio=GH5@RFZF{*`wY3eu^l*vKV%Ln`=IKc63)o%09D$PF-(|%87-=4ssOSm*QaY1fgr;L!s!47Ma(g1fn~y zn6(vMB^om?ux-A|@H&1Zx~piw)9M8*xn(68bTWcnZR#fecaO71qh8Wlbf2|lRuFdl zHoNaS0k#(@qg3G)QME5)@p}u&Xp@a>J?kNEx_4OrIepl@?H;?RJOt)8oML)5M@UTT zCv4g15ioe%1NKaB7g2aUhWN%O5skpxQzkyDkcqsg(J|C6gja9vfz| zg5*pv!pP6;NQ`2pOv}=eRM{jlwf-vEt?Q|zII2?iot_0*vZAMW$(=f3yUQk7eX5jP zlRlM6YK{s+Ch4+~wGV_%3U#teJ-3AFfky0t;cKCZp*G9t>k-1nC2Y3XhlIDRWCjBc zB;b-JJ5jn?xU$EI$&#~J?UflUX&{zF$K3wqaxKAW#zgcDumQt!7I<~BEo_}^i>jW^(6-hIy~lb%f9h;hQ2Z3) zYNe?6A`ol}mY~;_&p|qWDVij&g0%%J@#OhPc=P#JIIwdwoNwEN#f!Is{(x!P-MdVAL;%Fk^B7jA<=IE!!fvTPnk$?+;iY^ViA>zim)eGAUp)?;sDBWzN6 zh`m2NfuEJyvAnnooaECzO^-c$x=uzCii!p$O;i*a z*=w!e-i1U%A*76Mn(i_bB}zs&b2>BioKAE8IOln0JIj&-Hwq$?0qon@ASjiur z-d*(<4oJjeu_P?a@=s3!;XA+dT}*NRi`J+A`yY7qrQd)5Z!z%Q?ho}q)qm$7YX3r& z_cK}_&3@effAbG|AHTE1oY;PxQ1}6b@6QMRDf$ckod5DYe%k-LxonnTuHYkopXyIb zDU=s4C=}`F70Yef%SCG+otKZ=dqFf?d0E~pyD17&y(1TDw}?V&+vUYUPeo-D`sBJd z-ikUS)rq;^0EC2tNujDHl9g(b@k55A2AM7iyfgxhN;M+lH700HhZ%W0$P#Hsk0C(; zV^QsVMEZA#Q3W^};C=@A$?P!N_ggmEFe3*=X66##n?IuP-wQ~6 z6Gh5w1u+wzL^a06B;D>b@^U*vI+m8BlDKoEYHu}Kb)uH2m(-(3zJW}7aT$3HyG}yP zenIm5+oWIXKKiZn0kK`$hBBYElbz*H&{q!K|<3G5B1?e$8uPs5V_ zw7Li-{>ssBYVuJ}hBH-7%0;Ogd}wUu3H0pn6gsrI42?YGP1opMMNR{}Xj5rBT6xEv z-tbWa&n1(o{e#Z{A|$j`-2h5&JJD0eEdUA}=*#sE(3g&Ae&Hk#+l%RI&JzyzInhgM z?qI7sfm-2U2(HzqA}j+Ja~ImxF&i{@%%=Ur=7NLILfXYt(ka5&1 z#|D-GcmLiytj% z+`{jf-Gqa?b9h}RbFT5c6u)pYn8I$bJK9Y-wNL8TUXArRl0b2v*yrC}!H4*kGKY#7IvROoSq4?_40AxSL3WsKr=ZW33jxm>a9 zs-Kv=i5B0!YNTxYULuyYCn&F$nAza4eL-^c3lHpcR zGCyk8G0wI^ThUV)C>GUMiS5n~Q%Zxr7r))KUE!tUEw*dcQ$(eAb6Vc_c(2XpIQx)B z-Z(Fo6U>^&`=5;94oJV?lS*cDS*K-uL**ha`q&D7c3>zsZi@rI`B`S<8djg__`ZY z)fX(n-I+Up(zq2$89%OO14MnN!J-#u(8!DBq;TH^rV}f`yT%Hc{hf<2Auo@`6yC%} zO`c4Dr3jmEN@Ts$VsXy+47PAn240<2%Y=S+u~;GCyOw%$xw6T4>azk?hP z(yi(&A2-8|&Mvwpq8s|U9@AXIZrB#vPJ^^MVYvJujqiU9ytsv$jPC%)mS3ntLMu#n zyiK!S-G?{BZqq@o_n~Ry1A5!K4G!tH(c7Q@47J8rXj#$?@Gz*MFBV>bBW4x!yGgZB zcJ@3?9rZJ$FS<;Rr`(75ABMAK)pZb4TTTqKljKKj0&(L~Y0rsVxVj^a-pf7-aWguoTSF?miHfJ&cI}2uuDfYXQ5^J7N}wYa#lTAz zPaV(3K}JX%JuBM+i3LH_H*z_2Oj}Nodk9qZ&!<;QmcWs9OQ=y&5RAMTOuv)_gT~%K zn&RpYV=bpq!!aIUe!-Up`Z>et&t2)_^|o--+m;S>v4Fd7W>n`E0#s*BPjeg$QjVr3 z)`oCWdn66390pY>L+G||-=eoVZ6wsQ6Q!_)Hr^@v@qDN48|BA!WIqKLjexpeIbQL9R&JofX2qLAch zIn^r`)#2WP-i874YMU&(1s4wDnkPW)GU<)~L zWdpjlU@Mu@n2h{g^2nX0Jha~H6v;O#Lx5VI^6)md`zhH z)M(gs&YE^b09YLZs^8=ZekwERJhiFdTQ-MEZUnC~&HdxvwPvsN-JkF0`~B|E{_M53%0sCg z?}=mie*}-m)0VAtdJVj(dTvb5%1%dj$0`1|qS&_VFM2%Q?{|3|*y11*N*&LfwsEYf zUgr3_r`KCeb&8`+gy?M?8Vc~L1ef=ra|8q_8x6y<5$GqpC_;&yA=Aqs^Z{8dJ z-qyFy4_h40vt=eBs=WB1dDsaNy6$#U?Mq8W;P!Me!>`%9eWJ?6n}#3 za*G9-193$2yLfQ()Hg&iBMCT4lZlFRDhLfbPu}rQ1N(KQ#C#|N`0Ttyeon{+k^IYK zu^|Bo9gH{}zY6j?u91(Gu7k*=Lb4&X1l(R!POL@cV5sO0IghJ>(@-r5FRKS-22JEZ zNi%5O)=uhDJHW(0y2(k!9*|+%OQNnm10u`kq(pZV-1~ftw6AytOfa8Lte1zq>!;J4 zj})QEQJJdytHPzT)u^OW9X1ST(ajQVIHI9LyUpi8fByv(>|X>mjh4{!dPZ>Ed?^h` zHHV!!0;^b$dhxP@{!q9Jo!hw@_TMdqywo7sq?RU1c1WkU_q9RwUr*BZb<@zQz(^{y zPlv0cgGrXZ0o2azAqTrFVA^C6r7rrYbNC_kKdpm|y%$i+dSx%P0_zF7c`wQLT zorucoifGsNL?m{OqV~rRAtkApeNpL$eD}Dp^QARtLxKgF*nR~GqK#n0Q-8Fv-A*=P zV2;|mLfEBW=i}2~EkwQ*&iH0{F^lN8&Q(v;Mk+f4n8;s+i+WJQ)(0-(4*29T@dHh6 zNwNXDxV?`Ve(3={qTM;YSK%m7ZZCJM_#LkI=QEsqN(3%y{eoqrhoQ4+?yO3^9hJrE zGl{1@zVV(CtCIVaD>$~1OYTs|yPNl6Q)-XRm>5gHRYB?VY_N9TJfzSZ!U|1a|Ea`q zCMz`&Xt~%ipPS7H^0l~KcP3z*)lD>ZB89~&&Be_@=^Xy-7p8k{GhVmn43g+AMj3lY z;Pg_A#trL9kM4a zxK>3LS0B%AiJF0p-$U9k0l1csWUkeFI!7gAxM)!?yQ-;*GkY>HICvPU9rZs(`8%lk6xmp7m#*PFO4 zlj%5fS1b4K=uK`%r5#qT3g+zEZSZx!I_{*E3fI_of-~yO!L@68neDCnZ0ATW_33xx z!WXLLnO{4MXQuIhaD5RD3|7PWR?CITHF9|5+G1=!k%W3GjoIvF1=KXI4DE>xWgX?y z(4aUUy(~V<>F&FW%Qo%f7V^tE!~Tbu^(JvOs;6*bNCf9ph;U@>2%EpL4Hp(vamC-f zkMj(TIis@OI3X&Y+Z||(hXwC)4#A1I{j3F-H{%4t&qBG6`Y6_P8Q~P18*$5wF4T|{ z%~h_HA;XS)?AM+)G`eynvl{M#g#shma6pcQoma#eay{4q`1fV3Vi!T2q5&g`A9ma@Rn5)9il$DuT zO}9XQqdjPzX2=4w3!p=IKDG7UOsD6~f?31u78V`WY(WwOOX4SJ;`mh%#Z z>gIID_o_5_@(>V*RntP1c;Ij8MP;6qKqpCw@iPU;tEH9hR=6Z+@Es(^_kDrw5gjuB zvN+ej`Z$@y1n{I zYSSDzYqtU&I${RTtJus!nxuw4kT=B@!VLLLJ>#oU09-jTSS+yEcks3l)7S_8cdBS;9G4Nq?iB*SqI@Jjm; zQV_KPiZ;neiH|9qIhaTez3U90jMkIL5C_ocexyBrL<&hFit9%Y#9Fq`bD4Dn;0Ot3X5aXtEf~H0uk;I4u;nA3o3?^{O!{qtg zaKYKx#j>13M=TzuH_H5;^<5tBt(Toq(^}ROnJqg-eR8WmFp!z}ZI{)Y<#OvUY?k#M z9=2#3Op%415D1FHw#p*bKQ{LlX3NG)`eh~S@5_3S9=I`Y84>ZFfQOj}X|??j@P6XwONH9o&w*KY)v3I%5$wuaKr8tG&UG=S zZ8|m(e~oC#uoL`vx*gTow+?nqcco7bZ-rSJp0swUFN}ETMXPj0uy1}KeO_|_=3n@P z7H5lLLZFyx7$1cVUPowf+zB}3cY;QyOW-DX33Wf846mF|rkW0?VNCm3y4F|8`|TXl+zYd-vXfBhPE7Z%i%RE`C67^BbYrsa9%~ z*Z~K>@1VP7kD=V%UOJ)C2kW2wLLUS_hl^jnqzTe-XgO;d3sj$mMqLzG&<+JO$5e%- zUQ|KyU#YRJH)f#%U2WDkqJ_51(_yZwb)G48@rzteD>1U zmF@TFncH8p_pndn3b!xnOSXS`I@kW(p;CL@;%a-5yLI*k3=Q`3Kspl0XF6PJuXZNa temhXzZlHP-2M>E$eRKPKpgC5SeD;rl=6VF}uoJyoXIBWcZzjWmtpIKBQTPA= literal 0 HcmV?d00001 diff --git a/test/_correct_full_wave/_Ref/syn_ts/tR.sac b/test/_correct_full_wave/_Ref/syn_ts/tR.sac new file mode 100644 index 0000000000000000000000000000000000000000..f3121ef4b6996d76d99257235c954d192d8c6131 GIT binary patch literal 3032 zcmb_adsL3u8vjIhq;gkZ_nN+U@Ar)^6QK*0R7#&zN@1cir3*qvI!!K37e!K%NQKV3 zFK^1orP~?PWQK}_uAk6lOdS#!*g)<9CHgxUENXC{e6w{kI_@~+r0F5{J8&T^Jp)T=tF-W z>w86&O3v_U{4+}=`37h5pXcXFZXCYMi)SJUBQE9@bS_AOmGXE7v`V~lbNIVeJ(A4l zsr=BlgOcc`V*Yt`heW@sns>S)lRSBQm+$khm1wRS;A1|!E3wz^=CkG6CCgGD@gW&SITIO zqnL4Du8C;rbSCQQc(jvYnH}qAAr#6n`^|LFvmz@dW&omq6Pc=;C8&P6enq6D-z zY!~CA{}oCa+{3JJIf!iIlb9vtS!m7Q4l?Isi_qaoIZVKc)99GyDW+Reg3_bUGy6QN zkoAiiW`SQlGK_3t+;V$R?9PXb8Y71lIC=1R`xGoOGzSG^ZW0eIYwE9^MYkqgqFeSn zqx%Zv1=K4RdS@C5?$(QiISz9LzRg#N(^?@I`M3xTFPsIxVSB;7<`bbH-aweFJyB?z zr6@Gk^w7Mgm2|`6GFlgvLVYcE(G4eTs9Hn^*`<4uXlxz^oS~RCOJC3G>4mW;-6Gk( zCGISWHDOIPi{NR~BW%nsmFzd+f!uB#ea_}-30w61I;T2nx$;Jl7Q3pm zBF$~Icj_F9emp_z1J2S=k889$=nnO8?4o{RZNX1dDrk9|2^u=f1mT-dL91<*P&Di& zn2G{~F!fMjyGE*T+@@F<$Y>WT+pC1Qh$_K5Z>Qj4vQ`-A%Mdb#i-o589HHOUNs!Hk zg7MJ-+TSUr#&2Wkruykb^_3s;Gdl-Wls{!R$Jnq9rG;eLQhOmo&r)D3qgl=-nQiWU z#-d>@F6X;FoML?z=OZrW=AJ&t#hnf3JhBGZy86@XbcrK7)is%2J>>!G?0<%J^Y>M7+Udd@NSGb5DPjjD}oa7$NJ;s%lt>tdXhj1#5PF!AL40oZ|iY@gnp_>Bqgkb~fwB|j(%_U0XQj%5#3{=AvQ<(1K%rt>r>yoZj-eE)~_chPD6M)bh3 zXqxX9L6Q_yNVX&x8qIjalsvSD%laJ1lI|%)>Ev$7h;st1&$FOm%1i0eG*8-D#M830 zP4v~)xwPLuljw|BqObp2K!P?-gt|dz8K>9#p}X~1)L)^AH8=5aV^b*Xx*P$+2X$Gz z{4K7L<&y_LsZ*_Ksqng53sifu4K|rXLT`(5SofwMUOUnUdk>hg!6&am#bh4(`{jVS z73&!nM1~7qhaIrwOgXSN{2ot`T*j_|6*5A%mu(y`ZP1#_5qML-r)jT|#!u@FZ#%){v|Fh-V0-F7&>lI+6x4&{hB0unRSRh8ehNT)8xW1_0f!r_!Ntg% z!1l#S5IesG$Tb9mYvXSN=u!YwAKd|#+8V&);UQe-dllu9S=gt%5KXa7!J<{#I5oNj zo!O*-t#8<&#)Kl|E7d~gtt`Lu+;yg9@R`KNt&^#?8q2%<)WVeSuw|T{jxde#=8V3B z7<7FiF#0VE!PFBIKwotMKx5oNU`0Q0s3-*nN)w=)gA}-bF&VE*lk$c)gV4w&8DC?| z^RICx>S?M*#p7pUj|IiZwRkp;OOeNu)F7^A@=)TC74GaQMf#=ocr9}Yr8_%gSbe_wp=b{Bf9u@W1nVYIf)3m>RIg#vF$@tIr4P^>fnuY9{3)fu_sza%B1 zjk)Wv_17nmxq&~Hb)83XRy%NOZ7*tV*oVa{22tbDR2=GGfYOtbaMp%1DBmI;oAs?g zp6-!&OM@=j>$(mvi+jb(ngel_?mb>{#!AdPygKEcv=Dn#FJkUA8{oE^8I0PR@%X;Z z5ED4B2pPX%!FmTCW#@bfS~DDxMzbR@J-QZoJK2E+Q+xTB#?b&pE=AMU#{lDsZ78I0 zJp$I0si~CxPJ4X5%)~ RY)~VohpX^GV8dDBe*&_kg)9I7 literal 0 HcmV?d00001 diff --git a/test/_correct_full_wave/_Ref/syn_ts/tT.sac b/test/_correct_full_wave/_Ref/syn_ts/tT.sac new file mode 100644 index 0000000000000000000000000000000000000000..b8c7478eb0ba9d424d9eba9f7ab06a8b22e12536 GIT binary patch literal 3032 zcmb_ZdpK6v8vn?xkV|qag%MeMt^IxZ)aTZRQc4KrRx?qGRHtc_K}JkBBczgB$+dCW zYwf*tk%Wc{MMxx@BO{+yMF7po;9ws&`!|n z(C{B25D2t5d&l=6_)urRduV^x!Oh_V|64K6w*D`A0>SThc|`a~#o|JTLs`}iEtjFg z@1EXoEtj($#A2~SaIongo&<*P{g%BCm$`q@`t*JO18-~_9sO^i@c!-(bwAmE=O1eS zLYY5jls}sNxc~p=AM`#xhh|uaEoT}Uen7)N*9ZP7dIGC^bU;7qoDOfi|_RP#5=P>Ql6t0bW#VOr^gy#9Qp`sP7MlF=1KFl%8c2m3YDp zR(v}P8XP9DIJg7`naiNWWisgJP9f@7ut0+TKxFD2fK(ssMllvU(KzcY6blZb9{*w_ zzLSgC_A=D_Z4oNAKZUksl%g)#)5!LA2?CqSk&MDo#GGf*Sldjb9*~DjRQI4@n^bf{ zGaBu+PC(92gVFYp2sF7!ihMI2(YtU}l(h3Rv@+-o+__N)*q zMz_HjXH_)G_8v6Ro`hTk0wiaJkXe})N>p$~E|+XjYw9K>s~Chllv2=3&v+C$y%2>q zW~1hXXHeU%N_6ta9u$)P5XlN<>5lyhbhzO-+QM=Yy|^G6>6Y}e0^K!SkyRD9mw&>A zzE zk=R7o;iW7j`HDiLhdp%FyX$mH(?MF|Scm#hE;OT@qo{p^MCy%u#v(hswvwx{scKh}_Ug6n2{{Dt*e3tcn-nP7g zPko!kTYoo$f0eL=_x9;Sal16iP6jITl5`t-wQ)UtZ|^E$pL>8XeO7^xUVBXVO~a2u z@t()Rl?7^|%zS;(&f{jH+$qMQBXt@g?`UOFr|T^tIoK-9iokr;8yFF`tI(|$I?_ADwO?QHP(vyh~vgq?II3ArYwvm0|6)M3)Z zj=yV)>a5Dyo7cwBMJFq{Bat&XgWaucWcy7#@`WbuG%ZIj^%GHe`UKiz(`CdtPp60P z>(l!#%%z7%3sBf$M>>U8rR!@vXwf<$y{l*)y+%5TR<8=78%pX>#r)~CPM|&qQW!Yo za+)Nn6{4x@##5rl7ty(+MI_Ap9nzb0j0=(3#o1PxacxUfxwS51xV#+~*+%JMRn$xS8%dNMT<6I53x!P6jY_g{!x6!(SeJ=l)ts8%w-NLo9zHuAa-p$w8 zz_GFHW$zqz=Ax_E=(l%7)uj|jN7~2@3uU-awuqE{ErClmN=fACT1b8ALX3kGU1M?f$74yh@X#ku2jSdS+cu-rsVR`HDwsn)V$ z<<%Hs8{*4qW!ICh9pl;45ivM1Dvu~;5D+*Qlk+5osu(>)>gIG%g&ud-)$rByY6T1~^y-*@)Bn-MKIFNZq zAv|FuCjONG8i)5{yPznzofKoc@S|{Ec87W11{S&>6q%bAe-B${Wim@LbKuvS@r>Ew zJUG(W%p7=}1>MUn@QI6C;TDg4?5Ys}DW1cQR~AF{R*b_2m%_MjnsDsuz3{YoEaS5) z8Cv!(#3P0wuyx~A+<#>%y#B5Nzqq9d<<0W(aPwX2fLkI6_t!)<`5v%4PXSFk?hGea z%An3JFX(GC0;P5CFm1{+7;+ZCd9S*mRC$DnFx;@S_~YsT%aVW80Qs9 zphE3doLcpOigi*Y0)sWwRIfPVtXzU;OY(`K3nV>$=Scrg2S^jDCjJets7S#=((U#Z zv`KtOhMFstV`xso%~GiGL$aiA?pZ3c?^hfb{ftVAHXz+n6-rUnfSf)eq1ue(NX9QM zz|=yWkO5svko6jGU9V3~+%Sq)zZO$l-t^%qVeZs}glb%SIhBgAO2JXEh4L4ChIO@@ zD3R+)CThPo)yTwP^^}Q}$iD&qAecy%ojQV>TbwET&NO^(Unn)i#$%(Av%n@G8z&CQ z;EC14SaYTyV?4hV?;BdcqBj}G|6&8uO0HqiW*c04qYSHF zZDbY(Y{MomD;c$139f(pIpZ5u!8{TkWrn|8fxGt9^~s$MRXQPS1{1;!L`8I(}o|8N7z6Q@(FUWb5H#$Qo)zbgVM3fn2CD-B?( zYcg0mO^V$gb%T2%b8NZ8kV2+T_-DiE)b%$J@BhJsit e?<}MIgvD4XIgaY?%)tXKsZ_>Kh1hRGG4)R;yJ7kO literal 0 HcmV?d00001 diff --git a/test/_correct_full_wave/_Ref/syn_ts/tZ.sac b/test/_correct_full_wave/_Ref/syn_ts/tZ.sac new file mode 100644 index 0000000000000000000000000000000000000000..bab873d185aea322f0e665121f3511850c29d5e2 GIT binary patch literal 3032 zcmb_ddpMTa7XPR|rJ>uWD3@=H;_ z_kH(%UlqkfriL2#LQzgdkyN9baT%)9nVFu`Y0e+#JkRXsS^KwsYpuQ3Z~e8^8&6nE z+7ws*!z2=kOt7_k2a@-3y!WAPi(;wbJ^x!V!m8{q{Unm#XE_DC*_fN3P#jL3si>}^ z6~AkGcU4!VDa_5y6-uS?ADSfccY1T**`m3B>Gkn`{{z2T`}XaBGwHk8@8io=|Lwnz z{R>t8tQq-1?T7jQ*MHyd!*^h+nYqO@x%@rK|2*&c$LJ*aqu%*X{5b!2^DQeRk`MKL z++X;?f^B%)PdKEnWWWBolVGo3&1p9^CnM-d?3(VeBrKWB2AlT~Wp~oqe2=Gu`|B+B zm2W%ITz`VqTv|^Q&MRaSJO$#%v9+w?XM%Y1g9VqS`qbn|jyu<4^qX1m!$!{hZ5?{`^a3Wz*O$r%)ClKXS|*bN&kNwbDZHaTknSzw>ae%I+n`dikq(LU?Q^u|Z5 zp}RksvGP;yd5R|~^)=+s=l*2WdO2rva4qQ^;KXfG*+a@UGu%rMMP_u);X3tW$c;x` zIG^BH(#7X1F0wa)3>h8D(U&qwjo#gyJPVONFLramAtj{SgM-}J*X86q<1DVrt&H4L zQN*E{x5yQY#_*rDPcpkvj-CE+cTqGXF7AYF zc|NdTb1kets1D64WYA7w0_S$`hXZ%*Vd%C8c>9lSILymXAbjWHf=5r0O6W~0hApMN z{fu#cKshb#7l4-OL-2Ax1A*lt_#`_DT24-e63Jm$d2>IU%csJFDv7W;z7s0H_!W-n zH)MLJ7&B{*Su+kufsJXG&^PnB8*BlzXv^x_J!KYqXeI zr?0|n-{l9TovO$=?IzY&lL>v(f(4KCI-xvstk~o6ndmz`MwIoGiJo!S#hI&*i2GE} zh-NpFL_+t382C6$lrM-9g9~)U$lnshH1`Ve@#LrCwE9HRVwAV2EEp>M zuGD41z8TL*SJc6~>1B`}bcfp)hC`zzKG0b|2WIqMz_r5#C`ng`)3vU_^|R)~U#6WD zY|1?aS&=nV>2!s_qzra+{jFp)S-yj+tcE?zz?cPf8WjB8si+NnB; z$)gsFEe2Miythgmbi5{}j$bXll12$(ofYu1rpO$5Ho??!)L7G&SeYqrJ;HRDT*LJ3 zEvcy=e97Dp2bnp!gN%b_6{B`Om+@+Q%>;Us!u-!1nV1(hV2Xh)qe$onI+_E*h%u>R zQ2IKtGs#d)-@aAEdyb3MJ6Ms`CUYk!6}cN9P0b~o^>ULg49RtfND+ybXN8zu?n0Vy z6cA@NAZmxIpp@0okvXm)QE7lL`4F&mTPZHsvP2+Cr!b}s%fvMQ<$}F^iIBW`kzn&Y zSFo*35SpvY1U*A*cr8`}1F8zK`CO4cE7X9XS$jde_i%yHQWwH+?!{|lv*A3W4s`j9 zLP+($g-3P|;+{rAwET4><1ld?&N#A&|9<^M;Zl4&o}+z|(}{M#x!OxW<(+zrrZxd5 z7blQ`qL9Z=^8~-}R3S%*6egh_pr+#pj~pvT#H?z1chP+i?(Pj2_iw=WT4jhlOwwn= zHo#}F8@Q*u1+MRQ+p?e!ySPjc=d877j+PIHpBuzLo$0s1lpn&Nj+g*z1`A-F9uFII zkHIsy8=#YVEA+L>hn){9Ve04vxasp!7*QSzo!w7C$|ey8?uv(#CPu*WH7u+v_Jeoy zJ>k|aW7sW~!B4h~ggqzJVWmSCxZF_+JV#ssLpL7)ni36=tTqFD8{RV z+iB2Ua}BWv?P=%8W)%3)j5fP}81X>C^*0?tL+l4ptnOqSd-*hS$=-qEw++Ml)xX0T zWwUW$QxXn|^TWrdW4uxsjRQ~RVtwuazTtWfm&PSwR_hYh-I<0{+H3HgUrymt`T`Dg zO2TsS91gjdf-Rnm6vkG*#bX|62nslWZ^x<&7g8SMtOc*}^wuUUU3C{*ITYjkvPK;2 zl8FZ$i?CO18a`uj3Ky=5#C*yw?3g+mAGcnC^<&55jm>iG$dw?A$*<7hLUYtrn1&?& zLQcg%3scAmN+=NEH@<`$|2aI&2Y$)6oU7 z$#$TvV+=5Avj)X*9H{$tA=on!f-W^XFyXuoxbko!X!2i6-wJq28+z`d6Lnf?Hw%(J zf4!VG%>9ObzN&;)bzDWOmtUmYf3>AOeJ;_PrEWA)eVUH?(VUj_B+;=YQd-4vA06KJ z?gQBO4eh1aMokNxP4ma|$&N4+x-z~8K^LuK#wyK@@b(_eeel1;;tonHU%NOz4OtX zhddJASfkn?9u-$Qp-o@sqjaw?5lAaQbvIm4f!$dYm1mE}>13mmhiLSC)L|qY&_Q2B z9YSN0G*M+}6wdoBAU(>Wg*->-( zig^qQtsc+2XzQUjm*iaCBQ+$okdSd(y7*6mwvvXIpYU$C)5viTT6n$UAymj#fzNVq zruuO|$(8Nwr?0xG?_41lM%|-!71?s31#y&z+E#vY z^Kr^2YXRS_$x?d;9QY-7GpP!9TRyudi(+Oie<>p zE;Fc+dhzxkS0`(y(wVEQgdL>14%{QQoYbI?>6lusAfM2UlVwp`X`XChx7md literal 0 HcmV?d00001 diff --git a/test/_correct_full_wave/_Ref/syn_ts/zR.sac b/test/_correct_full_wave/_Ref/syn_ts/zR.sac new file mode 100644 index 0000000000000000000000000000000000000000..e39fed87ff9d8cc63fe648cabe5e46db251b9fc3 GIT binary patch literal 3032 zcmb_Xdsxlc8s6Pi}O-q#H$4sA*(7&WxVpnDfUu&ok?JzV&|Z_kHXA-lbN1&Ro!9UGyIz z5D2t*t9kE1Fj&v}uFa~^dXe=Y|64K8tne>-0>K}?-)<#A~*52OJ|MVmfzV}PthrP{Tv_Adc|Gce9AtXtuFXIEcazXOMr2PJ%z@T|VHa{lA;7H zw1vN$7RCjw3gDy19^ev|?d72~gtJr%z`o@uC;;TvSH{FYr6T;htl>sd*Sz z{QXIOmr^80&z_Ctagc0V}J<(FW6%%&sd1SM2d-QQEph z%WOYIX77mTvZEv-u#h!mzJP=SQMf9@oLpNR$kO5s$bvB<8WERI=21;L-ccUR1YV7!h zMI`5KAzS@2j7^?qES5zU;y&|D;4|A3xI-t!Zq42}SX@F{9p_6d*V&W2ma}yBdNUFm zT0#R}mk^iuqpZTFjBO7B@^XGYlV55n4*B&gsjb?`y2Ku&qFo=wcFka*Qb21jHBsUh zi*rv@(DbflaOB`oW)>-wgec4qyVrQLDY^25zN*HP^IS;olu}4rYfckN?^A(Tg}9Fi zWqMx)K&+jytbli-xv zP?DQQM{8V)0=_Y!H<&Uj>|t3X2hZ(^t9 z08+X{h0U=Ru~vQ=>3D0*?v8wdOJhP=dqgSi`RoCV?{7#UT?sdM#evCb1?;JEnNLs1K;@A~AU$zBp323-0(BwYT&6?4 z`lsP@0f%7WTLF&W=1zU>zQv7R-|}6$SK*L*GN_C@&uwk2g2(!L+~I4X@T=NozGu`O z7)Z}T)uq3K%Toy|`f?J@t{ytXe?JBy znG^{}eh){797g4Sr$K)5EwpBMJVYtyqmr-B!I8yfNb^Dl*m^V~x9b24O!81hY9oXm zxrEfOq(PhAE#xxvDom3VpiogBv<#)_h4=wbqjGd9E*DyIn$gJjsn`rR8`BVq>V|cWdFZ(FBM9cCD6-%w_>GsM&rR;b z=7qPBeQFNqMUBOqE{>qP_qL)Vxyp3(?mi@L?T3{qeQ4l^x8SEJkE7h0!0=iJva?mB z7ry8~6ZB*tZS6rBb6&&k)6K|!p$yz3pCa8i@-%eq8x+?tgpSVdLlIk)=`qb`C{>_D u?>v5h0(%Ew)-3_Pe@c^D_s_*S4`$G&W*uCSF_zkz55vhK9qLh|iT?@HOQ?we literal 0 HcmV?d00001 diff --git a/test/_correct_full_wave/_Ref/syn_ts/zT.sac b/test/_correct_full_wave/_Ref/syn_ts/zT.sac new file mode 100644 index 0000000000000000000000000000000000000000..e6b782d1ff4a7eeb88d6e36f75cb3cd2ef87c219 GIT binary patch literal 3032 zcmb_Xc|6wH8h`CsbCrx386{cXbI$LzWh{}UFhnVYy3@i)rDSYPC|jZ#Ba|hTM5*_j z^LuGVk)g;^j4WBw9#W*@W+}Senen;Z=KgW-=QHQ?dCvFwKF@Q$&-X}moHpY;vAq5t z!QpV`;$>DJfHPC)AM9K9z|zTbhX17)ZCdjuJr3u$T`r;PEcEqHTV`dNST^3&v;57| zhpX{sj-|f7zNN44;@>@S)Ia#8AHu=VEX-N~IT6U49Qq;oTu9^~sq`g5ymPVwIqq4x!aMrTDly;lct-Y>!0TNZ-s1-iHd>i~&v zef%&OfID0x?AXf%P3qdX{Lxb2z0nRks91w?n^ib*)mNZ$i!1I`@Bt6^Y{a>8UxW76 z%{X`adT_Bg2(NAr0KzG|@f^-RV6r|2%l;e#5?lgs-nwWoaA-g7fk~je_8fj4k_Utx zb8tm@F6do(603RV0e#0)*yLC)NOj4_o6ISoi3IyRV?nAb#l9m1?610ri>6CKyqnqUyc~49GbV ziFxY4j_*q9#gu^CmX&1JUSm+(%q2GMg}^6Vl}P{a04$43Hli?B){yLy5ruJkbV$t(duUm89d9Z$gy!YH;t_?_&^O{94j6?n)$1A# z-Kq?0lq#@KW)*nZX2MU4Zv!bj6=;+-U+=Y&3#=7B%jB(G&&1nqW0r4Q%!t>?GtmZ( zXhl&f@{zYe78}oDv&J_5MTIPibj!)%yIj`1bt}q_Un?kUR^+WwbYs<^y1+g{i{#c+ zP?yU|FvqG1`FyF2ytB8nwHb4Gy*Km)ineKVlAFgy|I$FhK8t2|MI|%BlYZ!o>N>`A zQ#|TwfsB~+2BZU(7!~ge=)PJEd@d`3`jxv-{;c&dm?kyX(*Y$4;%MDiV3^EpmT zR~uwkUpHrti81W(!9Yf(e+B!T7s0$pTEcE)^ckVxAU0u-4m0i^&1!Ak&lpsU&?THS z-u2DP*b|zmJTRsz(DZo36aDVHzRtyt_s;P`C$msq`Yu z-58<4{^Tnp=zA2;z?Mcjb@3v48uSJ4&Ia4!8ztZHzAL|t);_mFULL0zOXd`NXUPdB zVNo$#G-<*-o4c6}8>wMB+cMY{Y6;AM-$AygMwbb=X3geTcrfpZOxRevSZ1$OF&itC z!SF5@voXDq==qpbEGH3yR%RWB8dEn=oRcS-8hC{&3oTLZRgBEq$5FQX z7GzYNtk+d-#Q2@jph02t(A2V8B<9FRKl&|1Eh-1lqRC7&@-l#JT3XIl_BSv?6^Iq{ zPiG9A4A_GnY0RFnAKB(Za=ab&W&)b3!QLw^g$B<$>7gJu?xBY)t=@N^`d(fE`;-~l zHh+**gjgUs2Wx7rd4m3MK@y5;Cy=D(^Ej*O00>CXB}&)r(Iux?n(8=>VnPqFtp)Y$ zF6CoPV68Tr8EniLHx^o?92w zdCD=UOSh1kScjr7e?3OiBbTAt2tR66D~kHw9iVZ~-O=;NR5~Fpk303o!TPlA5C#s*1DI#PxDl_zI<_>96*P9&F96Ln{W#!N>$Aziz-C7v%JwwU{tD(>5Ur~E0MP$<^MTaKmA<<4*ntHkt2Gx`i zOX*v1|NC;{YjXtyQkQcAey{0FE{7#ThtTch zq(TrJ@-ij7d>@#sB26A<1Y?+3jzfYL!-@6hG1X26JL-yX@F_dsIjO^6Ih6!HbE4Sx zv?+EdeHNk<0jbgE? z#&Hli=zyEVRe|BrPQDfG$c;&yJVol Slk@nyx=}F9YsXwV3;q){7k91z literal 0 HcmV?d00001 diff --git a/test/_correct_full_wave/_Ref/syn_ts/zZ.sac b/test/_correct_full_wave/_Ref/syn_ts/zZ.sac new file mode 100644 index 0000000000000000000000000000000000000000..9c5bfde37330c576e3a680259f4a92d5a5e6fab0 GIT binary patch literal 3032 zcmb_bc~lMf9&aHjEn4N?a*-n4-*RtBQ*tvHBN9RrYQ|Pbno+7T(qv0}n8+wC+T8AM zxo9y_j5HA`C5a?kX(7w-Cd}*2jPtyi_s4tZ%?X$!heX0 zii#1L>N)@_@BNDChf~{Rb7k-CzZD~$D*wW#qViiUkI=;~QfapAP^zQsa*eI*x1I*B z%Qc5(QmIt7W{uVFJ*kKW`eg$$`{TcG{dvCs!LF+P_1FI<8UwZ8`@Piut-tsE3swKf z(f%O!!}$NJzvuh#ofu1`cJ?CCdlda~-rIjh50&5Z7yf~Nj{n{E*WN1LDj)Rs7yb$i z8Z)J&zbFOt&5R^>{P%#IciV`~pwGd$+(_bDv;n}Wk>snqQ1I}=mt@_h2yo4CJsI}Z zm!K)dkG#p)3<@fPh|j|SkT7E*c^AQR)(#r(m{Cv_Y`VRE0%p^{h zsbJ-cOyY4T4al$NlIjUrAZ8$*n|vHJJ}f3D>$Abd;iaVQRz67gSCXTS31C5L1F>KU^hR3B5pW+(!n# zc>zvW_K})*ufclTSLB^hI}kooqnE-p2F6mO#rf)Rv9BgoL})>`DlH1;YQyc%G^u|4 z8(ty`(uF_Ox2k%RLiek(l|u~`gF$9@aJ_K{CvT0ysBo#!*?D7Z#LzB+AI3WukMF7!g;{;5pSS&$||`ZX!@H{u7ol zk3k?RV1@y?bmy`k(FVh8m{&d;dp*0)bdOY1Sg&Ao(PV5PFhi^p@zkCwWE8xLUkUfZ zd(CBhli7Kswrw^q8{dY_Tmd61reKqAa<~v}#`WEYxsO2!oia0ueWv;&>1{JbTS{_J zLtP;34T*(Y2dVME!BW(3upL>49mGyw-(c&{o@8^j?BmBe7V{~KHZiMHCfMDsj)!Jy z;g^FlID4Chs~kG{`XmR5`n*WVzLi{3w&rVPMwO;=IkDvPGS=`whh#1{e87X&&*Rl^ zG4Efy4SD|*%_1@kB(dj=`L{N1+_~@yuQ(@U>r8pCwLf2f zzyTXPS-~~Rh9jF?AGkDt;G+j~n8UnXxZO)3iFzI~= zztyec+A+Sg*=i$y{_ZAT)~Ck>_s2`-7V6T+acOkdg+a;@=^+SdMj%_!2Huw4jD@#X zacPtfwhVUWt~=t$$+Kya%$iT}iCX~2IBikhGCG5Mggdx-Y9A7=kSQY)OsQFewNkr# z4vN1poR3jkh&P9d@M=AOWw-xk{-`|yH7!Wv!TQF$#iR#zr|yCGE(Ky)y%rxb(n^vV zu$_NmZX_ulW+Zt&*c3awIFG!FG}!cEL;~0M;;zsl+_^Ab`B2b~KXwwKFAf%AT@6*` zp=c9aT-nbB%4vA{J1;!eZWJz1zD7?Ygvx}I_WZQ*US51Dp0Dl+=bsGh*e{@wJ*d>j z!tw*$9%Ny~lw_`L^9Y5{7%AymSIEuO67b$mQ(p9BKemQBe5NMCDsm+^wd;k=`dZw% z&Vp`kjbl^95p-nBP~5(49a%nU2`>pkIKyN=-y-tB?z2PjopEvWcFs!{qT9uLtfbfv z?1U|`T~IhbfrX;)VZzfk5|o13Ou=lde7TPq@nLvEh7-G?+J>z?68ZTU9eiI}!4_?t z&5dq$q1cLPY<1@lAPUjpf^j>sd1^N69I~C|bjwjmxF>V%J4>>z=U_KzjYcmhWoAu% z++&0rHV`DSg#A0vTi4rc?BF0|V=6+AmYiWPhh?FTt=Vj~&`-R*%#_s)J%VbUhB1o^ z`$)IR3(|JflzCpcM4i>^Vct)(SaaBN*5!5@1?qBIpr44m?rN~Km+28Y!H!hjb2F5@l+&;iov>2t6cw_5;8-6`#4=UrmsCktjd6yB$L`YU zu6a-rew0cdZG-01blAzn6A)bQqK`e7LjQ?b)IR7UkZEeOF|z{Tr{~q#uKh|l{rnC3 zMe;}R=gIqM?cG4=5Us^>H$8#(qGGAgbP_$1yN8A!mK8qC^Q7TCTO86fnZ_jDRwTxS z($mQgz>_&i)Fo*iw2d0_-G3h+L0sY|)t39+gIsm?296)Q^^`6e<#u z>A6L>pzZQBS|XT*S_c01lChrf@cewbtgHjZ549rlXmz?{sW%A`To=<{UX$n58-Ny5 zqpKF5R_NYKC%WRxAf|mLndsgPq-*<$hlLPsSW!hTNN$5QQe!ev!y1n6wjoAciJ<9p zAi0*_0~(wS=@YXws8garzYhk;v<%bwibha7-IwUi$^jbDFKoAE>;lRXck!)QS75iY zPodvIK;yGXiig_^3)SVXjv3a@9@zVE@fIos0ySIlPM_7b#Yu<6%id}!UgtXjKhskR zW%@iY;!6p+v7Cz~$Db;~EIY)#zS?AHj1QQ8=00sLF@&+Y*+joJ2UOUPA_2Q2f&ArK z@|jK^2*2P!=Ei;xHgqHqxnm1xZfzjPzPkz@9aWHq6FP8z;7>#mW(1pd4W(~>)`0~D zTC^j5Bs?I}qwacUu(?2m4!daphdec=+8T~<%;>4qJZ~{v;$qI;xbK5e-HEj9un$~m z7fuH~^@X>jYw5ZSFE}hLkaoGa!U>hf>8|}LFyAGCR=KCZ+}cEHd3Qg|IvzzI-;9HH zXI9gTRwEH`cta<|G{899>vYVFbFftNGMy|bg`lmTR?8b;XLLI~TUrYTo7|(<>MlXm z(t3JgX+6Ah^ENH3tcTin>u7>T3GAPAh>G&<(H`yjtUg45S_5p@svKRUG2NUk5U8WE H$vW(xOmTtz literal 0 HcmV?d00001 diff --git a/test/_correct_full_wave/correctness.sh b/test/_correct_full_wave/correctness.sh index 73ea4d80..5143740d 100755 --- a/test/_correct_full_wave/correctness.sh +++ b/test/_correct_full_wave/correctness.sh @@ -24,12 +24,14 @@ python ../compare_sac.py "GRN/seafloor_5_1.1_100/*" "_Ref/seafloor_5_1.1_100 grt syn -GGRN/milrow_2_3_10 -A22 -S1e20 -Dt/0.2/0.2/0.4 -e -Osyn_ex grt syn -GGRN/milrow_2_3_10 -A22 -S1e20 -F2/-1/4 -Dt/0.1/0.3/0.6 -e -Osyn_sf grt syn -GGRN/milrow_2_3_10 -A22 -S1e20 -M77/88/111 -Dp/0.6 -e -Osyn_dc +grt syn -GGRN/milrow_2_3_10 -A22 -S1e20 -M77/88 -Dp/0.6 -e -Osyn_ts grt syn -GGRN/milrow_2_3_10 -A22 -S1e20 -T1/-2/-5/0.5/3/1.2 -Dr/3 -e -Osyn_mt # also check the time functions python ../compare_sac.py "syn_ex/*" "_Ref/syn_ex/*" python ../compare_sac.py "syn_sf/*" "_Ref/syn_sf/*" python ../compare_sac.py "syn_dc/*" "_Ref/syn_dc/*" +python ../compare_sac.py "syn_ts/*" "_Ref/syn_ts/*" python ../compare_sac.py "syn_mt/*" "_Ref/syn_mt/*" # rotate to ZNE diff --git a/test/_correct_static_disp/_Ref/stsyn_ts.nc b/test/_correct_static_disp/_Ref/stsyn_ts.nc new file mode 100644 index 0000000000000000000000000000000000000000..07165bfc429ead903b18bb8bd6a1c07c81d7fa9f GIT binary patch literal 68880 zcmb@tXH*qUw=TNKIp>^n&atb=AW1+#Oe7Ie0TDqF5D=AQ08tSU5fn2biX^eCkt`}A z2$EE?8@v1&7NJMm#LWr8R>7J-%iN4BURw>;1IvR zb)wq|2@mC?H}^+zra6y{)GP*oWMW0;QxqI`v({DA90p{aAZ#ZH1|L8JMfp^ zU*P}1@8n;8{}Gq?m*0QH-TQ~1^MB-5^q1dX;Qzqy<6nON5!d^d-+#m{u~yLjweO05 zZL;ED`>Xh4-P@sp{I4BXkpAL`&A(-0?XL|J%YSji{J-*Zuuk(Zd-Ut4BU--A3M)CM&>cu$GnrDR+XksCMAuw{EFjr3$;1Cv7+B%V>rFrZ1dQH`3Yp>e2d|T*ZFQ$e zXfR}2lq|jpEj%u!&zdEn&AAqOrMEr=1s&Xbw?mhpa?Ux`a8!|?wzI2u8Db}>DNJ^g zU?)MD9Tk3Njh!I(tmfz$H$hvhRMqFpQ_%A31#7M&^JrFV>BCrRH~JkB!nLQ;6Mffc zYxaCl3fVV4?C{mO1KEoA&h+MAgDk(;$c3pMLk9C2D}^d*NcG5(@^F4O*u?BblEsw4 z)We%*%_+LzaWZ4njnx=5Zd>$qnKKNnb&qo7r~4BWZ?ne5z0wF;jup(E>Jvd3)LYo`C(EK_7y?^J0+My*7e2Up@M&RyIKuZ02}mZj&I3KWWU}7KK)?Pm7xQ8=+a& zufy?=+))2IdwlvAPgLW-kiAu30lBs2nvm=|p@*=-KDj!Wh9}1|}6kmxN~PKY!GJ-b0Y{XP(KsJV4MLH=*u* zTuCrDTjg52I1{W33{c@ZKf&gJ#4UD?5UkgOj~(#1LNNWL2q+aFC+KzrI%I296BKha zM)on{X#G|6rK~tfH1nRVB>O}V`k9%1DV~u8-tf+W^h_hjx5ys4obe3u7AG4wvv))8 zSLT)MdlDeWr+{D7r#?fLMZK44d=ij4>Dwc=rZli*R1qn0S`cOJ=;&tvd0`7Ik>F?)tPgRqCAhNY^4DZz2(D%`fvnfd1ZOhWilk!@!7e`K z_{3w8VEVDBn*N@Vpng_uQ!7LpT&N4OiH5c?AWGc!FNgO+j9}v@w*HV47fG9k;I_xL7@cpJqNM_&muubIOVdfs13WCu3#_fso%{ zM|fQczQEL0FU|7=*SJws+Nni?Eox-xR9y~1chL840p~i}>R>v?C(Db*2XsXwT8+RX zBd7Br4kR!dw;?$4Y#54_=5j?Jrh#JnqgkBSc0$pU>Y;w}d!f*qr%mZ{qL6pQhby=5 z3}p9R3aENh0x9DK@Ikku;A6K}OVOwiK`x{4-DowNV9((0uBa;~`0pusha8C`M0pF2 zoi9`%#5|);?>F2e#6&;)xBMn zJJU)w_cH1a_I+x0QWW$AetYkGbs8+*Sj8SytU}2=SJW|63MKE2eWyFo3nk`Dk_EUc zp!n&So7zkMP-L9qb6Zn9bJD=`nNJ6$-?{WRMba8>b9?z+g2|jMX8wWm6J?}=gr!d{Prb;T^Iv*WHl1p z`uCgCIGG4W715;HYew788?}nKsug-QG*QSRRsrUMBD{?cc0;nvP|mlv?V&8aaQka% zAt)P>sS~$11Ik)Fpv+FRf-(}7vq@&*P_pXUyY9LGDEh<4l-DO5a$k?+YgT21q}wgJ zQs=)BB%v8&@$wjg_txj)DAJWA2iL3Kr@zClZXH%$4G2 z87_OFSkjx_C)4L4yUXQl^B*dxr(Kp2}>k{NhZ=JZC7Hi*6^x%x`E`_K6c*vK&XF%GwCB0mdP`mlL#T9f(u2 zEP&*^!0zrXU&s{2ylN{;583az(nLq4Mq1smx4`MqE)eQmK zpU(i*ru;=GKI=f$%hJVjm0?i9JIDS@YCaUV(9Ljrpbcql$ab-4FB8nQTG~Ok{Rm0- zjPF)ZBkPCR zNoIejWstJ(+V#^=J1v~ab8H#v>`Qp}J3SEUFkcKx661qf?9DCtYhh5eZS<7+G8vTg z`}9!fp$^#05W0OLoJ8=h>QT%^auTv4L)`;i_X$OF{=T%UBtlN(ruCO+orKtxSV7HW zJ_KuXg-!%rIJit&$=(W9ysnS;$xOkd}P z;N5bl@iSoE`nx&Q))`J@pi_mqmAYfJwk^Qx_aNhA8|vrOh?i6fgpm;&iqAHzbC;5e^`+okqFN8 z-E{}g@}i}JoZn`1L6GWq)s8+_2gqk$$o9kL7!-T6%b}zuA4=BhCBNAb3#HA+M59F- zq4bk>rzr1aD7~|BtDxf})H?Aau6FJK)IAm2-99q}4OZ*0gGz-EZVav7QrLhdDw;-0 zi&oHV|FcuOz%ghc+|PKxunAgxaIIj$rlEPt>L-5%7ihw@SoPh35$fwUA8yj~g=%D< zR1V@HQ07vXz0PM{C{7@)+4{h{hb+d=2+3)34vco8Vs1IBB}F2;E#BdN1!7gf0?=&)Wpkp(C+y#v7zVyGWlR z!`>Wdbzabxr3OH=b94b3@18+pyT@r4x~d@-%^{n+y9Mg`JbTR%xe2vd*N>jwp93{q zrRSz33!v8I>d(HN=b>I^8MIKn41t|XA^i&nXrikA(@by`S|~*6?8r}nHh#azGVbR? zN5N?RNBx=5WxuzD6?PKtRjOOPc~&2~KV0pfU+RONS)@r(<|yd7pM6wN$^#y7ze2&q zA_U!@x_pvRN=ShfJnl!@Vw6n-!1wp-Zj44LxNEbiA=6u%hiMv^7uLw|OBCTKq2F z3U*V6CZFu+4@F&t(7w*HtbQ8mpXh0daJmS!v%eQxQ1Uj?YS)BC2%uqo6L zI2P!;a1H98BqBZ2YalqKt++1A3XOBBb7vC;w#P48 zbZEDfd6wG>ca1D@)$Q?tt|8qvQJ<^f?i~zM{!$6hm1)hTyf+>0I_(rLZJGidPCM;; z^_T|QIR6-_(k$M#yXsP-U>ur!45=h}9)a-6kpp#oIZ(g$Dn6w}4Yik*p7LC%h00*s zK-sqys-Dl;S+m{%)n>F0wE!GyI66HRf6D}G$*vX)n23% z)d1~`eNF|(SwZVw{SUvBsi9@w!CiEjx1gDSL%f$!9yE^b{d`$r2m*3dk)_l%sMkq* z&C$&cY9IB;2<-m~WwF-HP^qg>PF6eUOV1FL3lV0A_GwVIutqP?a@+seEn4ieyb8sx zOfN(!v_g)dWFyMiWpu_)+x4ykiQtf*y2+?!MesY9V^=O+ASm=-sSCpjD9X3RBOsRr zwS?F27`YlkIBZehZ6gj%hKI>mMRlR+wd9h^Q=g$pne2o5SbJzJdE~^iCS3@gs~SB~ zbc6;HU8ErEQ&49z|5=9kG}IC@4g5}R1;tNAbB}G?fmq8?pZ%*2P*gOwC3__m@|9cc z6t6XdOmC&@J7D>c44Md@N>h(&U2#Ho)mQSK%Vk3S@-9COgHl54j+eFxH#?!B z<|q8}y!-aNrs3d>y)&3zNLsx|Zv^G+Zg-toeG9eJiassYUx0dZKMW=%jzNPP@q61j zo1lSod2g@8wqLtV&XlYKpw1P$^B=yHLM`vSUMu-WP{aCd#(9>jkpCXG)*sOgc|Nm&@GXj4@HPwL@ zfqN$*J4tZ=t1eAQ@6mR${Am$bs`7E)RjQ3%uH~_`z?5MeEe4yOWr|f`~*HHP0DPP#dCs0+%Ms)0XI8-$= zT<7Ai9Y?m|AYkb0N&wVB8MU`diDDS4m`b(Q@ZB^~f2XyT8#S3Qy? zcq3Cl+Qbz?;_HC?p=)!5qCm@RWnDU<&WMX$N6iV|N z^*YinLOF(LiISs_p!}Y*Vwa;oL;0xBqxSQwP_9>xKO}oUSel5Fx*N+1D(5t3PKQ*W zRTZnU*vI<`7L7Dxk4M}DfBeP!T+g@X;h{-KleEGJxnWLWF)?pKr3rmZUcX1Evs+K~ z%SRFF!&}9&Z=VzD%w9$J#(p!sU##- zzqv`TF$1|a-oM$wPzyztxSv_UcTmEr{KuY3rx=zc38gt?GMy+( zL2rJ3^$cSe8hU%%yE^^}K^|GF3ffW$HtLr1&CT5eKTWAdMfrO|>~~6l2`)%T*G{|( z7xW|)7+69dRlOpVHn*Gu##spE#n(GGWjzR`^D8OTG>DLUQ2J{Bi$p@Af#%m;0~Lbr z;+wvTZ^Z=FcMhvNed1tn&xIAr)B7Rgu?Nhn`ywC@_Dt&N`gbT)5k2{+C>M&6tKVPb z_Jra(D@sS~r=fU>_i+{#0rbvtLFKDLHk$iUr)!CSBPiFDY(AWoBv=>icw{Wr2wo*| z)M{Enh?Fk&k&yxkiSR=#NFbDuKJWVc^#(H`8^tO--;+$pZrUa0?TjL1s<0mopFR?j zUqa8E5nLif1KqCG{n$@%o2!_JdixR79Qu7C=bBJ~MU;4nN(>~w8aFSUdkC@#zU5jU zT7x`7JCYZJmZ5;>z_e_X6BHKxlxs;SLXijOYz_+SMST=^JH2?Hpyi?<(zT3Xf+|P4 z;0c>Q!OBj2ceD;CxKKk=$ai@{AmX^;K2v%^^pFtk*i1Jeo*=)MS9gk#Xzjmts4$q2 zq-(Hgk~>O>Z$ubroW$+wM@LY~VNPGxQcwr~%emA{ zBd}cC-ydV_2I(Fut4RD5h3t~nuRblEq5L*CgF&#KDEAb-(RD|d$l8hiXP zhUCSH)+AcV0NUoO`M`wbbQ{5Rmtvug>ng#a<-B^*Ws%^CJ0nu2@QM&PXab!PAcWnY zinnIGBt)s@M(*(lv62jI|k|VW;2>uc+uSb^o+xsP-qnyu@6U>9ne0e+P z39>WqMAWa1q8_Pzi9UoUm|?bgSQ0Y;Dbg~zALuqg#*eS;y{MfaJ5_7(__`zHEK+__ zSOy{Y<>1n%NE}VCA3D6Jz5uOrHHFno7!j1%9f#n^qXd0lPy)CO2-cQsre73z2~H-O z4CQ_F1o!u+^|43n2;QNg0G>D{f^TNVq%7?CKFh@ai~sDOsr!TpsDsy$YHH zo1cdMV!kgyXLy#)$l)8>h~BxA`KT=FITZ5hP@ydt*2@Tg8|@4>gBsr2gsMOq&ck>7 zhgY}PnFkBqWZ58_{`$24r%lK~H|nPNsTs{N);=bG=Y!VOPWb(heL_%hlpX9SXd!5= zp16>Q4-rf!zg3)zttQxn8|JnK1_%z3SA)2QC&3xlnf=w%lHihH7t-X>B{;vd)4F6g z5gd{DxmY@Vf=z9vM)dA$f>E?zDS_n(L3xMz@Q1tIX!+TJTF*%DY4fX?1I`G5I#;{+S&unGA<0KxJm_D=HsOoCxS?~DNdB0)WC7*2g{jUX+CCXCYNp!s*P z`42@KQ1_myH1$GN&|yG(w3*Kj%yu@I3eOdTO`%vLRyJEmF`P5!eJ=;nj6HeS{>lI{ z@H!=4G2VssG1xq5v?gJ{h2kkwU4LQYDFJ31ejl(IaY_-j;j7p@?U%1-%ObJG-_7}o zx)57=u0`=cg%MjfA=RCUEW|eN>G4{zVSu#1@nz@acR*&g6zd3f70`fE! z&QI>XfP9{}ao6=*fPyEEk@5OVK)yN~!_6`Q$nWGZa zGqKqXt|h?*FKi*N#HZ@gT5)A{2HQOH{LImLJwRf< zbb3^d4Uk&JnwZAZ0NKSfsYsIGT&QASSF2-id=k7n4Jb=x|lTC7RSz`-Z-4c$P3$y=pI}-ScR?A$Bs1p zYQ$Pj%ErA{yMT3V+`qXeCmQS1C62pm#$!VaWU~pW0oVvTtpH~eFs$Rw zjeI$&uUM;mOL`?&JXU9wMqcD#jpy_Vn;P&u!JCy*v`(}Y<2}X~nofKf!UtlDAGxYd z;{!6yQ?g}y@gLmrSzdJ)@%IPIf6LD%VPCzD>xU~MY-oSc%O%%NY-U%`X)U!WY>^9Z z^HDa(R;20!)n*~Knt1r$iPz)UYGk;<=q`0^WtY>V0Wod-;c;f(;yzF8me;eRv&$sx z+3$-6@AKNRlm~Acl!s+;BG{(1ZRs`Mk@)i(xgR$^nmS#Wnbe3c9be0*BT=EnA zqe~H$z*9KE)PQJxt4U@4wusImdH>>3H$=0$;hNl;NJKtCVq=g~!dJej4~QzL<0JW% z$^zDf_;Xgl3ikCLthe@rK~48}Z0foE@MUd(Y@W#=WY+r{wlE*jAy+klE%FK7pN_kM zHwoyMli`~9=c!58h5MZN*9PqggUfsIu5YiMo3Hrb!}pV|Hw)_Vh5I$#j8*iAtnJj< z6i-&fAY^F$m?{pj_Xs)k788hP-Rfxu8-WCHvdiNiX_1hqFD$%6i3CbC^xwzwBA!;# zqDZwTVn6sT#q^FhqFw9Mt6zD8uj|{enWP5bgPa%nrw$8Yl@kFlb|?lLOl1i(xv@iqm+1qWhP#UGw;ObAQyC8!zysJJtmH zegaXvH!!F$&Oxjsq3fMDG!U;xWtBr+G!j$d_j$>~kEG?)A}h~#Be@NgyH}4(BZXTX zy}kN@NIv?iVEI8eB*S%Cym7cRjI_1LAULxm__IJEzH#a-dsQYPg?OqO2VaoeqBzyu%8;>f)NWVcs z`V6)=$PXd5lecU&zJEjHduFhyNd|m!fApN~Q+xb<8&l@aQO9k+u$PY7Rums`diiMl ziyS`MB44~eL=K-Cv$UQU`HU~#7l;&~VL{}d+b=~)3L(~?cgcO8?nT1t0pCl7A0xT! zqcol82&A#Lc;@l>Wn@SrytYFn85tibPqHl^LZ;M@(>I=}AQM~e<7K(w2twZ|4p4z~ zsV;s`SzJY`J+RY4!5v9WY&=vK)McUIO`rGD0NO7Y@!!p@RPr?`7rtRhzQV{zUPb1M`BP6Rmwa?w~IMOSapnFfzhb(-Jr)5PV(9RWs z^FD6u$SGG$FTis*a-7Ypu%jM8c9sc0U(7clv$>%%WAjF2V5(LicM2iZ%)LPezvLr{ z1C`g#+;c`;J=zS+pR*Aym)Gs4o80(H&)FQtAa{IBvix=smk|E@=Lh@XCv12ZGv|6P zg&$TWIegSgwE~;2vMQ8HoCjnZEy4#5i2??;eM!Uq1%R&4q2xK+GCt0C{gID#50Wcf zU2x|~LdFiyMz4q%A$!~9tz+B%W7pBxOS2;1&@QDz{1Mgnkv&fSB2N+_bLn4Kp-2zX z-{%!ov?zg8D1`}sT0SK9o)X=}S`eo~Q;jvnAfjQuRVSvDim(2pVQ1Fgk56kimnJ(B z_<;YD7MreVyvkecWusRNHnKC-w0U9~P>kMZzU!(8xO~kpt+(?)B2u;OG$$`mUKBll zKAZ$Jb2eV)|NI0r)wQ0seYgpPC)l*ShSHI0T$b+G2LWVS!6zDiJO)`NDLgG?4?w1! z0s?dQZXpA`i-Vo2ZAkUS(?kP;6G8YCr^wT^}!}2jem@Jg1 z^zu3&EN_w4zz;T!rH*W>F{BlP+qX>b2{IV;xgcbFE~Q>3}Aix zsMztJYrycsi_^ABG(aw$x0;C+WG-cw_j@sbdCs)Q@q*d%Y& z=euN{fL8x?)FBOLAQ(#W*rHelO7EW1XH!c9! zSomf7DkB{nlr6iO(f%n`tT}D1^4&cLK29J6uW${5R zrHDi5C*IkeF819i9edw(rFE)28(W|%boq2+88BJjdST}r14IuohEz5Q0@eIBO^0$` z0Q=els*3o5ZG4yX)^#heTkwH7d94}n7;iY_>O}>--ma<)G86z0LH;LZ)KtI~s4%Ql zcmf-5+PN{eV*s2u;!P1D2~<4u%vb|sfzTx3^L3FI&?#m~ESltElN3|m72k{Em40p` z0gs#TflCItIaLn$2%T$)c*{O~TqMBw+M*~vt?F)}l&~M4yPI25m|u^t%Gf1xZF@qB zANuah?N1PEx_~2vdN<-PDkx?lt4E@;ReGXn?uhJkqeV>yjQ6#)yp~Xjv@*O_rze4XSCXUbF>OL;DzZ{X?vpw8r`x-GZ ziEqg*t0Mj#Bgu4huaPWFJfNv{7HJC#k_P13kf~tai0`#jWOqqqxjW)MavB@S9ZCyC zPCMgG!A2Xhd`SIM?x7}9_{J2(Wvl{36%@FC@Z4qlXlw0I>eXIBQ^gAeQuxKtyljE~CF+R2|)#Ft(@ zW_bQ?6H!e!&S!Uq_%K8Qc;X~!B zk?bB?_<(x&^C?h*FG|)iYY(;{>OdBzGl$s`*EJ6x_Kv4W;&*p&abq@8yHC4;Dx(lq z$R!oe*@G+{cSTcv_CXG5XNKQr(TsL;?1oI)F)rYw_-NN)J|sOJ>EP#8;qWK6K&pai-jJ%L-yNE z&!8r}eO-H~ERoX3vhR{oZb+9mC2zjh1(`5?7SSp7Lsn_7H?&@5BfEErCXMAqXea!Y!r-Hb zY=r~TZiQMP^RypEPmu{SjIEM-(yNG6Pg;_G7G)sOSy3s~uk?tO=efx1Bn^C1JLF2W z*E~M$Hv6NZVi+G52A4PPq~gC7ZzuIrTj0Z|t#cOdN@At?_cPp46#ld9Qr@9jetdcG zbqi2EfEXNH@6^eiM0|P8Gi~mENa{BSNALPxq^3hMFIja&y8V+k*8(GvkyqL@Wu_f6 z-RssXGqHrsKr``b7T;jU`gnfwhOC>INOdr;H%D9oBaSb!r*NLPVhOm}16-dR5 zNvhQ^18GV&x!knXL^?sEB7HynknSmA^{O{LNW1<<{McAEQlCmqrTmbM6m_Qg>(|~R z2`a%mm5C{cCrb0>_}~}B!1K%QM`@n=rFk1~dO z%w-twmDWDh@?zW1IVp9FN?gW$g++tc37^LnRkZ-tEV?PUua&*Q`a% zPrl0W<(VPg#?#MlS7{=#5|bKQW^N=+WaQUY#US|u?~-9|22w~-zP0gX9my?}PK`1T zAZdE4%E^dzBs$6;bmO1_;z_$?$JcCtn0X@0BH> z@2wRFyytEF)ZwrF_$SyYziVz0FYPUw@bO*8i{kU$Z+Ecb710z=Wj)fd6@@Y(MXi(A zO8ARgP50WcWwsv7Ig1>dKXCWzwsN`4+Z5M+k|?o@3=!k5nmM6U+hTyN*%}giq1tP zKJUbO7%n9&+ssO3KwqVFAYkS{C zEF`e+Hz4DU&8vU+8jw7vzkjU;*rsc4b(vcqwr1=sCrN9DEv>JOCLOqqEwIJ-u-1lQ zv*Wr{XN`QZnRlBT@6wyF89yn@W0beCnHyId$gIP$nZ_|or>IZZOu?#NRO|LUdp%e2 z{z4HpWho-hzi6A2*u}H$W*ZSe79S z0OdYBXUOI=ptMjv-Lhp2D5kQ{?D4z+C=O-pXzktJw<}_>n+@C#$cB0>v*fn**bIef za{AcT{DNP~+61=oGU1fbj2yNmXLEAvzALuu`+I3T*#TQPl}Q$D9*@liny_9-Ovff0 zuPL0>r^SXB0>n~w?!fvZ#V51$-#Ts~l1%#>8?{EOK74yiK1**I*)S%oOevcwim zRLYa>1F<>T&p+l~c4E`Q2Is_k9%7>bsnzQHo@0a5l+qJsd-3^`&Jak);j{ZJEG%B? z;FHtFCyfkr@!uEx5<<@?;XkPPzaGw-#w)HHt|YyDf_)_QHFMP# z=9^*o(!=*Vf{$Op7pz^sx?4TLXGR1kx zjut+@i0i?|#N3}A*r>)9Sk?9_hv;Byy#;o;V(Wl}cazO?qyqBd_;>gANdbzIsOBJN z>+N)`%i2L2kbk)CWgE;0$T+30dKXk+8@ESt%1>=z%low-1c&otv)#`z$8lJ=AErOMrvaZ9dRJw-Ll^HW z1vaU!TzH}Hpel7j3f9A2{(XW&6`Ps(^mwPOHnvtfsj&bU0r}w`s^P>|KrOnr)Kpm< z(3a(A-u?L%(34oDEF!)Ey5Gg@3m&3?_HuKS`LEx&BVT6+>uuqLES)4E_=BJJ;D zRp_v#;prq?;1V{y@YFiPaSg z$lP(Z8K3t~`neI$j1O0jo_0K&j90lh^_6{0!TN%^%QSL!06UuCH~?)VEOoa@PK#^U@hH!+2<=OV7(LAa-psbuvp8g&pAH>jFbuY z@?;JIn)WyBEf=@X2@LYBlW(5JHl&|hS~{M<=J!TS)Ed@c!*gcJxS9#pHd3KxKzf5H zM$YeMRa-~oEo&e7_XZ+TjrXxr`Frs-4to7!cm|)D&3`ap%Zj%<&7K=Y%UI7vtPKM0kMnIWxf4Sy>4PZP*Kj#s06|lu>_d2k60gf+zCBo7^fV2PCY{1Anz*#s< zA$wO0aCpxA*6#=aY)Xset0$=eW7qQ_^@R{XJ#aPW08j*^t?;~4f`izS;Dy2uew5gl z{H=HC>?Ex1IDT1V-3rl~RCmRO8X}s^F3qiMe?-;(B=TG?H6ni`Ez(z_iLXU-{QgEp zgO4q|kQE(Y#_E@J^xRgmu_gK|x!$~-fLi*iQK``=VAaChKFePPoY&c(CKhP|p2s)a zl_b0XpWw2~y37T@+jH?#f;R`?5ih|n6%+$bbDEgPwHkmG-|*5K-U4*%_9MLo4uJe! z*#|#fNo@IW^-bC+6KwQ+V$OQYJFMyD;ZJOy1&HZwwbM`zCt@0wEA^osLW~2P=SN!G z5IyFw`!%kIsA4&&LtP)@>!pfe^h6fk*B5)iM3luN@11l;-u@ z^U5~rEPW%yt(J9vxXlgmtdVEN8ml5c&nhA8WGUij+H?K#U z-4Ms4XM|aR4`L_?>0`GU##b3iIDb>;VVyL$)Y*eP0fmDIo2ZNtU`@CDR^{OjxWeAh zg%{)jo~xA!S?q3ryIA6$P|0P$adzhUqb_p5RIgO*kQN6hzP0Ts5z50Blq1y70a>i= zQ9FHa>@wbwnUbLN%MhRRS$5;6zmNC}hNwUAy+cC3Vwf4i7Ln+A*T}GxE+jTO*pk@A zi^M0w27;aqAqng4Z?G>|Bk@mtYe!hEkk}EPctdm=2}h^x<2vq+c*iRC^0p;?3)CGA}Fe)faIRyH|AW=Bjw9q$GWE4kUGae_zfj@q^08;%-oiabl<0M zzaAVB#&9<5u+Il%?DqA2)6MO5X2bY;=lLDT^1SsKyd@l2VzHXH7Wk2=bI?_e05N26 z$0;nplLe{s3~Q^MNk=lxCvNzEc0gj3J8c-hN+JQ;+BY$AR}hcIWT$u03&c4wmBYZ~ zfjFe>Y8_^z5L-G`hR;=NB*|!0{&9QUWy7X>8q@C~W%gYkcZaDW^`r^FKqMjUu6Z9- zV{&B3ar|CUV@dgLl~U@9!a9PJiA z&FQ8Xk9I|SEsEvJBAZu78&r+IAQO=k#>e+hAf2wZh?byWr1~s`?EDThBuCDy`g5ol zN%TxmzlfbgqJi<(vg>G&Ku1NJlReM&`O`kePlGl{%5$GrCQv}~HzIxBUXw)1Y|-Vz>+*F+-4PPZ421KuF%Ch7OQN%crPVuYIMLIo1Ol>O9rI1x!zIj&Z- z4=(K$fIb9Xl`RfyH&pUOqp~d z2S+kT>MJJ5lG8oAhgS|6Jyy0V2)0GK>FSwSOC?fI)AZ~Sn?Q zF(Yxa&4o5*UL-}g<5%3=JtQ9|<&tB26RB?9d-Lj1Inw5Hyezn6jtp0iUMblUL?(0M zMYtJp!|ri=s>27>Vjry9U106VeFKOL|RV^r`fmn z+f_4RV$`V;Bsb=sR=J^wByWeDEcV?%5*ldkdOs7ArMq}fXmtfCKOhr+`n(Tm{Cq}9 z{!ItzE0B%_ehfp#hpofcGDeZ5^vyhl9c9SgcqQw6=rD3>mk2Sps#}L(oC$fKd3>`Xi1KFLj+B{dRjI8>9vYll-gG`@ktE|4? zfq;%%u+70Rq{IElW#W7^Qmba{a5+ql6pFaQ&lzv;FU>!DTM@~OL_^KdmBh6t>(s&Prw zTd9cOrr~&Tyf_k(CC9X1ZTsQ(0i(`#<4A_z@8>~E2c#6iJN8Te1JdaKeovqy9AQlN zTNfO(kR`KUto({3a+=>GmAPJs4v^Q|hnsrbm^SigrzW!Y24_ew9)50}{`%^8R+tg3x#Te^@#f6?1hn%jP$+3r|RG!5Fd zMj5FaJB4<}S)RB3VuvjKrF>slK1U{Iy;w_M*#v z5n(w*{W+Yb#YY-Zh3F0HKW{@cL{Qmyh$Ld)@1Z}HU5CVK4z=I)nWP3fRe3tUVWfGuNxGCkQZvp7`wvdtgMgl4Z z`mErx0_^*a9W^dmqKN(~#lrxdTqGzH!*`pL6^Xa*xND%#gCy&@gV)ckA*nZbO>A5P zk}mY+FkG@lGO=05E!BdN%$53AR*!{{^y94OD#U#xMJ*7X#P|$}Z<*{i3RXm-KjcGN z#s-j3Hf0}k)Ej)TJC}7#of3(QPQO`4#*GTKlUB*1S}RW;!-+70e5x2 zeek|XAasoTo%!RlKy)HwPCrZth^ZUAOj&FI!bcdbx5mN&cVDH9A@LS4@ZRa_Rb6Bxwy8?hhO( z>URK~_uu|XQ|ShR`p*?z_Z9-N-iePjc>_Rl^H?ywe>9MOPy4>4a0p00^sWCXCjlfA zLE(gVA`pEoy(y@w3wYdiXgS%n0@{p>@P^w#Sm%hA+$X?iJ3>YO!4oMoGbx?vgO3eZqYcyt7DF|;kYIs_q3GJoYDQYtp^u&C&a)AoAlV`xL1 z6E~pDm^?9KnhaQcPF0(6Gyraz=1WsS?||UK2S=gcdq8ZobC7TS0g&w4Go3%;38X9S zrWQCAfsEX>pEEIBK&q|bysKO-5bx~@Egx0}LY{+%OER|rM^Wm)d))v)$;N*~m_8qC zyH#@kblW?8k&0C^$bSM+QCYQuiax|}I(S?wISVl_FJ~;?XhN(Kitjoqs1Yl5jUtmp zFgB^i7|v<8h^@88za70~0I2%Dk%6e!faPSANo2@6;8Mo*O2$Y)fI^w&?AIY6x@UJ` ze%}fZuilj(IWq(#4~#tNUbzV*8$)BNpEpPB0DLjxfAeO%aQc^q&F z9(%KKTNBV5Wds*=Z|`^hmf_y}Q5>)6x?69=;*2kxGvhpDP=UxEG3(~9jw5Q0mbZ;- zv)j)ssB z_*r~W$;;p6aWlTM@#E{U$HedQ)6zg?en_5 ziLD}qS%7Y-M%p;jV0#^@zWpJlDd5y_O|zsE13b6F#~8cr0e<&`A_X=pKv23%BlE2n z5X=_gj2XazKvr11jA|0#lS^9k6Yc<9XE)xS*|`_6k&)~V)Tsj6T?195TMO96lg#8x zi(J@nOK*W|KpXx(cTyp~V-X*FGORJ+}Cl_bf^EZHmZyT89af8G1|-1p~x-shaxd7bC`2+@DHsu_Fr7YR`5Yg?;6i~{lN?@twPc%0*0(PCdue`EYfYVSq{87CcaHS;^7C0{gS6>2OkP0Jk(fEp75{jw3DN#Z_QXuF9joxeOG2o3cL>>~VWP942l4dE*w1=u>k{{cv?%_b1=+ z_TeWBRYD$dr{Gu8{w>bcH{rJw=(N+jM)7;3E!12=w{X8oH9Eo#QgEAW(gGg(KDhmt zY~&hJtw5IMF#_#*GCJs@Qfmmm9Ggxh*kr>aUp zhMSwg_20Q0for49Z0i5Pf&Vgh(j{AX62G9MS1A=_fM2cF=|bvX@SC2)6Cd($89%WotQBbB=!Yfe1E(!3uH$wU%6N^36u#*0sKUN28|y7vy114 zfX+o|t|4a!=(EYoM_9FizL9Rw+k*t?URsTHx}E}>dpOg5duO1$N_O#3`x}tH86u6U z@WK6MaH+$kcj4BY%}6_ze&MD&#J@XXzPLAzQ~|!nZ{P>iGF$HtTj6JPgMYqYB>Gqm z>U+xzS@EkKR^v$*5N<~7d1)`l4DQ$aN{07}S-5rK6P#WjL~%PD`MKr1-*Nw_|1~sR zSq73}O2hX!fe_2Wo@b7i!7gU2jay?z7h38lKAwBaYqv%?fS0% z>AgLWEqe3?o$vzkR&J|)d^Is1GdsJ4aRLQn#fgwbJEG1IVbrFkP1K`xH+{a%;0~)k z{dt{r6t{OQskoI>6}O>P;H#O0!!3zxRg+o=;J!sgX8k4^#g!y0-ANKRz<(x9xbgXd zDSk%lc;YAZPNLrRnCg^F6n^Eq+5Cw%4*c4X5rtK11b+R1d`-yZGkznnPU!YXI*?N2 zy0Muk0ohE~$qdc_pm?d2*NxEv)zQ1t8EcC`{cMqI(~S6h6O-Go%zXr!;*@g_{q}%* zTFO`A>=g06w6&q84+8QBl6NC$l!3$`Af%(?2X5y??b9;%W!#_GPM*$)5!|;uZ5Q$D zFYpZuPYBVaUb#9qaP+tSX$$E z-k!8NSFMKMnQ>fSY#9M^b;E)`*TaFbUee-(j5N?_^VYDK`vG0Sti$0wV(DJ$+HYOI z0Zeq=_1}`#fa#jR31ORmz__+wDeT@043wvPEO`NFTGw88`QHKxYxZT!7Yew4ld>&& zer~um^Lyk|&H1>|iF}hAY^C^5!H3H)FGt`PDOYP^+~x5bN>0HM%q#d^ub|sv#`*aD zhf(^Ib>8>`@4ISW5?S#F?~k!I=yT%_tlh6UNB#uL+m^iiZ;AIgRrOa$`Vbi0Cay=H zZUJVD-SX@b1F(k0=83me0Nbs#i*DQpz^)$Zc!@p@*lB{k1@&A4*5wE+erE`n>o!V% znn(jJ;iMLccm$AL#Q3eYuKZtLlho=|#|;-%uh5j7#!tNFn9x4`7r*{n_u-Sb1Ni+y zdnJ>p5TY)Z-kw8iiIMe^%TA08V&s1APakM%VdM(-(x3SHF|w1z%eb9>pb?r3y6@Bg z43TX|Oah01<-fGALN}{`ZTs9vW10|fnoK{RS&aj3CFbBljzi$_XDd>ZNC%!Yasrvq z58$@HbnjEMAaLAQ8`y}?1(rhki{)7oKrPwi{`!^yQ5WswacQu{b$x&9bd{5U|07a* z#lGVr{=h-5y~Cgkqj1d+mK5Z{s3Yw&NK=I{n)btm^ne76=3V6{?*U4Tdh*oS*+Sc&&GiwtEA+ft>y_y%2Dfh`bxS?gf10Ltad+jv&Zn_@mSH z7YG}uQKYq-f$;9mNJ-`g5Ng{T=o%XYeuo##r)x5RlOR7tP9h78=0Q2TF%r1HnMcTD z@6qF1|9zp7;-16*)h|tTy4a6V?RPJ3{xQH94suT%0~4a}RjR7=nij_5%A9-mIt|AB z1j9e&yMi%RXFm$j zgj*~aqMYu77}5Exp8X3X#A})b_VYoKGo?s??+Qp1t`{55`hr-%b_m7lFJd0Y9`7lm z2W(Xn&ZQQ5xPKinzeDOC;O8UvUBBqcW7M3PZrsjqG1lWWW+5{-F^;CwkG`Kji*b&R zS3Z6`ig6~EecOCDjB&`u;5{dqFg7=(FI^Kyh`;0N`acn3yn13b6WjR?*el*QeZF-X zcsBOu4<3C5VKb*EK`FfI+G2Hy7GN+PPGI0Sll&J)ccG!V2=6c(G4b zbs#HyNaUx?g80phACyT&!2K)mv3?2#?jPc^u)r1L4=I%X&7LU2*mV{LS!gCOK3~r# zm&dMRf_eM2pGmYa!5*=kM3p{FK-5*jYwj_|8y9BTT@{3Jofs9_a=wbOM{=gn`~ScH z5nbxrMAY$HKbQImE(8COxCW-_Nf7N%rs7)ZCg!)pir%@bprF0L#-$SmDg>J1UUOMc zt6h~n6ifwm*{^wvL!+P?zoExIJpppxM}qSXO+e@_ZtE-2-?&95adPl%7{=n}S-Ndr zgbCjd;7_$3BkE*f7dvy;F^PQ&|HJHPOhQ;1B zF^WOQ^g&bW6zO3wCFpd;+>yOD3wpk;MTw2!pp$mvzc;iYpmtA*En_7N#9a=cEN~1L z&!}-_w!0k@ZkXg#_{)aLHm>-cS(wKZpANZdQZQhO!CPlei(bOyUa2vUl#*c592VL& zrB^U~2k_C4ADetQ5)4>1x)3orp6p&zY zH{x7x1BH;wrg`})pnm5m@7bs<(A_#seVp|S7;;4aEjjZRjQ(2vmQk|@BfWinQ1%Di z$qyc(BXJ;|uA6e_^$Cnl;)G%fZ!D&0)HzeQc@@)0CY|V#KEyPi>!?h8sKqq6uBqQ~ zX2(M`q zl@0(Np6m2Kc1l3(l8(=5=Os{JP2Xb~^a726t3h<=2jE`#hS_?Jg2@P}M^sWhSiHHu zddVsgEQWG6nz%oM(Hp-M!GfD0XLyE*@t7AT;d}iJb($TfOIXPI(auJ!@8929K8VEN zzeVNs;u=iPk!R!j+$N@JznbT$^chom$sl#_@qJ8o+~9TZ^d(Gu1^+EtjS3SQD&IVk zXNOzwUY1%6CkHa`+LFu4RKTKkeRPy@3;5sMY3@Il2GWL{dLJ=iP& z4+(EDr+w)t7=9M)$U8`?A}@l2Rgy+tHWOHQv}86>5kGgkCv0Y`5Ytvw3d3KH#LQHR z%daYRV&=We9^uCxV5U4J-{h!`Fhd@}nnx3Qn9jXRVJsolnA*8(6RxA5iN5%LOCp;F zlc;^X$M8G@6SmG2mg{_o|3*C#@_@OVwoEYnw zZ(f1|(im6Ltrqp-bxhD`q2HjN6B8Hlc3SY=!DKchqU5f;#Z*W`zFip_!qkh3N&H+H zFwLvOS|aY?Cs^IJ|q6y#LyQQtpT&wNDGCW$YQzEB6$$@N*aTR=SQ^M!n2b z*s;S*+fHpL#&%%nkK+4Tp=eC|{GizT3PVhlr0zSJ_%J3@;bkt`c@h(&aM(E4V@~W} zy-%P_^~E??b-EnXuVTc`h1Paa41uEn&?T^#xH&&SY%H6l2KxRC|A10uv!M>T=Y7h>3EZ-(O7BqV*WeFRjCs4B(t3s*qxHELmpiI}3He@qG1{qtNu2F{j;;R06xhRR^F(o& z`nA~Fr4O$${d<*7K2$E4$?XgCULlV$oAbUkc9M>m^XU{u`J}g)yX^E!viD4w%aW6S zQ@A9ixRKzTK?ne6@|*Qz{7c{_`%T>2j1BAreuwYbF@T{_5a*>Lze_S^)$8NJr zJb?MVPXUW0RC>0N0jaM_}n9pU^$BlznPuD={J7(ZM z)zHH!mIjQ<>%v3ALqH}>7qilP9uxU>$tQl;9FrV5(vzgef+<)@HK;Z=U}{EYf}a#F zVA`z6xK*KD+UE2_o@-rp>t!xGgTNJO! z(;~p6HkG_~r=2(_z{Dr=MHH-?b0&U>tAKe3bX4N`g6>$C0@fslpJ#PQZ~;BatDw24mvdRN2D+gJQB;SkfK(=t1lI6g0{WxNKH zX`TGxtT2Zuy?FeEqRIeMD{F7jc_E2uKHGge4i&B31{i$__0Ydn062yn zSvCqq(31Q9V@I0+Dn^x6y(CT`?J&h7VRHZ?85D^FR}_KU=YhZe@mOF<&+1WFvcLq` zBWp-r`C+0miVDT#l$b=1OS0x36DISXnIQqOVhT?Je;%aX!<19UC7R>dFy-Dyo!b{# zG5M<{!E{G{Vxpw!kAM0!;LZ(?%QyWjBI69#&m?dT4665T zKR)gX`aJz14S9Z`?QzpF%fczIU2Bw#s8QTlB$8%v;-ewB4aKd(#ODA2_KL zG|quq=V!w9Zz<5IQmC_!GzSe_MZ|h13#jQ?R+f5^fzlkhrTe!7Wcic1`hNWc@c~ui zg(XrDYG#V&Nn8f*yi;qwey54`#KxK2^(Ks?ZB=?NFAd{feG^n}!jJJw<~=%}lZpv$ zsQmkQF%1*07n>^)*T94$$Aqkel^YJzMY%6nE2T>A1gwk_;roy;cAP;pPjx8df0FRhvCf z^AIG~llgA`(*se}jOGO1kHG(MpT_JeBXF)&i=T7D1FJQ^U*d~yj5UaGZiP7C#eS*6 z=X*3O#u-H}uvLB?;~J=Fde@YNajTzGb(a5!ac*+5)1*ycY_SUp37dx)%^ALMW8;VT zg{d=LHx&GFhk2J+Z@6m#dtYv(&bt;6>e8_I-Eo6hZ@4<1ua^VKLyqlPwS17e5nU8$ zkOGnqB7XNrr-69weX1thJcw%kv=<)V1VK%E4`x&jJjNTW)bF-{{oBxm-PhNFg-`Yc zE~5}*_&mky^UDEaPB1R@-}J&*y2sb2?#W=RwgOF}`g<753cv2%9VU#adaggMGZ3S3 zcy;P!ZX|xcCm{QLh6(=7PN%flf)ehZ_=@Cv8Ba>%zKjHvY=7y( zZTVY0Gm38o@+lr-;pbpo3Fzu6!KU`4h z{eG?xBhm4uUbiLsyZ4?BDETDe50pY#{-u)O_Zs4jev-`LH;rQK|I+2*m%={weJdo+ zb!0w}=DOdF8*`R@^Q_n!w{B_|6sg9EJG`7^U+KR9l##jL2PwmWPG0)eUnyE(tWuB4 zO)LiHHbwC$4_ja{==pN~k_9mD%V^)U&`!|G7bZ+rEjF0(hSgXbDOul(*ruPt}BKuOhBt{A0)1>M)W7D z)|x&i5Pg^%78|v8K%w2_;q)g9$RNZn=0OU6^8(iCT&IHH*gk85*8=?d8=g87+8z9A z{hOTv2~zxTue9SAIL+}tb;-gXc&6fqnelH(0yS`@lqKI9zrDeY4K@CyM*+BBnMcpM zvW4O{O+T3Lsz1m5eS3C%JMI>c)GBRN>Kg;u^n*necRW!a{t?y}Vhog;T}<_)H9+Ze zvDiu23nCXFf4Kf087^tc=gkY*q11%A8=BX9kUZ@t$KwUF7_!v4OydC0u9&7N{`8Ck@s51Ln&A&ceh8AowSWJzPrrx+rJEZ*kZ z>>B(+7RrR_ih&Ab?$Y?gEK?Gh8J0H1UXMhkRF2(e2bhp?!~36af&`I~yW6xZ^(TZ@ zCuYK*6d^h0zSp;$$dK~a-ECd>9;B{TEYcLxh_nW!l}DXJklq^}*0Fnz2!3YbTbF{7 zQSK*2+u{~vQuZK`>6H>P%isKON#7G$SU6@xc5WlfN2CoSWA(_YjQ09j`dDNwTflNv zkO5iipZbwMG>fdp*0J#*DP+ahNyx;VLY7q3q0=r*$h>h)ZokeRndz~8iM-E*Odd@t zDyPySW7cD9xgXPzVKQ#xtgkhaQI5&JBsPbXj!El@ypTX@tw$f2Kl^~R-Zhv0nwmrU zx|#kxcbtj5L?1M3`tuhUvR*Q|RVvVo zY}t;8H(g3VwpAbA966JWY-?xTZ{j17t-^i~V=5D}Q3aMwrcz`zZ1Ltb)gH3ccKO%8 z(}v8WX1Pg5`H<<22nHsGJ!ITGuSsG04H>Rfmc9@bLz1N=VF}8)NI^V$TO?^5sm_l? zQhSpiEzcAe=fr5F?+MyZ%Hxn>==kQD8B1hh^;X@_*BhCuF-&}+PDPeGVd+l2S;!`O zCW=I31KCA>&%PMahU`ClOKBX{K@Mjxw83VOuu?jh9@{;end_;~f zxflwqw2(tL$I|UqL1b64ClHd9fowE>4vE*VAJ)rMR*7a=0*KqRp|0Y*XuzLB}|ek2^%> zse&(LDrb-hCcAhcq8AzI{CCW>?FM2kQGckp`4S1;6!`v!xZcu@B)69DSs_(rlE3v1 z`AGZC@>0LFBSN3v>0B z|M)B9YDUZ4emsrHLfXpr&g??23CH?)0=<#TEk)z|Yb(fEOZ3V4mutu|+q3Oj^Dwfz zRCT3*g92Ht1UHnE*dg=0G@}x3eq@pm6fC;#g$xM~+?upBVAJL!BIb06m34!azCi?u zWIUMC@W@B+(XE2t5>~M42b*3r;Kj{#P@@j82Ofo<7{zN9yaeK-;(Hn-0o_v^Eh8|=o zc+Fom)CxJUUG~PF3?S#ZA56<@smSel+$YcTrpP1RF+z_m8F>oYQoAJ&BafsYS(+RU z_lt75u&gk#9;dL(N#*&ebh=ZYRAT|~|k=4-AI3CP_gbzYn#4teHm=P0D< zAuqLE$z&=qxs2wGx-^j69c#9DQ9k6Hc13T{PXO6tPEm`D z>d5lID@t#B9+|RiS!Q|aBSTUz3GIIWJPoobN#q+e%0>l*U7tNJ%vb`QCW@u@Tt7Lkih4fhnj z6*=q?ETRp`h;h*3c+SrjWGYsv*p6pL$l}pC_N4PjN8Z@$W)v}gJ#Y+9mux`_-3jeQ zZv+X-ymL<~Vmb*bJGyssPw^7eRH!av+E)o`<1fCHwbl~Uq9c&5yB0xtMmp9jMt~r} z=jdwUN<2rv-%I6hZZcARxoxBWlp2}5Jvi^|V~FgZ!RLt&zmfAVVM&_Eo5`#w@42 zzjGn0m#LL*k;FK?Eqm&n-alkE*KRwwdw>kPx_`XA-H3Gb73~KviV(TJUq=+Y`jFI7 zQ-#WVDoEftzeD{cSHybavS*%M3nG`*2ocM1gq4S84SchP1eFtoHUru71hwrxzmLj1g4atfP7Y>=V>+wEYw{wgh#$tRn`m5J7q9usdJQm>~3E&gBg8d<3B@ z*8&tr)d=#s!q(r)D*0woP zYg!&UjQNNRhDdeOr@DyyJKOk5a{*GnfD?FVxJ2M1>Chc`kg^-nCDeAQ-Z$8yM4R5pWb-H#c=k2)gO+(?j?42wJKd8J!Qd2$1eK~) z{)-tih;$-f@D*<e7bpEBg#2Q!mCV4XP3w|y%Kff9!Cn&GD%P%lEo29 z+#?Z=vO{Je*Ngb-A0nHRkxjMJ?#PjGQKjxkA#zwdpX5`dgDn2I)CQmJAaaqvGd=ht z0i72e8q9By5EOMPIEsol2wE1b^VH5A1YI&cNV{J`&>0_kV-&bY(6*0HzTLG<&{Wj@ z8c;4mP+M^|t^Tivpz>)ny7Kir5;GwoF^`c$vQ!TlI&`~`Du*$h$+`v7*0Mjw=E04S zM>ua=yaq8pSC$+5N8}-UQBK);yhl#v5o7#K<;eB%8C{>H737X*oV+1cirf{v$6O@* zkaIxODI2LsWVuVGx?VqpbSE8VT{PwqKlg_iKb{H#_khX8;sPFm+=UZVC(=$5R2quv zRnj;Js$GF$&wrmJsHU|RPr{2W5~P-N~N72-WeMs<@j$&^z%zd;}kh}_mB$G z%ja8DdUA=FAB)KUH^7C=jem%lkiS5-#~HO^sEP6ZEpY-^n+&;KVz;14>qDM+Y4OTa zv&cL7a@?WbYveT=qdpr9#VF2;Esk5~okI zYq5-CpYJj5wFV{=_@jn8$y-_p5flUyP4E&I2m7P<%` zj$YbMFWwWx)W6kED%6JN6)$jTr?vyak!b-~D#!r0#X%oMV`yc68Pfg!t(h6<0q z5jo~#cQb`I7m@3>SJMTudgL)V8_rNShP*}6p4<;FLPw-$j<%StBJZmq89UzT#QUY5 zEMM3|u4R`kOTRlKM_1YNVFB^Tl3sR^)U23TUyU+XTB;z8<2~Ze6K|0e6F!<}I1RB1 zCYC)I%)y4oT@OuTNC@->{Rb=7D+Eq`3KWq35`p+eiI>q&B5$4gj+@^FWcHJUCHd<- zvetV1L_H=4IdDJvN4711oE1%?Wdfp*YyS#U!z>f>P>nQ>TirvRzpJQbUiu)fhZ4#3 z;>E~|`^pQk;ZfwFm*&R!ycW6I(`hI-G9brx6)xUl7G$j;cQCAFj7;1#bN$7L_122x zUd>22Qcd);DU)DB5{Wt*(cxW){r14s&OfoRY7^yQbn_&(^M^+R_f`ZM`qK_6UDiS- z;rTCz1*nkuQp;l**BWF!)jZcE^A|a2zdnY~eM?+l?7a5#N63}#m7r^XA#z))l>C~K zLcFire^vUtk;n6$iMv6h$V1)oz0}bPN^cCdF4`Q`qvL^>9b^#{J0 z3k{vf^!q_wQHBjd!g=k(wth%Umfzu8#c`xK8lN4Yo`}S5D=sQ(7a$Jl-*pPRI&cuH zSK3;?f{fJO8Zz7yA!Ga3jv zFLCKS^tVQ?{-5XWrVJoA?*c4lMhCe)mi-uO7>(T22u|sNr;*Ec4X=`hDRSm)zkYG^ z9CAo$>5bE3L)Irv?`KlY5cfmvjGOZ3kYRJ9xSRbl(l!3OLu2=Vxf(1EZZ*AUr;H!g-?r=RlNjHg1b(%h77hWL?-oT?&Tsp{Fi7sICs|2z; zV%vRvzYjU;=iM3lz<``896MMbSmg4giM+D<5^^!JmXLmt zL9CAiBX=d)kYlr0e8!bvWcO4v*Xbg$?n^eD7jdmeW=8X&$91fbk=CGIRE7@H=T|JH z*{?<#cs$wWadV{9sYVzg$B+~$lO4IQ0TSHEthh>k6&XDfZG9g62AR^19bXDjL+1HM zw}O3lkd;mx*IBh7WHWW$clC-qvX5Mxse9Fm9NiRi3~N6kC-s95$J1Sq)2s38+>$!T zY2fas=Hu(gDdyaY+58A{EU_JJ5bH6t2ugmx~~?UYP>hNL5A#s>9-n(=|w&p}jijRxwh$5OFe&P=FK}(}f!^ zT|`nIm-lavxg#TK;~lf4=g9cq;Ggo|Dr7qIMsh&U9GQQ-@Vc0RxE>W*6=APWB5S9g zrF|cBhNb%l6woMrVkB+P(o1Dry>N#Tl z@MGhtbCe~rc;C*u)-j39PK@P5;?5!CrAf|c?>S@`OjlP%wD{`J8_(+*3LtH-=_@tk z)<|8)WKM;cA1I5EkH*FwBAIu7!SRX)$gpbpZO73jWNh~Nz`=6`namCTz52rdnWg^y zv~|r2SK`bOik=7cq}33T%Y3wQKgt~$@(Oo0q_ZQVM5DS( zcMXy8^nBwVGvfItFm>km-LJ^(Qdz{_W->A_?x8Bad=XhBE9>ljqCl1c>i!RNNRVY< zuJhAXIb>R6)G#PFZ_kW?p3wgug6k1I%4V**ruV*4o zj1c?Em;U=$QXES3{XdeL-}MBMA8ke>9CRSQYk7G$jMyIysGs+m4gh(oqm>dr&x2As zPX3N_AgF4!acUhFf!h4L*!uKY(C|rZnCi>`&G>J>Htkf1eO50<&Cz?H6~f}U7f1$L z@)G|z*!@5=@(1~$bDP+j|`}*ysR}5xeKb#g~QJf>mC(RlEC6maiGY;b|op; z5@gpiJr?h~fHYZSt*d4)@Z3v=OLZY2q;N7vqR$7!e2ij`%1(fkxY|#J5*Lv3*X@Y= z_!*S`y`Fk%-UF(0;ap0msX)WFgpJLo6SOk_#*h_Sg3j3`JBmvKpt~xWu2MDxdLwy} ztG9bWUr30X$H5Ww-dGy!W+a1d-=%>j-*wOlUJlxN*$i4`aa|q{8$shyq01MO2~aB# z(36qu1C?y6C1(6fP?Y~q_s3=g$YyFJt!C^1?Xr2Vb`!eMS?i! z$S;n56<`LH|Q-GIN7Idu-yWkwp8T{Z)HK#-+|$>ZYii)E3g;Lbb>NhlL7b1tH2u? z+^h4f3PgINz49ktffQSskNxXopzzev=g6cMsQJ_-a>goyHdn^b{Ffpy2=tz+zZV6D zho`P`DWrf2@hGz`W&to)t4RAZlMa^b3Jo?g?O-L(xw-PF1gw1XOGbY`0?Qh!y1AdD zU;%eDGvqSB%;eaM<)9KU-tG_&NumeCYlr5yMTNj1ez)-NQKH{bT1SX4YX>y!o-VzT z5(1U&bN}t>=Kv4&^J}XyydZM^t&67e07ymB_8#6o1xnIyr8X*3Kz*x=!M`&Hbmw>Z zHtA=;FytWj+mt()%4{oNh!p^fRR(+CiGHy0abavL1F$?@y!D|e5BZI<&lOH3N9nqJ&tPvQde3<6t*;c+l=@VRlvjOgz%aUdI{JOq7) z+pgCfCqT>YZ|U&+R#3Ck!!cQq0jJkycIc%V5UP4kKczzflD*HQE&TF8sc1npuBiz$ z`5k^_*>{1#x!87Tra~~TS^W9qSqNA>KWF%|{U6x=I@y2g2tPRfaG?yv6TrFk@`&@j zEpX|(vzlpM4=xMqA8)+!0T<4$|cpvbLLGZqF^~ z6`29XM;V46IjIqSg!yc;qB@`>I>hJuD-hI=(?qXo^8hPnm;0^Zm%yhrqWyRA7l?nj zMN74i1B%CvH8_-(f~JR=lftd{0D32w9NKHa6sxQCT^$6Qkk{S0I)30}O;)~EJ^-$b zYLhnUQQ#3~o_tb#13Z^*GLiRdf#>0oH#eGRz@zDw$7GE+xP@s4oW;9?bLKIfC=UWS zFlij>##@4wAJeTZmPIi26B+edqX0vt(9tzo1JM01;PPC}2hfm+Ics*>8mNTBn`7^q z0=x8}+-=qs5dNW>aYH5#WUFEw&fqdZgP!rBW=$A?O`npUs5zL^&-kX&+kt)6m%@CJ zi{PrIn5uut9=vFH$scaafY0)bPydP%9OFOr;R|ml91}mjK_Wx|Uxj47;rto!>OHse z=O7Z?x{@9yI9U=onx!W!X%)d%zHeZ8vk}ZmjK>_RQoxAH!k}C6B_s54(~F?A>cukKS(uu|~gc$K7R6IY*QJT80AjV|_z=ZB@an z;rGWbkv?!#=aM>e#SGlb1AEf>ap3E0p~e>(2mUeJTf33<5J++R{k1USl_c;Tef=}= z&+5GV`&1SP~qb089v;xnIYc6*75&YnJPedi>8x&6V(i|a&xc>HuQ z<@-*sk_rD9Q7QwTRG};{8koUKWq8{GP5%6$spz#b9r0T4d+3b@896YFeyl# z3HqZHO@{I774`kvi^aq)-k+E(Ho#;lPdzWV!GNi}c%RFmaBIhxhneF?lg*qu!`rn7Y%M?rUB< zn7&q{Q}MHoN6=94DqQL-9Tr_B>~h79p|G z3)u%#o!8un8>V2KEzjHC(gphUXWCv%kAeD=B7b4=Adn~FdpzzWL!8?XXcRc-2)yZ^ z!@MnYfX%Io`VCeJ^e)W>)6^_L$&@4aE`0?^8pAG!UW&&!dbnClT}v?`O+LA{m1Imx z_;yhLBnPHq=bry_j0e-9F26UKNSrTE5g2|-NWd&>yDnZ|ufQC&!+RbMUBaAGQ(Gpt z+lbuisq&a8OH3zob3WNJ0u!>=AHSMoi5Y)e{JTd{jG1?(39vawVJ6oR%cIY>n8Dc1 zg2y)!F!hF{-=V2InCwW#bvxNwOr$8Ppjy2V<7D5iFq25e=n^k0;V=@6B*L4;E71=l zE#J_8TEK(R`Q0&>QJKKl>+q&|4J(+SuYsR>>^zabv6<6$bOKZEel7pnzZBEDaM~~` z_7{fxNb3~i6@(e7lHTk3z=|1AsJ}>xQ^f!4&X`>+Q^7R0dxpMg`Ctw_M~c30-^ZMr z-p@TT_rn}JTw?yR+hSI|Ey3#3`TgXZ_OUzk+JI||lx zTTBQpMT=3AU>rdU$nI4E?%%JS2PMHLfnwADJ+@B*On0?o3w#WLTmQ?OKWTj+rZ}zJ zlR*LslUAJ$c^aVANN=lU&IU#_fBI&Ux547v(#@9pi(tE}a#q;xE7)BSC3mJu2K{I= z~UE4XOhz(3V zbm`6Oeg#aPSYkbV>w}52d%KPZ(_nmER_9DJZvb!etYuE}K8X2ucb)a=0J$G>r6U6o zpg|oh@N8@q!1&X%Id1}(FH|RRyWzlzwUnP>d>)AXzwy1OdLmCb(WYm&90Jb2F=nj! z0sixZKPg}SfcGrrs&=^qF)Rkly}fk;q%@d3MGf;X#UAx)`pq3oC%YB3P!X>U7i?_2 z`Y?@uJsH+2&6rA}`HEDv4<_@6W4-Hf5GIy<`P5NQQ;Z+EY*<$8gJjiw$+x&EP^uS- zusq)jS^_HXC4E)FP;BW$uHigbynx@{EYsk0l#lZSAsswd@M+mH> zc_PH?4B>+P&!%mrAZ(BO`qXG91pV19?fS0{d{*z<@iy;)BZKOn+e(w5M-g>FmW1f* zy*ssT(3FhvHImL&mKtIbQswktF&#`o^GjjoCIco~{-L@%!W0wuK0|x#8*z?kgqO#m zs}vNuHB_v*g+aq2<-rHp0Wc`8NgDL622)RxD!++pu-%FJG9}LiuCe9WS0u;52S;^I zq`LqDo?U|19wiVeasBhaNE1YeVleY^5+Y0%eD1DsLnxoD92wIv_&-qeoe7~K#=-A= z-2LytdO4IuT%QE=jy2D`*mnl;d3h1S84e;>4=0WVRy_Ds|6vw*5mvu~Fc8))C`?8(;s35GbY@kax!U>=>sm3AQo?0-5wmAGLA zZe$-EBiHYq9*S52~-f;rV}2_Hw(e>ouse$ zXu#hmpG%;b6}+^swtr;k1;?8*DRmF}!DRZsYY{TUeD!7K#jTP9kmMkpRE#kMmP_*` zIZtIqGdZPXou9>AH|W88A=kt*2bj0t1z}wQpaw zK;^MoXj!H$h@Gu@)LT>vY_G-cxUE?NiRkN@@UAZ0#9Vpa(-aXN;dun{Ea1YF|dJ96O1lTewJXtl10;e8+Z&JnwMBaL#uc=oncx$++J{eO3-*&T^ zy+4+4Y_j>c_rpFo_8VRe%g2Imvm0}}2n%>UzmsTYIRmaQLaFm_7J!5O?~i@<7?C&M zel_!CJV5LgdkNDCP@nmsdNgYfWT~#-&+P98Av^c-_8|&j&Gy7`ktq?+0jjz>Hx&lz zKEL(9w6A~;A<_GI-65cByL#DOonT7t6e&w-3s#&l#z)^xfc>PlsknC$IA8I=J-bv5 zZi!{_vwfW4@w4lM?}J$IGS0qh7E}XXy7aO~#tgut;jxZiJn=r)MEGxTRe)2B2Krth z4mR1{d^iS6Fnyh7?Df$RaND1cUr6!+jpgB}z?33T2=8WG5cmmVz1+zzDi*-wBzV0r zrxBPM-m&Cv?GWeHAkAVa611&HM`MrW0FH9JBkk`97_IzT&r89;Z0>fB@|F)+o8;0d z#6Ji7*Fis0g)_h@ivID1m|k$9owuS@>;#v&Z>%!6j)TkKnoEgjDLA*kYx?>q0UW>k z(dotKf$f*5x(t~Ou=ww^g2MG@G_*87V`5bH>%#IVUIIdC{0mQUk$6dVlnpH7`p z2D|1UjjoH1V8fyLmO<|pSd@>PU1(qhlVQo9W8OSq$h_o#FZ%#=wZr~n$nybpi+eNS zg>#@VLbxIK`720j$X)Bo)c_&SU&gS8>`Ekd8C3ZooA-sMHi^I2%C~9r z#W*K2ztM92?)p1qasxj((=U(o-7hU%+Wr91?!;9!EMG(Ft{*brdSxKR&ow*mO?etC@(rs&*2xs;`VM{^>=`719lU z9!^9S1-rIqkTekXqR8-vO+PZJer37!b2~CfKmAL4ycf|^R(d$&-XQH~Pe0qeUqh<) zUhYo(!bo92(XCr^3zF?G%gXS6kHp?AU5g_4xN7a;@qQ6uNSn4Wi|!Nt0;b7CcN2T0 zdsAg6_N$G!-(vSKTGkRiM>p{U{<+B1`azdM0SO_K9(sy0H^MS|>ni1%vk>{lTEub2CM0uzhnR^V8xlM0p`D<39I0s2SGuFdka|t$ z_j_17(!The&Q}wRXiXP)9SvGQ^iPX{i?Q)Y54U*b65NRl;@19LZT)}@hX>CJ+XWNz zJL!+df}4@CM&{=JGZDmGrJ3JKy$%_%Qa;HY?mz}t?{x?!P9VKUwo?b`T8Ms3Q-YT< z8_}fkp`P2_NHg)&iN6A}NM-%c@Wv%Wf;+=6V@2=>zlTf_{j9{;D&sb! zBerLEz^%uKx_?vdLk;0?u&F#+5Vaer7BH~2uhk=pkO1{v@G>ITyEoO&#v!?+^NTd; zpGeYBd>2nhG7`Br=Tnoo3|pS5-gx_&1^B+!w9h`;0zwg;tzV)QK+N7Q`?{kVNF6oF zIh5`La{XI7p6{^+vhw;+F@@kPc{4)aF)0MyxY{{&g(f;x}d^hRVYsO9+NmW2>p!&mp-ye4uURFW*VeASr;Www*vT;Zai zv^}=;`1BGe+)&*<&ae&SKYo@vx2X;!4!b2j@1ElA=4v8;(19`4U#)+dXpitIIbG^X}6sMjHsq91(!;tVNg<%@^ChO2r5=l?X(&$ zP}eV#-fve7nx7xpSjPMS?WsBs=hmNu4|7XM#I9bTZ3&B+-5&=u=4BbB@kF2&)-d^A z69ethEn|6GdO?e^lTUs)4b+r zB87mZv}OtGSpn|75mxH220^eNow#)64~V;cPuI?^1KIQ-Yh-T=Gt{cKs=g zxSt8q$i$sxa4(SeoPDMJ?>s0=#MuWH5xhl~1%{i4ih#Or={Py`0qERjE2k6vt8UcC zyN$K8gin=&{n{sc&=1JH(S0Th^nIg+jy#tnd~}P?_s1}SZo-MQdY^NkGwz&vt#k%x ze2m|N6qiBEcpb$7mAApdRw z)&$jz?a{WtKmIT9@{=VHKhrN``nmu}ub(KNd%FTkt?Z7X4h*2KdVgdpDh{Z~@Orgh z5a?8>O06Cz`mi_3*Tts=2tO~uZJ!AP12K!`o!8yKK>G+E=UFY#PsdB<78F6x0hb%j zW&`@gYWfc8=Rix!I+e&13|g?j!7W7Ou%(`e_a#$+lEd(pNFoCTUYg6e>@@ zR4J=#cc1@wGU;{QpX6YzhgS@NO>JX9U$d0Rm+NO-F+)VZ_m1CcxEIJz!$uYo`QHCo~^H>n0DIm&w6TS>9IV5w^@9faLi`m@Cy!StRi$RF?U_Rr37|3rm zk*eR|Bz)R}U-NxzK+(h6|K^e@kiGl@{w+6voSyo2t#UDWl8K#*U;tSS8NyA@;PJv^E;kRwR>Vjdz^?gh*PM&o0kb zAti^r&A+ZMi9YdeMAnfcWbBg0vh|2RviNBe=urC@S%v8Nm7n=a%rSI@dnqhPS1yz? z2+t9zUf1-5VZT{t-h`E!9Y4XnV{q*x{Iw#-+HYTu{s$t-!TX*rE`#8Fqkz<0Ch)1J z3_JcD2aa3f1sU?zNYV|dSe@xcvQn{(P4D&~(%*jrPlA0B`Q^@V?q`G#Q=)SkX8e(o z_75y3xB@Bb8tI>;yhoZpKE_%!CL==$NwZDY1fn0K-HoJ^k;`b+vmMQ6kb7Cz!Px9% zj`)2}TzhO20vVW?sC~*DjrAN2x=8o#- z9V%t0Ta@#EquQxbH>d3tv7t0nH|65=T)Hw5Nt?!fo9S#J{JucuOGw0({odEuFX(2!ndUZ{~OU9!&)?8O-uI{wykGRg1EuY=Tko=(A!q(lRXiI9#W~N*; zuEEM>5N>5@!lh;YetP&qIp%dKU zV|i!zy2rutZ|iXubE57k@HDe#KTOP(>P7Dk+(vRW^Ez_{mc;v@_$$Hw1yWvP`u+MZ zj+BNyx|}b!Bl3`rvf0m4B=;+5B2Fr!i z9&IO?^}S=8e0Mj{tku(Ss36IM9@oXMMQ`jt+=TT$LJ`T90+&ocgLXn0$< zVtkn0NF)IXe2h$}VIlh3hCRjFwKQa@ZP=o)`yj#35P2!tbPu@)y*sSJ`vAGPJ}F2! zPMjNeV?6f#tU|gQiZ4vF<&kob#jWAi8`L znb+0-5*&Pq(8mI{BGiube8k*lTg#AgyJBgz&K#0xTll7URSXICKD&pu9!Em^dWGok zQ<2DLgD0$~+K^cIyT){5f)iWFna}aW6p6mFGX4DZ9uiSvnyB@B2omJ+Q;>Zcl97o zR~{H<0n(ZcPH0v zUhVJ$5#v4Dfb-a+Rh+wZ5B{ej+n zH@{kf*w3Ahd@$E5545(Zoc4Kb&`Nk(T6er1REsnZ+>`JEO0>Uptp&j+x)I5h<+=&v z649V!&jpZEmza2~!VaSUZGU%Gkqx*KymS_X(U=Pcv)dQVWpY^V)3bgdc5} z?GR!V-}2RW{WZal9Ov?hhy}8fO5K5Qg734(NqUf3Fp##kPNallU}F*{)5@x)_fm; zq4%I<(fKQNqZ~iIbDw6LeKU;;yEK_-Hb217eABHH6uwID&}38rx zpW8s~pToJ@4*>n-D%Eg}U>xOl)6TjcjFEHj z@>De#X-Djesi_5nz=NmU+_r)4VvUq}Tqcp55c3V~7y`{N$7{V(B0xD_{cYSdv9C_d z+O_ApDTtp8;{D#U0en{l-+F5kd%+23vXj`kx|tH6=WIhYb;b0rd&>+HLGivqvF4l( zsAPH13->G${oc#S{LkJ%9S|&S^4$k?SGhF%E)CFo=gunm}Pq{>WCb-_ab^7KB zz82kaky26!#?}UAvF3ST6sKK$<_;MQ#+B$#N{)lx-QOaH=PN*G>h8W-c0oR1 zrS=+Bp9h29be{_?>0tQyYUJa|7hw45_vs4{k8~iXFwL=l~Gv;2Fy)6d)@&ySnb_RXgD#Mix3(ylA z^StuT7wE@U{(7wW5jhI=;}&h7L2K8e|4zC;1l6a9Y%UgBfs%CNkLi-XL?1CS!?EW( zNa}nz&0re^LN7~p^Fu!dvfMR+9i4ST=eMG6GdFHREkTCTQ#&!hS48gGx9{gg8YA6i?0t zm2UVF=X=?R%8ne6?6oLKN%R8YdH%X^CmG=DlH~sBJOk{W59m0Y|iyE7bg zB)ruvyNJI-`BK->#$uqgI);26H3J&IS#bQtP|%+9y*QZP0h;s?@n3glKuzTN)ABv- zKxvYNveso#IHF?o;wk{i#wB)-urf&AvvIV)Ndb{v-+U>X3Bcc!xMGmT0$b{Yd&(0X zLB6H*^fOHa3M+c=^Vk@GqJNRlQym1ABdM%Ed)|XuUekqMF$vJP7S@q(#|c^nLqFof zh+OD>-P`JH-GIvVt*cyqfw+IId{VEuK>HHiS5xdBXnD7-hpm4FjdbI@6&qzxb1hu? z97*s5anV*6&;TU{O`nVpoIutaZQ_0T1Z3B0w`;I{2Z`YUAw1U>M2x?fM1ODs0f}H1 zH9}g?#n#xm_2EY#m6VUwa$N!P*_%zfXUae^taUy02Ep60oiIJCejHS|nq2QvPlM{6 zgIfdIxL;7ZdnF_B&cIhRR^<2Cv~8m<1rW^4)a8_Z(8F-<@oO83Jr zG=ah)3Camo;{NBZZ)KnS4GQ)=7&)+$m`};gS?326xth7NBvm@d)E^}o?fDPHtCdT_ zgnohG%x});S-G&~#NXZ8=ltsa4t}{|Jt$OHyYjaCOA^&=QU14 zLGifEA$ws~Pxu&+J)vgP;pqxx13`ECHunK=fAIl{LbgUikp`~;_6?!dv?V8y?gzSiq6ZrWpmXp z?_;FQeje5}>+m#0>;igOQhSafLCsGpvRV6(80Jc`H@Jo*o8~8`n+B2e?^1Ui-OEJo zbk3Gb-%^k)6Yno~#=oFMUfe6J=nKm43gQp#>H)Qo{koZtia?_~G1xTu3ut{CJ#yEI z$jygLeH^v-o^t&^*1V#!>17svT$lS(^%g;-H@ILVX;_Ebo@K zuOs+U7S%2HAI8_sU*=E@_DnYW*Ez0YI7aZc`YMA&PlO}!2U)*Pa5W-X$F+4^mKsD_ z)VAq<(SzhiE$Mqb(hylLrC>bs7LxC}VE%Em6qGs6Rettl0JQ)xSV9y+|IsjcS1=x^ z8G&I3zsZ0OkKnpyO9ap-_mOLr720G$T+%Agnf%Ze)3A>l2 zpz7X=4pf$d{N0)-HlHXUSQWy(&~UkK)XM1UUwwin9slqA-lYpj+%tBSv@VKBl_C7w z0tn6(cxRjTSsrc|yb}yQJ~rG)eFz5i13`iep`ai3 z`MR*zS)k`{XgPVw5uC!Pilf(vy5t+PmhFK=5PJP9+a~X&S#SB{g*$jE5>~zcw;{0+ z$?tc(_bnMilE(ZuHyisg0eVOGX!h{Mf?@RR0ot=JFyR}EeP)~wri1#) zc54LJ&+kZPOWJ>6N*viKAtW%&v=`d(P!8zYyB^k%2<}_W5A@pA5QGLCm@j`^Mhusj z-u=zLk4Oyc*qhmxNTutu*Y@zeNbAgeU|`}UL@n=&%-)zqv@J2>mn@nQ^%iTl4?hiQ z?`%HBVnpz~Hd6I62%VKW^Q-f`wjw}1XVd(C-VAivBkb)`j)B2!ac0i}>g#f%)6OQk#5`=fgIzxO4GhvcxyQBqtsfycGq5h4UYI!%3iNR>HN? zc@YG;7~{`w`-+54Md+S2>Od+T1@wj521L6h?2F2!41)MZc8+lZEh<_NUKlsD#oR~mQ?pKl3)nG8@ z_q*f1x&hc_pYvbThrlxH_7|^dFRxA-%+`S3~E#kj~QF?D5k^pcDT!~crJ(L-7H3#-zQzO21jez6rLaCh|eBeY+ z3OGEg01k&+&IPmUgB82Y!+FvshrJ3=*q2)`mfI- z3p3$r0qGyel1caG>8W4H!u|9w)elDzR{G#{aY_#|l{nDraAyn|9w|SOuIP{S_?>UR z>3Rzk+6WUakxJ~r_JkYxdX1@CtUgw6H4U=JQ4 z6E2fG;wjIN{+kr8r2~vW8oy%x@rM?uC}o_)Ud{ncHK*jv?I|!cE;@NX?Jby3yqOv( z-~~Gw4bQG&V{kbgzH9e%3V1B%Q&Yp(!8`cCyDOY?;FYSysjE`}&SMQ4m?;KzcO$1u*FMe0L?HX|C+XW6ijfU# z^_HT6{?hFm4zl!X)8}uQfN}wE)Q9VB zM182x>>hCtOb*qBOwx-Wc%vt$lHw zBNaSUT^qi?&IcV=$LJm5SCOS}79Xup1G%4kz5i~L4f5ilxbfP*M4kc{6Nll@PZc&<$aY_0a7O=gWGTVG`{+52%ygvgnQO`;gQ8=(bc1BGYnO}04{x;9 zt&lWT{sw;n{@`OK=L@}n)Up@*aJC+_dWU|t4JU(fPGo(9z9-nUDKwOPj{vv+BQZDL zc!6)tmTHAt*I;{Zcigah4|qHOOS&k15+sKrmw!A_LvAUxOZH_M$d9hE-!JC}^4lTR zHu0($d3ywIPH{>g_oNp=x~0j;apgoK>l8Dxm6ktqh`SD1oX~Vj(jP)5^=AuOa-5L9 z(Q|PV?gb=LRJNKDcNvMEzh>!YVS_|^#I63ZP?5l-$EAF;Mzf8@a1C22qW^2O{&z-b zH>f|qkn`-V4VX(d{%1NW4K6QDj1+z+fnWZ=t*Z5#;1?#OcFsl=>GswCo5|%u{$5l5 zr%2b(_LjwkJNK24U*TSgbO|%^(iIC`6*`057#ewBJlKsKm}6;6JNF=~WY;p`U>SsN z?=O4NyBisrG~Kg29*pR_WIje*-H()h#Il(!$RjmlsXv+~wn#ftz41RCS)^0OZ z9i-1ItQxW8gN#(RE)6AzBI7}>zVbpqhJCwxYhUINdx4n7Yotz~{F&ZoEItC>EOF)q zTBXP*%;ThShbHp9h@xUXdLS?7C8M2_^T^$?m6=>_iJY#gIW?S(Ms{12YHkyJLdyX2 zw$yv0$V`{R{}7>nGH`uy?bU%v;yNOvWne$j+EX3zxik(@`JXa`3={nCb%AfSzD>wT z?n>B2Lhq0DAC!~W{U5S+36}h56pb8)an1_Dui@s|jgE>`A}^U?=jkR7OkwAMca8**9O3+0rwSYPRCTQ)=x<$ErQt z;uWFSGt6spDf^0aSU*xO#$Q2HChpmfi3Ui$e1652s2`Lc32H_}vJzbV->*#r{E)7; zZ_l0*XJoV|@z53CAIMxu#1*tJBkLO7I>}evXlq^TOc6^uayNC!zLRW%yr;UdG_P5s zZG&P7zbsFpo#9Xat~2dN+bbEX2I>iYaJxawP--)5myQ3MQTPV}!hSc<$wJ`KHo?G_ zKL$3JD;lql;-J@E)lU;8bo258Km9(PtJ}Ej=yqSC&TKuscjKD#9O6NL-_+Mt6Z*ap zk=bEOB&fRLHM`dt=`X%euF$6;)AxKD0g~rWJ zA@U4cW9q!UV4MAUSEGlK;K{l*^Y&^aILP~d*%)RAEJEE6dwRqFNF$m=IOerxL>8;uC^mX>x~0KM9wSXk2)y%Oa@t&zqoh8dMC0o zlN&fp_ye6DFbA@|S|^^nt!+0yb0N=#toWagaOAy3N1@ARC-Mm1%+6SfK~^4P=L!`* z68mGdy8Iagp60vKh*UP%mf zUO-mwCE|CGL098*q0zNN=t7+t%R^czeMMAEd(RJyohRBt1SUo^w9ofzQ zo=8#PBlxZR?i!N6f`K64@qV=$aP<<*i#hZPyfa~L;n)QDT-b9-IJO+Tc+Rx2OL~Kw z`sde2%2L7p>b)@>_`#A{*-d$Q7K{cb;#ro(Kxb#0VC#@Es7`i$&rNOzvbMs*;~Mt~ zJ=oC_%f@TSBu9V}cy13e_YSd0*(gF*!ik61_nRP_H*ZzKYvqvL2Qe~mvLJhUgwgR% zHe^$%nBeU94w)rP@ij}ZA&I2WJ4N;epr4($>%22F*#FWR8lq}~+k#i%>41Fj*fs05 zjd>2-o&MXAS|kO|mDnvKpIDJ%pFVbuOmhikb1sR6n3%hm*AmfLW z_2k!7WGceds-`4@OlCUGws{F7BljzkKY|w!JxB;|JW+@gLZTnOZ|{UHx6B+$KNGo* z*$12k3hM#=0-60rb6~wkh_cb0N7PHr$IWoU7eLw7!28ApYzi&8sdt&dN@mvax&J1Z zo2bi0+;jrtPbGGObM=HiWt+2OrW(+K*n-#74 zxE9?%>-QYmwy9t`8&&0Q~uo%~J zV&gpy=03wtaiRubS{`^b`0O_@Jo%J%C(;{qW3BXjvfY6?azT7flGsaRowvUHeh<>L zZzS(#@kQFz^X5Xx|vK&pl}coKWI zB9i6Mi(fa#iQFfH#9I>5AZ(!F^&*Y~R5h$ddwhinor|=a{GaoL4}^j3@Rt=Zz0g&= zX&?(`KSYu)c}IZRfipkId}?5l73TD|<0KeP@_CMKOoQG4KldAYAhCat%s(4k589u@ zE3O>b0JY+=FHb&sBUQe-Lh1sM|8VicB#X)uq#mu2ln#{_W}j=VMNG5^4rv?pM!sY`8#)8k$WXS_AFV#V;rW z#60!)vBPJ^317^+{F2_8TrgD1F427@42D`34}Qlug8}8um7Q1gK`#V-1wTnMTpo+TRpF{EUjl=)&w6;b%qGP^t!5XGia>PgC8r1U^R zFFyJaQkZkr=&~zCq+y*f-YE(azu6lpG%1ESnyG1Gwc~Z8rxqA|EBAunoZoR`*bDMQ zO&9)SqJioI>mj*0jF<=hOg+aR4LZIbPiyx$5Ph_5#z(Rj(0>+71uRtnJ&eon=D0EF zIJ_o*b;t&4^L5XRC>PLD>z7EVx({lP-Rk?rML~J1wbeH7J|sUEn8Q<9h{$NTgF`D5 zky&M=cZBL9dC|$_zS00Bdnqxuc^4BBACWr}vwH*yFsshrSurqMlfDsNt+KCf{`E^P zvxmbV5ED6djWZQwo&3fG6QV&$D#7J8p^H?LW4Z|+c7kSal|k6EQK0hG`l{`?2-J9v zl`og?0(I-NV7($n(AKla_s_KlO~nUuOct7;_E6L!JJJGF7`(GtVu`#`&$05L6NGPG zWuUzFf&h}0=vbe*bQsCV6z4{laU#hj%5?VtI}(f7Tj?yNfdt(jvsa5}A$AKC&T&53 z>}MGt)!cAz-SXk?S&_IU;Cm7<5?Lz-693Udenq?h`Be7R={?Dy6wf?qUXunY$%{MP zM4o|K(ydcTuha;=-3_tAkUyXvQxKEEp-()YGh~i4+MsG)Cti5F2b7O#Ds|$7-|8^f z9i6@i3OkM+Z?SxXBpiN@r&zv5;w}2;_YwYCv9mlPZYJ7D#4hJX>~(_cys#iXY}JZ5 zr|5QHX%S|dTKJQsKiAAYU)kz?{7yvO-_)4*wnv@CCIgKxx@xuHmr=O2Z{n z5gI0-cq`_#xP=ucNTx2{J17M5aohcZ3|ff%N9khW?lmO5FLaE4K@bTsv{lEwo0^AD@TujBI!C6BESvqwfl{?^MAA zjvsk;Y5X!S)fee_HB}auxn=h}d;BaeCoyCpGWH86NqS$Di8+mv_6YF3Jaq^s`T6W) zyzGz5wSI~iTvNqmxphW$uN31_fdk|>Q>?gHkIb`z)*)OdU43ZitR2o1+DgMDAK}cr zNAeu~dSFf~@U`6BU(__kFntMCp+9ND1@D?(VNAc{ePr|au`M6duIi1v!T2149q;=B zCQx&G%(-F-6Vj@iyjbxW6PBn*R_k}bWy;Lsnbv}FxryWN8)zFiX;#N>`u=`g{)O9M z!kI>#Y%p)Nl1BXWET8W9cn_EVQ+anqHwu^M_Db0@FOSQ~?n^x4E`m$3`lm+^^5P=% z?OT*oHE@2)8mTTb5oi9fTP#>-7Jj{z)H`W-9u2otvt10{f;J*-T1ytcV;qN{TzvGG z3FBuFYP|nM6BC?h-3!o&pT`{zPtak&FkYEpgwaCz|s zu8R)2I627A`i;$fT;bRnlZ=w*Tte3MIIB6YB}7wyn58bNcQFzuR$^tn?b+Ehk_+Kt$*9(?@8IIZwm7 zPZeX?6#XMq&xi3H*yq8RP=yJ#J)~qLrC_3M3+ESJ118=!elJYN9Ftf)#e2?03X|O3 zSiP($fRj#&FT?ypoU9$RgU^xySLocI-JPp~E4kSIqb11Tlw;fRWIk?OS%J}uHoSl< z%Rd>S9Mr@q_QBsuuzp-Iwy1m5oDr9=dh9$Pr-Dl>hQB%jw{W3c_vp;&5}aisGd1Z# z0dzH-}#)L_B*8C1En0Ug1?bARxOe&7sxcvAXOva*M zE%tUdCde|5F-9zKJnqPThUd$d`^}(sdI!{kr`Q?GH6n=GFF<~YlZ9)SliP8#% zQ!e8oMST&G94FAv3e|q%B0rchMUG!}LjJv&LcPMT`OgJ#mk5MY&Dygzc)hvWc!?Me=)GEtVWt+4U9IES$5S%#*>X8uhs~uIF$X zZH}|F1zzNj>?z*S=MpRN!O~s#Sz5eC7(*Idm3RKJ1YeuNsW2A9^wVuY(=e z&K~7bkO{+SAZj)Gj*8Pw{>-s_y@k`e0}l9pXTa&gE&uJ)CpeY2C(xMEV_ap~-9Y*l z6D}QjPpstLEG$li8LMPQVEmfp{c0Z zE2tz)J@o94z17K>+JO^C|J+%`RHof#B&ph%GG)cXzf%sEjn8em$)16epTzpEoo2x) z_1hOOjNZW2<2?K5;uM_fA0%b-KnbTaetF4hYK-f7wccPd7RB|S%jt0)tHSj|^Qg-Y zGjRGV3lH1bIb2hYamKsgDK1NGI8l4V9}YR3a}M2`iAgrB*|C2mVak5?$_bOanC8)Z zkrPus7*)ZtqhjP0M)mDJXy#smX=fxgolbp=Y33QfVJVNo)NN0cmbXb@Dm@?HcKj3~ zc$V7z!g_~s(sBRg<56T>vAKu0+#J}z;>c(X&)cK-58e@!k*};FQ8JOZ4 z$2EFiIY!+gPqr4pF`YjDD-%~2G5U3<@TV6TF*Vh#x-I)J$+|J%cvi4=v!* zXS)p@-2P%>EWJEe6^t?MkhF_4=N2)&Z333q(h_FiQ~l;?rX!}mnAGC*^f{*IA_-L?@3=U$7w(Kx6MOyr-dU%8(Mq$-%el5G=P7A zcG3Vtlh3%>V*)WVBXRMBrYg+nr)I&Lx&)^0YZbM@(16i-By!Ia{A{XWLtj=M9aDeK z=F32u!IWEfn1;P!!U`+FeWU45Y7{jq8>~wz@TAI4Hl`_$V+8I+7l5-Bj zzpx{Y4`}AN(2CBECDsvKk(p;%Xf*(*ma56Nsvp7i@BMUNdcYfI#jLG5&tQ^nP2aln z<1rcDY_@9SOib?jb12}|!AN-n1++OajO1ruAQ|C{NxeIDgS7Gq3e7BTM6R*X-Q&L_wvDX^Sg2%N5o>fQHoa#zFa4A6_|ylBfKzb zx|-~+!+w~?^R>J0a-#^2+}9cwO+8F$@1feDtD3}jq<;zLiN{E^u@|9jv6$TZ%ayF! zM==Faxy^ z?oa^<%rGuy`ZlWork|c?H`42bDZ@zLh2fvL>L>AgZu`e@Y1+UC2Zc`PmiH>Swpoku zOH`g#MCUOH!@$dvo){+m&TvkwAQF>N8S)eky@W|E37k0o=L#k%#ChUD{y$9YuOV{- z=VeS}Z1>^L!JU|*0}%y!X9!dNWb0kZZHj5ARY(h-@x`bgc%S#u%Q0Pv*&`1pRWL*5 z%maOG8JMZ1%Yoz0vKThwA@J|aDa<_hLw1S8cg$#={78{W2@_Oa8w%LA4cEKF@d7LR zfE(BQHSv{f;=0C$0nUxu~6=$5=gUx)kf|>b#l>x~dsM|c|#*%d$Q{5{al3W{w zX`UF<&~Qt~Xo{|mCps!H`jp6Z?Wr|P|MzP1q;UskG^%PI=Xeq`S#ulH3OtLMygu*2 z-6f6b&lkSEEc_f(=KOd4YO4al#~yx=`bh!T_Z+?7c8%aY-(zRIR(*iEmni1UR zji-$IYgci7cB7{viuZBq{zJVN{rPa^SJsRxocD23Q>Mvk@KapOIl)7i&kpCSNW4gT zt$?#k53ry1P$PJhT7riALNU#j7LJCkn;3O-52e`jHbyr$vyp9Q$MlZeOlaLK!3=bM zR}Y!~#SAL>*V7ksF}eAKa;%! zH(8DD3+}%{4CwTCRtmJ@rhGK_W5NUhO74D~k{=uw zeqi)-jleh_uw=unOecG3B+|BJYqtY#f=l2p8VeN3^y7Ravc%O!VOHta=M)diMk{C zjx@tZoT~Y$GHz4?S6$qlq*joMD_#iU8;c9UNpjIwWEGC#5}J;}4^y)+m2I_QXqOWganxR+l$JaYB?Px2}#jW-j)uZyRJi9zbRogFgdMpXk%%1Zz=;p+f&$RZ( zpWwlySh}Px+z*FCiEJZ7MVh$mrH&iROis9F{1#XHBNjOQuje!Qj~H=2Cuy*;+{g8& zouB=9{uS3drubxU=>we38+sx{^*2s+S6LHJ3B}dzW{*A$>A@+7EkgXizsJc=Z%ZPk zlyRB6WWzyiL5$Kn$6B5|g|FiK8Bn#iGbj8w^aPUo#D#*f**`DU~O{;<7m?rM61OW$M;OHRtiRY?)EC-21L z)R#%k%32q2x+e23rDrNQef&eh=qC@HZmxR7Ad-aBsthwn9q!}W;iDby^%=OjrknY+ zDTXUcy%rxDFToWe+hitGzT+gFL&N#arkLXJMsvMQ3q}zMuRc}0fKg^Dw#?=dbsts0 zLQgmiqoiKX>=dcR6!TqMTLyMv@d@;f?|C`waVG;#r#By?o`zqM{-?TUv znsu=8`{;8PXSUk;;+$j;F7}DjnEpTqC!er7c3kQOuH5JQ(@VG!SI_ty^ImHh*GxKb zu9$T%uGtZ=Z8Nb0*XXD+B!bZjt*IIYhM ztEVv`X>wqPK?%m1jr2cilF|ERG4hkBe4OP>^3cOO*|tG6o)!pw6qCec^2ifIDYcl)ow{>nSOg|L%KybOJokCt zpZELqj>&?3IUqf#`aGEiz7u`bQb^{wepz1D(w@vewEy${QY#8P7xs$R!8Zi#~rUZvfenL$MJxV&vTiq z$IfM@GtB*3kEW$M{mi{Uur$E1s%5DQr^NfdxgMS+)10e!^(#3>#{ZX1_nH_@CMM#9 zFZx|33&rxZh@BQBi{yMZrWSfc7SoDQA*MEw#d&Wpf>?U8xb+6{XI4I0Y#G*i4-JyV zBuQ7ghokHAT5fvW^(W=Xe&c%FbX}vqXRqjSR>!{x zNmSM&ls2^6-H6v?&u$Y5>@e12kVzb>{Fe+4a4LMosdvF>c1Bfg`XZda;!oL7twN@L zv-R_-0*TDN+DpucS0i&(-_07E{7UAl(LNjU(3>nomQ!!m4JDKMq~ql8`jJK6Ra5@B zluZ_8Z#$Yy7ep3WROR!u-yw_mEi%2aN+gpC{ZFU=Gb9U}mRw>Xej^Lw!nQIsuMPyc%E3NVTAIKawLaAG(5oDgkYR36jJD9%q z?*KGA^d{E{rBa=v<-IR{clZ8vU3CeLd$wCR2htmSE zb(fJx9gU7B^PApaW)?oHhu3t!x&NbF52u#SpZ+&Pk2SC&(0UFn7k@jj>;!f{;B>RGk%WDbmLXba}h-{ z-ehQNsOue>b3}=C=T;nMTntUs*%f z+P`=lxgUaV@H~E<_`n|Bq**<#f7uw_4A4vX#iWC7YQMOgbHfweu-ZDwixAPZ#|3s; z10Cp!;;^CMBX4v8e9u>wc1QopY*(FQa6td??us7t$D?1b*R(SfrJ`xE(gD;aO-P&U zs#>AaUF3)5W%2`xE69k*UU&VzJ2D>nPD(1b6q&3u3q2}agZ!n^IW|1Afcy<}TY8Pw zp=&fVqfyWo-Oxe;1U)y=t>bM;OlL;X9kFmmwq`MOH<35#rsD;4&(@0R)hBy&&zg_C zaia;{wGR>Pa(Rtz+m&j?kfza%;{}{+aUST(`A=HowEO66>5YIOyIgeiTj~9tlzZq` z^tQZ1A`9}gYVZVCQvr4?S`!oe^8pz-CAWTE_AD})RT=YAl^dBc>k4zP>_=wPmIrHM zl99Qmr@4dh-N+m|GI_Ua0o`cenvneMfo^9h{*2}}M)%I{YmE*tqX(;3_mXJmAjO@2 zr44*Kq*N83D{+~Gl(b<|QyWK%POFOn3NU50KyO&m7=%#QBc-EO~`_6<_R zPtyu>+d?Wk@1g4WU`Rb0FSy?-3~4E!4~sQ?g|s(2c(S|CKw4pD}o?6na@N-4hoyoINM{|*lBUIEE%6gFvo^acyi;Bp)t6|jGYg=@$ ztfkZIJVx}ZJdTJ9&_im4z8B%wKS0`-JmNXHYmk1cruP0#KFIjFTOioe3o?DbtjR1` z3z=$yek3FtLPpgx=0mAONT+Cmodd2z$_}lMV&mu0tv2&ZB+5_dKWQ6B_mV&8`-aTf zc|-~s`F`bUi-Qlc6p5o9%s!25*!FPT)Hsi9Z@288ICCG_{V+7t+)ah-5vadv^U6UOuehSj@+y8uS(gYckoqEDPzK1MxVl-bk zrXYK7<1(4R2XPap!{T;wAzmQjc5$9E#M$L^9QAt%S-uBL_0I%Ay5DI!r91!6Eta&t zH4#SVR=#tUw9*PV^+o4c1Flnn$m1V4Ks3RHrPL-BL^sh?-8tr1_2cl z+xlmEdO%r6w|a-m2T*kD5TC8-A%|k=QKxFdA*GVfEX|z=NHeEB?(p`kc#~b|;g$oruB9YOXJGyB8pj_z*V)7J;{)px|9yrW#F9HT z9T2i^4rNhh`9f-ogJy<{w&+3!{)^O_8Zzj(+Vk+gFl39+Ri2-96;N5bf6y7u0<=4C zb69P+0R7Du1D6?N0R4qqUxH>c0iDV1b;|MgfF|XoNwvHwpvv?qi>qCPRJv_Bk=<^P zt}{s`pko5Em_DL4K4$^(wnK_lU=(up9v$)vtA<3krxvcKuR#9uF&-BfvZ0{-qKgVM z84CPXdc9n$4te{$_?MOpASc6wwxBo;GFekpH zaTo@ePo6jsKhOhg^D=vism}p>l9}7qvqqKke)!RgLa64N&* z68SXFPxA#7vo&3D3dZJqCGV7|&%b~o-? zY66OE&madP9^gb6blbfG0q40+k4NHtfa~q{pQ+h+z{U0;*^y-!aGdx7Nd}VuFVMGC z$Yu*z-!RY=McM$ytCiID*EJxOhlVW;UKTQhZM<_2eg|=uN(X-ieIQqQvt@-|8Wd1} z-8?>~4@GI-5xJ2@C@K00?e52TIS(eI$;uH>jw+<`a$_x&-u>{ertdx!oxV{fTR9Fn z9@MLKvlOBiPZb14&8Pt!`|p#Q)qH?Yz)<8ybu-o{K3?Q`^%oFQ|3j|fqy_wB`8Bp> z6c8;UQ|wUz!1+3OkuOaf;N2dI63Pq#%YyUWVx9jW#g|7}a_zK`-rT|c=H5++lQ%n+ z7n}>ZgQ|uFl&?UcMN7?;$&*lg+Oc3)&Kb(-v6XH;OThR#OhHn6KcPw_Q~d3m`%wA6 zk(rMQoKW_ZLa!XN0_2&%d7O4&Mt1Hn$GKTF0KtmA?f+7q1M%LNB$pa1ASvR|5#E#n z#2>AzjOsrEB7+KFs{OtI!KVLCnompuqQehL5Bdba0iOgdJUb6?6jfA$E8XbswSkPt zH?)v?%O@cM%bBvOmMGA0+=5&_hZhnEIZ)V{)J#t-fs%zwkLxT8pra(V2<^`ACb^w?(2tqt@4JigysMz*>5*1&>YAjx8Q^Y#-kD> zAg?PUfTH10#vE+|kYiD*-kC))j`}0UUk`5rk=W5Gg;)K6KklEi&u1dwD)iA>$5~^X z?a`ADi!0IPd69d^y-uKeDt3+6edi#ZygYT&SAU4VnaWmB8VdOz)ckr%s|O`A;yz{? zpN2|-$nMj>wotQ@1tvPkLQq(2$@hi|>RWJDD<334El#$|8q60=)Rcv-Xb_NU$$0g8 zG9Tl6HOf{z!F-(8lQupK1p>ADi1zW~UqGo`b3)S}^Vw9~P-7411Y$S1UF>t}fZ(0& zxCmAh8yy=CU+qQL8}$m8WI`bg?S5Q>a3RDW z_)x@Ag+Sqk;OrF1Q7EVRoINJ18)}eLB;L#SL4D=7aZ!z8@R;7QzNSKTh=^rzax3E? zm5rMd<)H!4NIRZ&dF&(5=gYU{8=?dt9XF;_o(Z&Ss*@{KG5=|v#Y|fD7LcxkLSuGBLeTe za*5ZG*4#H*j|tJ}IL*s-R;v?`D!PcxP4o!laacLEbM_zBPv6+ia|PqowEug(|8yMc zkL8sB54 z3{aj_i4)hm3S{E5RFS|7 zYUYwV4=8C54GTT-BW6uzJNmpdH@5CNu_-uIjVg27)|DYk)@ zw~%qnc`Klj+hid5TO7zC9zO?$!-3e_-j{wW3qYV^bfxQ_1mHR%Vo*CP3D|q%w`fc2 z0fXm*$BJ|-faZuCr?h%5V9E{Z?GitS`AI^w7*3b~;)na*&SGDH$m!Q&ieqjV&n>}? zk>w&#KK4y9tX=_VI-g3}P2>c+3%A~sFE#;WxnUu#traLm+qTPo!t$E^Z|6AV3n1V2 z>4Bo+BCOZD{rS%V3BYIOVoadj1%f8%)uebVk7X?F_j~3y=JWpW)#|n?!0Xt5o-^A8 z%!<~raX5cK8|fOXw9Nr1Mc-VfPoM)#tPd#FTdo4UhcJDtT0h`<<&{2oYZ;KDY=4@K zWnq0CWBw^v-d^c--;lj5Gtf-BC(2-s$*OPu5pUbw1ji^TFE@%60sZyaf>vP*z^frc zY?pckkv3J)!0$6qf9H&o*5w(fA$q3d)5Iw#X``P|_v``0|Cbr*^O6N!uo>{m^(;rm z)k1#OL4IUqsXHH$X+hQ{RHv;qxsWwmRaGOANo2)SvSs`FSHLEIE1im&$v9JblyYDV_}Ua(&g7 z|9b~CPVF^0)Nw&7#P2~E)dPr7B{e-Em&0RJBDL1NTu|>a=#c0Mglhd^x_Kw@Q2K_P z4l~Bn;>-HfuF*0Hnf5&?LM;EGThoI++O<#6v6OqV`MPOnWxn7lM`I__CN_0lTB`}L zr=A=it1AVZsSm#_j|c#M$ryWFf&~!e3T9kbXaLf@&D16?8bGmN_-{kOEKvQOwRU39 z8)#HLRk)a_1k^R=G@l>-2c+yTKM@|xfjqsBEKa|{a$IFTV(rS5P_M7UYK_7QB3Z)c zLR_bz#trcc`F&@g!X=9!9n&%>W|W+FYIYv-+`WnN;^rX>2ccZe+6q$mAO36o?uGfq zy1O{X^w5cG@*hoy>jAsDLE`4UXu$C=g8PyA4ZwS6vL;k}4v;()1W#mF0`bq>yTqHx zK&J9M=YN=ATz=u9Qxo@R9dloURw^6aJ-ODP6UGLZhX27}Y#)a^ zQf1h8)dFyw0>fe0{*Jo`TIP*a0X{Xi*$O`#5Gap5|KRl=5Ev2vyylw*h=?MLGz>~D*1Y3xE`~AN~-Ev7IxwHM%6JvqzbasWa8jKb8yXTUM>&d>PkB;d@tSHD8V064t$ z#An*D^}9p##28aQa@d4VE$Wy=Cw-$jEX+cd%Xnx^LHAqaeOtaisig7G!RSslID; z71Fx0K6b!KKnk}t-E@y-K$Fn(Xgud7pm%l2@3VXXn1st8^4R16me6K@)wzFw^=o?0 zW7|EzdfdMA!*~Z^Tz%dy@Aw)}^_Fs6^cg}HES9cbH5@?aV*OAGjxE*>U&w#!zPSyj3gj)ddL8cF2DwWFjd_gkLJp0vW#aiN z$X@ZTOr;zT84L6$%9)lRO=waPodzAG&`zHle)kJdhD|AS+B^i*u8x^1;mm;MAoEYa zDjCpmr!Tb+UIBFKfk&eZX93OkBiljUf`HP4@iJYH39>QbI6CL^4jJIHRhXQ#MdvN| zEsizbg;XwN$5#x~kR?=9delFi*bqBj#z+|?{>`e}gd7hB<|pj-UumvMi(p9`$C%C ze=)0D6Og`r6|qP&hm4Y6Zocr`gNzE6^Pk&@kU^|NcaN1H(z$Z$Fh;0A8VP}>gIY;Q z*&MTQXEh2v=+SKu{_KhFT=0_$2%13FPIO+o|2hF#J8-1js(6pAvlHH$%LpQC#ynxl zsrQkUGi8qjWUY_|+fqiC@z2P_cAV|^Ph3cM6gs_M=7F}xFK|Z-;?XHe7^^a8iLOdl z6Gxe!qdPBiCWYSmK#IWR&To%y(C4kb36xec5g~NTc!8@in{vscX=T$|vxyczwDJYjmIN5_9d|j~%k&-}njp@$1R`sBpUmeGqKKLO^jW;GKPNB%s&b_y1 zW{Z%e=h|JLeo-Tf#x;%9Ub)C@B+X#<0xdGxm|wiP=ZXxC-PU$AHbfd5?nJy-%w%|weGc6|Wc;UDcN5)u7T{2~_#WMu z=3I0)4Mo=$C|_&asiG@GQr7}RtB|P&;p$Er9>^3gyOJ^Jm6f3zY^OUJdKJ+Pgiz737l*F* zd~MU1a75RZ&YGH$ucE7Owu*OeSE4K0fswOnz35Wml+oKke&k?sQl|Wo4WLM(SAG3` z2~gg;6s#4+2&j%e-%-Dm1gJ_L)5~?=1yof+>ax1>fXY0(*nmA7P`d95)x0Z54u2T@ z2>CXF>@Evny;~c|<`Lh3{?F6Ma$~RN$TJ$`?@Tx7;W>#6IFcjpv_&Afj*UL2kk@Dr z$L)JQhWhAOadMfPPXqc7De9aRb3~U&E>&WUtLTdB(+Ikk7U=58DT%--A#|;|zg3b` z3Q+7QoJZnU0F~$P=KJwJ=s`lA~cm(smd z&Ctz0rH1>}8t7K6r2m(p zy^*dTVB{@1ln*)w7#8v})<2K|{US+)+;0fz80P3gosNNX9(q7OBvF~m2<&dgBfk7^wN1kNJl40eS0PObLi4XiT2CAzUX?kpkq`S z1G>#^RiCQ;1>JR4m|s+{L-&%5^hl~V(Y-Mii58P4KqLP_jda%z&|T%DrD&A{3{N4G zv92m$ddl^yepnK)Adk4;A0q=+^^hd}j1jae4YEb~0d4CsPL2W;IdDyA_1C_HERnMOoOl4z_^QJVR;r_u zu8USfbG7KY7Z>O2P)~Hfgz)vER4$~fIl?3&$PTHm3sRJ9-GMZQnRChHdywYj8STP* zNszjq^x&(49-!yCVu_P$1xzXBsor>Fz@{$xqFc!u;HJNKZs$t?LJ}oP(WnbJKG*ex zw}k^v-!|@IkrRMZ?v3olMvU|KGoAnS4jJIQnWc(uv;xK-Z6<(g7Eqp@%*qvOLY89g zzKr_&4y|r^f+}Iiv1P#%`*d4Y=>A{*JJdPckUFoTv}NKXqtv4yWWDvw{;}Q*$YOcgnJ*AK-`$IQ!;0l=7-Cz0I~=qE)_~;hi}F?g zAH`8#)p7-JM!9YcRheP=ku!XcsWh;CyX>*2DP2H-Q~RYV_BjeYHBUY3kvw>3zLcVF}0?{DAs+R*&-+df*e$fK5nEM=B zeY+2F1HP%w_1%!bIVCy;_#i<=K3?D38shKKGz46shq!=OatWF8kWE91qvyzffPt$? zt-Ck@u$I)v0$DO35Z6^OxE$bWW0p&)5CFVo{FZX48xVAkO}M*~07zW(KPNiyfb=5R zJM|_tAkj<<%|_(_{@ze;FIWXQf0!p9hC2g#)YqPq*9IAlcg~{t4bfe_#~x|9GmwSP zuVQe984~tI8>pqzA=fjfBP*{7kUQY=P)>|2G)dBRn1>9&j&p zXw(YQe&w&?)(L?e-idc=eqcIy>qY;ovi*=>hACW1FbeV?wP#D9%7nZf{%09CI3Z8! z$Iiy1B9Mzkt@=432oe&1IA*3818UEse*fYs0He>}H&XSO{+(n`V)S$YTzcn5RexOp zd@o&=HwUo(Na;4_fGSZSn*7y&;-*(>Y z?*%jpchntT2%yu3Q`d%t$06Iy3t5>yEDt65p8dxqLnt)SdR)b43X(3qt#zn81BDld z#J?C2un?sPEiqOa@@X>qpZ9nGxhGy$l&AVZj#dE*gYR<4zKl>>vtk&aA*|OjNk{`` z*BYUOg%m(=nz!P#)d4)K!UKETwm@*v_IvfG@f+!z%386*D2tPzq3#WCWT{w+Xl;-C5X3gp?d+pyHzgq)`tKeNvBBTIPG zx5IW>$W{Ui!P@38pi({>PWpBpFv*|j=b(EC2t1D)%fJrcxn&i8!+!_x`{&dC5tadb zo?rQoR{jLI(H1F@#46-h)?5Vjm11;@FWVrPrxE>eN4zMq$rkP6bxonDQNyCZPYfE0 z7@z=!q6e4hFDUw$!wKWXK(SA->jh>k5q#jxSi4gR zc7 z$gTH#km2b>g33#bkJkPOJtmqCS*J#_mGCHJJwL$du5u4DC|~|`bI1YRJ8!4kV4Z?) zsJRWKeddND>Z@pUh66OIToc)GeGW`wZcPupyN314jue1nofsd1r+PWu2G9!I=nL~? zLAobALN=I=&+*_Tg~@ey$c59!C%oy09M5ADkCZn;{E3d|wZTo0wXTx&_7WXrG;=xS z`CbSz{V3t4whx0irKp;*lwXh|c-oP18Cy@7e17^~^@4l`MmFX8IgsD=SZK)XN5~y2 zInTMCgDy>o@CEo`^9CbL)*jvhfN&~}*^2%F#?!TPA>pgg`2Dde=N~Acj7wtbyQWZxf&Q~h4S*tN-!l0I ze}SUOk{3xFQjjFwUNN(D2Qp~g5py*?2ZW>5YK>2&0ySxke}f1w;C0K}@d=NDNOR7` zOec&F@7G50sAmo`r_$HaGB^Wj#)K(%rOSY>pratll?~b5XcN&na6%@WP?NMrC(x!a zOZ{WxztC|%qR`nBq3A3QM@g`+Ko_(>_t91gKn|OPOBpO|kmuH_r_@Y56o@corrdc2 zMHp-Rk7l2RVs%OP%Czo5iIXV~bT4sGZ064KXG!^x_s?sJ+a0fw=!iQMb8jgCN9^H* z>fcUu)-g19$`I=r$hcLb-!p~HX`U4wZ>K`WM~&Q6yCC3lcy>zVVItsJiIj3FcLN+A zK8YuuX#=JSTKN6;89*__Z~eo}3|Ym0UTp?R$e4+L-O$BIZW!T_P+&ktFT09?HK^oU3D|O7Lfv8-r+BCQUVmDc-UzbH3LbX3nWiX4MWkv@xka%(omeq)8>}r zG!zR)zx#XWLXsYTaZR2i^_4lZ4aFRVX&+EQ1%88isIE# z)m_5)*U7KK#>Ihz#py6vu`nPud*I8L(GN(!dU@KrOaNcCdh>tBUI0%0-$gymPJneh zKqc+cc|cpLQrAIb#rSz&!g@Yox)y_LALXYTkjLpwSRN7#`DSLKKN35k(B=Gz?rVmS z^yGG%>R(qVg1k!qFCJSr&vMAk5+6Xmja^FaL5wufLb&(b)D`{rT#-MfGZ)Z#dp^U{ zm;)mDif~JkDG+HrX=CJk6ppaADqJaG?uRSfteE-QrHxATfJ1nhJ5%81FL!t^dYi;@270dnS4 z3t_p_5HLD&G^YM;0mRUqKnhDfK$6Vu!PSp{15Ubn3QW;9MXe zH6~4^QA`2>$EmPq&yoR;a>*<^1s5QsQrFJ9#sO9{jhHXz6#(5h3*TkgWQg0BrtiWy zjfCyU&yLCtkn=~H*dRX>kBk9QM>|fccnXMI zc{BATa1fB%9J_LloB~48x@SDCv3+2=m>z#F##Q$(E2^i5S3 z^S!LspZvz625DGpC0;cIq6-JfM_k6!kvWHZTwlwF0o}KZH_j3#0HLa2!MEWeAckp- znVCNW0;Ix__S<+M80)s~^BmLb7oXgJWvdPN=GtD~DryEqj;eFIe=uI&6w{qIe%}C| z3~l#(&H~m~nq68w1ArmL`#-G+YRGhpNA`y)4P>!5KB*swhpdGpOM0wtgN^2wk(yWn zWMw?xPQdo3Oo{~|ZGW90E#rJ%;qO1_PAmL{|1Sg`^8I%3@H-WleaXkHHHlL zD(a;!WIzU;OR8N$tdOoz`quU5{E)gs?yB+ceDpxOEH3lZ4!Z2{vf!U;D*Be=Tji0I zyU0?Dt`+qn6QC>~9lFI71DMRTWF@y=0Nl%8cMfR#0f%(=cjF^Hfb)iVnhe(K$64tV zolWluIH@6J>sJ{-cwHV*?f5cpR9NgQAbTd|7Begf$J)-BgQ5`ol} z?oBWBg&~cWp5zPOqmV{@E$jw6hi z-f@QzWDximY3z&W(s{9hESWRcoXrhG4$3bFI_+EqbW}H@;piM-&UzvsTv`m+UntJ6 z(kKHweM{`tw_SHfO_iESLQk$bT7q6-BC#n-7{efb&tP??uz7U&`g=5 zTMLEHD(F_w^~j6;)t-vz(p&YVh9m+y$+g@5Tk!y%JVJ@-cHg|L&ndqFG(?fpM>x#^m2sD_#G_SoOHCrb$0z{Z zOl2j1={G_*6z4h`EE3SwWuYOvX-vm_*UJ8W_zyb$&)PbmY9IYoXIeCU`wjY9K=|X} z2W#YqvsZO@r8V-8iS0cmy&-Gj&3spd#F3r$;1}M!J%GYY;mMOT7XVdO;ryQ)#{tdJ zYHybc2xyOR((LVH{l+1UrG@lK4`?NfY8W1@KIO`h4~cAMBh(#ap3P*qIF| z${Vz`RBh3bp?q|4y6PyCQ#LxYw-aio#E6cC6Z=g^ z-O%sz9M)>_PDrk}mv>@~D>5L@#rtS#1^L_bi&50q2w5%;e*T`G5!s~ctCuN#hWQI) zRdS1qk;BLL)PJtL1(e!1pt>^fXZs>Hf7*>X@rNLS z{au#X_fwXD!thpQcoYqCVEr=EU-KsNaqG4_FY6v!R=*Rdyzv9=*>CQVu;)O>zp8v} z&?`k3@!xJY7h0j~P1m{tJ6zE{_aL4u?Dc~iZkHU}Bk z;=WQ}^n%PqzRDV^cQNjU;;*|=-#n%|N5K)O9hN(H5zotNU_DmONwSAGktFhOl8XAy#jZ^%i zl26dq;+X5veWmEA4O2hyAQfF8_aqI~6r&qEXO;g(7@>#izT2$J&Y1o-vU)Z!7}Bo@ zK2NWE1(~0obh1a!LAJ8f_zPVQ5O;oFHoY_$;#Rnw*Sd=#ZsQ~d+1G@)n=FRe)xMBj z;W^!_YF5bd>Gb@&MOw(XeoP@PXdBY4e!cP^hcKj`-;X@z?gA-kH>MayW6-^)i^_@} z)X11q{ypmVH<3Q6h=z-C*N`Z_0n$!mCOSw~I`)B28l4-Qz2w~yk8akj-(++lVgIjd z7Fm`I>HOC<&3x(6Mn{8G52zi@3GMIt5`JJ5@Gb;$@~CxI{b?X zuU{g%`IuHbpiUQ3Sve)jf4PF4r*z^a7XIJ++Y-EGXaPA~7M$^%qL6Ew6!f$-9`ZB} zmM4F}IQvP;{iMmukZ4b=rd1J!JemzuEUHqF>yk%_=Z8qhA(_{un-L9hWoftBH*_%` zvxK+Ke>Wg~8m@|eI|fo8z0YrRwi`X@Oip##!Xaygv$rpaVsosYZ-2cR`ix8o3e9G- zcp)thBN;5|?a&chxc{fzA6-9C{Cq7j1X8O|hP`M6m|mCabKx%u61og7IA$M*JQBA; zb^S2j`mDRDQ8LBMGeCBhEoZufuO5LaIi->7O1nRpwn%N2BB*GnAo{3MI+ zU+mB9JBs<1lBKUwOvC}I)S7TRF<(F-MFCpU_7Y;?8F&~<(^?%NN zgOKt|c7-E865>AQnQ#tLL0<33juV9+AZdDQ3>gcA;$EE{MQ`L&?OxT@?=F)NMOC&*-Ebiu5Ln^m%6Iv zf`0JGx5IjI>J5xB#wvK&O;Kn(Rx8C2PE1> zyMIwhhxk{V)gJuGkfG+sUtLK8q`0-hEnjsRT`TV~YKV^mxXf7ED;?$lf9YkT{zpH6 zkIVZK7?%xj7Z*po!%qUH@_F#)CL6NsMTk>?*o@yrj7Z zY4QtAv?{Tl!0~x6*T4m+08U3HWllqN77O>W7$v9`wN#@{`xR>C+5R;M!ulG1y0}cc zwP3w}Kc9>nVf?NhkI9e5LQtgM87yA?47t(!l{)O0uV5^b@Kg(j?#|r^*#69oP8$8U z5h~Y#6o#qh+*VZuyh;xnw4@(l=RKC4p|~f2#Bf{HzH|?WoT0G3WT^{8J_g5jWzqwZ zG6O}P8zu^3ra%Tl$y)CC07Z#v}3_5zWM0%sC>WP!L-U_YfE<|p&8 zm#ik>fXtk`lYcWSknKvvsZ44Csjm(Ejc>;Qze-E;Qf@3{DQ0F|dyj`I2OZtbN8_M@ zzoG9k?;bSnu+S}-4TVNk>bk-w5}@8U(Ieq1(NL}SZ`*R~Hk7jGGU5-ifn0+rc?Ff3 z=#ImLc4G=28Tot(VF^w^cJf(OKYzYpZ|y$9tU7Xjt$@8_h;nSmO~x1YM#7^n+eIo8~1f^qyFgzJZF06rSi zQ7?{7s5o5?@oJy~AHcI0WJi6I1h}U|t8e=y0o&i5x6`}| zK$L#avKxC}($@Ok$(nH>_j~P|a*RDtjsFwS*oPOY>hJgQVY z5=J6Ju2_3t!;x>8@AYK#=X_xxeEehG6T2ZGqHvp1neH1PS;8uwcN|ESOgzeG-U9_jSe__yQ>{ncRqYHXwh8N?@$e5U5g(W?s^b09qgTs={nk zfxe2@#svEVU`+E_^+$p;u+TmqBW9`Jlq$yQ(cX(zjs9(mDoI zs*=cO1&o1K5pVmI*EvA{N>IjZF~&dUex+sg>^ZQ$U-AqxMFKm+H|c|pallD!lEVAp z3t+ME+(dOx7V775JdHLT0wxPbT={yXfd%pY7S9o^w><3cgqn8;Fc>ys<5Ly_I{&R7 z3@j}Ib>tfXwZ8$Bu5cSAC6oi1>YC~@7bzg_RvuwIv;jy?`0~}?KLE+SKQKE5^XZ?l zGAbDSjrk-`hH*N61xgo4KOGK|fx4u^jk%U@KMF3a+Obj>7Xj8~5=nl? zATW^~-S_D<1^VFx?e@|HpoRG$zg{;1Dzy=Q6sGy$$d8rqT8g_sGULMGk@gjgGf0Rs z6s^MMd^yJl<_Lf^@e2%SQ~>c@$GOI~5FjIFrZZi;1{6Fl#f_UR1Ld~$QbD0~psEzW z#+^zBlvI-TB9bYP<-+u39;Pv9q)IjUPSOfG-~M=fJNz~@9k(^QFsFp|jT;isW0!%I z)#Q;QmO{XS{+4uLZvZgb$dd^d-vIiGSLghbTCw|1;B?>JP@qcvqPx~q3cD};GeXU3 zfRv_?-Rb9HKs2hzUVHvw2CZJs4k2SGl|NO!UG7(YB>m|LMp1X|}cKX#7xgLWIt(pe?a(4ysC zowcDS#2MbEdtclKj^@i4Gld2M{hf+%>S|`7&lZtc&$Iw^NJS?kd*1=|WyPT*aSlNF zR`j`R1=QH@el%CeOcq-Q28P1>F@G-ugWXG+-+*1T-ph0i^LwH`N09F>0MApfeYGwG zu&d_$yNu}>H02fRd7sMBQGQ0hjiY^#-$R})I&~4MaKBHdio632q9}D0(wm`0g>Ule z^l5n9JakP!y9Ap2K4E$%wizlnJ92DG*aNQ5Nq?!-Yk{&?GuA2MHEIQv zpM0Y3?@}W>b2Gkc!31O^o7r@c#|7OznnT*oRfOFBA+v*D-b1N5Gc`}11E{V)&3y3i z5!AQaWZZi65E`4QiKvEiL6Z=ZZEmbD!I1Fqpw(j!svf9V)lCIL+}zqp?qkz{xAy#^ zR8Sd^3Mh1|@4p6Qi$@jd&6t49Z!z(GhGih-wjkU)`WA@08f(f^Mqs@Y*G}|C;DF%c z(X`5#Kge?T!}0nyEo9vEUC_V4P_$R7XW*k=A-dltHo%=&2l4I%7fn%Yzo=zg%E!bDOAFGn6++;mcRQ3b-V_QD0nC_ZmWmH;nZiSpI=xK$TAIedT;O6Zihe>g^&}+ z*Q=rMjc)bEJR;JIKt0CjXu$O^uuaXi1CWT*We3kgfC#!=uYaKph@8tOH>pbl68raR zcD-~UQ14>*<%$*H)$3Jz#h-~(`OG_=yfBT9Pe=QBmbIh17O6L5h#8QnIC|3DUI=nm z_JwhX;h{(iJEd@Z6O`E$$tsQei}hZgqPcfJ5Nfa|P##=bhg#zkN7#bQpw?H;A|?@g zsP0s)(E14Diq`vA{tRYR- zFV$&vS%^E+(pMvegT(Lz-6th6P{jSwsQMZnN}X;uZjF|J@+Ao&uXy^Pa{OVlTc0{q z#Xa&6yN2~Pf3jZSw=9E7H~fv(mw!VUlcZZ+cgc`sb75rOmIAW-Plr^W8bC)n9h;;d zX(Ic(f5(>1PXMN$GIcr{nE)3n97^kQ2H>AKpte4kt{Ri^La_G%z?B8fUue9Eu5K84 z$4Cr73W58VUj;Kj1{)cY%vpD=AM8*6g&sXf6f2LpKa2UGso0i3(9U8#y7DqrRo^kb zXjmH&d)~_Zx$k647*FrDUV}$|jZziVIzyS-&%4*-8liZ>^HV~(@lcSrGxxOP zCx~Y~NRHgSg!v5_ibZ!!kbX9a8a0M)WcPPiSdH^rKr_wr=)_-4pL+7Z$A;?wFmc%b z6H0jv7!O&;J|u9U+bLX0jwM(Q;g*V8zYXS3`H=ZuM%oVIUA17-;wI!SEwstL77hg- zmA*>o?Svv84~z!|iBMcg{&Zay)`MkUwCR6I5=to@9SE8xL&;=9$kiSPDDEt`X+$j# zNkfX|vnOdGujvn)bxtpc%P$%H{q+^3j5JhA#XF*7tZile?*7Q6MxE!oR06VPP~u!& z*bFEFE(fgkUIJ8I7ar+0zrp-9U!D{{J&Eqs^`tI{VY%GU;iJcz%OPWBPA2|bCB(%# z`*1sjLe4V%uJ=D~L0*3LpzOL+P%y0HX=|@0Bt1<@i~l_ZMRww)9n!ILWzgoslyz*r z_?Oyct`7%E$FFGVvNAyas>|DLAF3f&Ep;k0w=iVuvwWghY5}QcWDPtcW6@=)r*E`e zywSpQlAE6<5|F>$CuQ2Jl907WV_w`MyvUC9_4v*s*~q>bHTSO%O6UPm=zg5;2S^=3 zcF3^eg$&h0=cGfnAX}lzu~#QiNJt1s@UVLfxdt=`G_R2%aht3C|GEKx=DF5{2tFui zAva**EdRd%C&>E@4z5>nMCX068z4bxs*Fo3&z1{6(tL3xwJ^`qrLF9rJ?$ud+!7~Ez5Q>k%jfwmz5(i!M#rdoz6Em< zOj!-jz6m(1v5CR|z6}zs*Lim7z7S)Km)aAsz7hsR-@&v9z7t4I5{nCAz7zq8JwP+- zz7!HS>Y_Bqz7uLFM_Sf$z7h?V*g!`Wz7Q40Cn|sQz6)8C*+Hk9z6J|_B$7O2z5u7s zKO|p`z4r{A%zKkHz3ZG!0LH;Iz1M(Fy?i}?J){^qy$uSZJ=X4E$KS6TJ>wqH>@G|S zz5F9rP!dbWz5u=_G7K0rz5}1S-SYj=z6PIM+^&huz6qNdzSo6Wz6@qAuk@e-z7ES! zNy#XCz7R;zc04-qz7WD09+bS=;D9z6>p(4~vx6z6(Ss zbiw{kz6qD9)(%^4z6U@=_~Y3jz64~J0V%v&z5z5S;XoF5z5Bt!4+F?xz4bA(5S;2o rz3vB_?d-Vvz2e;eelXD(z0yf^Y{jK5y`Lh_3F8@1J;td-W*O5w-VePN literal 0 HcmV?d00001 diff --git a/test/_correct_static_disp/correctness.sh b/test/_correct_static_disp/correctness.sh index 1862674a..1ec3c7a4 100755 --- a/test/_correct_static_disp/correctness.sh +++ b/test/_correct_static_disp/correctness.sh @@ -12,6 +12,7 @@ grt static greenfn -M../milrow -D2/3 -X-3.1/3.1/0.2 -Y-2.1/2.1/0.2 -Cp -e -Ostgr grt static syn -S1e20 -e -Gstgrn.nc -Ostsyn_ex.nc grt static syn -S1e20 -e -F2/-1/4 -Gstgrn.nc -Ostsyn_sf.nc grt static syn -S1e20 -e -M77/88/111 -Gstgrn.nc -Ostsyn_dc.nc +grt static syn -S1e20 -e -M77/88 -Gstgrn.nc -Ostsyn_ts.nc grt static syn -S1e20 -e -T1/-2/-5/0.5/3/1.2 -Gstgrn.nc -Ostsyn_mt.nc # rotate to ZNE @@ -32,6 +33,7 @@ python ../compare_nc.py stgrn.nc _Ref/stgrn.nc python ../compare_nc.py stsyn_ex.nc _Ref/stsyn_ex.nc python ../compare_nc.py stsyn_sf.nc _Ref/stsyn_sf.nc python ../compare_nc.py stsyn_dc.nc _Ref/stsyn_dc.nc +python ../compare_nc.py stsyn_ts.nc _Ref/stsyn_ts.nc python ../compare_nc.py stsyn_mt.nc _Ref/stsyn_mt.nc python ../compare_nc.py stsyn_mt_ZNE.nc _Ref/stsyn_mt_ZNE.nc From d7bc395fdf1f08a739d386d78a8fc14549f1ec2e Mon Sep 17 00:00:00 2001 From: Dengda98 Date: Wed, 31 Dec 2025 20:24:25 +0800 Subject: [PATCH 10/12] TEST: update _compare_c_py --- test/_compare_c_py/compare.py | 20 ++++++++++++++++++++ test/_compare_c_py/test_compare_c_py.sh | 10 ++++++++++ 2 files changed, 30 insertions(+) diff --git a/test/_compare_c_py/compare.py b/test/_compare_c_py/compare.py index 7efbc583..a2bd130a 100644 --- a/test/_compare_c_py/compare.py +++ b/test/_compare_c_py/compare.py @@ -198,6 +198,17 @@ def plot_static(resDct:dict, c_prefix:str): AVGRERR.append(compare3(strotation, f"syn_dc{suffix}/rotation_", ZNE=ZNE, dim2=True)) AVGRERR.append(compare3(ststress, f"syn_dc{suffix}/stress_", ZNE=ZNE, dim2=True)) + st = pygrt.utils.gen_syn_from_gf_TS(st_grn, S, stk, dip, az, ZNE=ZNE, calc_upar=True) + sigs = pygrt.sigs.gen_parabola_wave(0.6, dt) + pygrt.utils.stream_convolve(st, sigs) + AVGRERR.append(compare3(st, f"syn_ts{suffix}/", ZNE=ZNE)) + ststrain = pygrt.utils.compute_strain(st) + strotation = pygrt.utils.compute_rotation(st) + ststress = pygrt.utils.compute_stress(st) + AVGRERR.append(compare3(ststrain, f"syn_ts{suffix}/strain_", ZNE=ZNE, dim2=True)) + AVGRERR.append(compare3(strotation, f"syn_ts{suffix}/rotation_", ZNE=ZNE, dim2=True)) + AVGRERR.append(compare3(ststress, f"syn_ts{suffix}/stress_", ZNE=ZNE, dim2=True)) + st = pygrt.utils.gen_syn_from_gf_MT(st_grn, S, [M11,M12,M13,M22,M23,M33], az, ZNE=ZNE, calc_upar=True) sigs = pygrt.sigs.gen_ricker_wave(3, dt) pygrt.utils.stream_convolve(st, sigs) @@ -248,6 +259,15 @@ def plot_static(resDct:dict, c_prefix:str): 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) + ststrain = pygrt.utils.compute_strain(static_syn) + strotation = pygrt.utils.compute_rotation(static_syn) + ststress = pygrt.utils.compute_stress(static_syn) + update_dict(static_syn, ststrain, "strain_") + update_dict(static_syn, ststress, "stress_") + 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) ststrain = pygrt.utils.compute_strain(static_syn) strotation = pygrt.utils.compute_rotation(static_syn) diff --git a/test/_compare_c_py/test_compare_c_py.sh b/test/_compare_c_py/test_compare_c_py.sh index e74e72ed..c489b8d8 100755 --- a/test/_compare_c_py/test_compare_c_py.sh +++ b/test/_compare_c_py/test_compare_c_py.sh @@ -45,6 +45,11 @@ grt strain syn_dc$N grt rotation syn_dc$N grt stress syn_dc$N +grt syn -G$G -Osyn_ts$N -A$az -S$S -M$stk/$dip -Dp/0.6 -e $N +grt strain syn_ts$N +grt rotation syn_ts$N +grt stress syn_ts$N + M11=1 M12=-2 M13=-5 @@ -89,6 +94,11 @@ grt static strain stsyn_dc$N.nc grt static rotation stsyn_dc$N.nc grt static stress stsyn_dc$N.nc +grt static syn -S$S -M$stk/$dip -e $N -Gstgrn.nc -Ostsyn_ts$N.nc +grt static strain stsyn_ts$N.nc +grt static rotation stsyn_ts$N.nc +grt static stress stsyn_ts$N.nc + grt static syn -S$S -T$M11/$M12/$M13/$M22/$M23/$M33 -e $N -Gstgrn.nc -Ostsyn_mt$N.nc grt static strain stsyn_mt$N.nc grt static rotation stsyn_mt$N.nc From 7c5f1d648f0896e44a84b10fb93c82263f41a6fb Mon Sep 17 00:00:00 2001 From: Dengda98 Date: Thu, 1 Jan 2026 01:58:28 +0800 Subject: [PATCH 11/12] rename --- docs/source/Tutorial/dynamic/syn.rst | 2 +- docs/source/Tutorial/static/run/run.sh | 4 ++-- docs/source/Tutorial/static/run/tension2mt.py | 2 +- docs/source/Tutorial/static/static_syn.rst | 6 +++--- pygrt/C_extension/include/grt/common/radiation.h | 4 ++-- pygrt/C_extension/src/dynamic/grt_syn.c | 2 +- pygrt/C_extension/src/static/grt_static_syn.c | 2 +- pygrt/utils.py | 6 +++--- 8 files changed, 14 insertions(+), 14 deletions(-) diff --git a/docs/source/Tutorial/dynamic/syn.rst b/docs/source/Tutorial/dynamic/syn.rst index 0b261c40..0c1f833c 100644 --- a/docs/source/Tutorial/dynamic/syn.rst +++ b/docs/source/Tutorial/dynamic/syn.rst @@ -126,7 +126,7 @@ Python中合成动态位移的主函数为 :func:`gen_syn_from_gf_*() M.active = true; { diff --git a/pygrt/C_extension/src/static/grt_static_syn.c b/pygrt/C_extension/src/static/grt_static_syn.c index 5b39b466..66fd3510 100644 --- a/pygrt/C_extension/src/static/grt_static_syn.c +++ b/pygrt/C_extension/src/static/grt_static_syn.c @@ -209,7 +209,7 @@ static void getopt_from_command(GRT_MODULE_CTRL *Ctrl, int argc, char **argv){ }; break; - // 剪切震源， 张位错源 + // 剪切震源， 张裂源 case 'M': Ctrl->M.active = true; { diff --git a/pygrt/utils.py b/pygrt/utils.py index ad2afbe2..65cc6359 100755 --- a/pygrt/utils.py +++ b/pygrt/utils.py @@ -4,7 +4,7 @@ :date: 2024-07-24 该文件包含一些数据处理操作上的补充: - 1、剪切源、张位错源、单力源、爆炸源、矩张量源 通过格林函数合成理论地震图的函数\n + 1、剪切源、张裂源、单力源、爆炸源、矩张量源 通过格林函数合成理论地震图的函数\n 2、Stream类型的时域卷积、微分、积分 (基于numpy和scipy) \n 3、Stream类型写到本地sac文件，自定义名称 \n 4、读取波数积分和峰谷平均法过程文件 \n @@ -96,7 +96,7 @@ def _gen_syn_from_gf(st:Stream, calc_upar:bool, compute_type:GRT_SYN_TYPE, M0:fl srcName = ["EX", "VF", "HF", "DD", "DS", "SS"] allchs = [tr.stats.channel for tr in st] - # 为张位错计算 Vp/Vs + # 为张裂计算 Vp/Vs src_va = st[0].stats.sac['user6'] src_vb = st[0].stats.sac['user7'] VpVs_ratio = src_va / src_vb @@ -168,7 +168,7 @@ def _gen_syn_from_static_gf(grn:dict, calc_upar:bool, compute_type:GRT_SYN_TYPE, srcName = ["EX", "VF", "HF", "DD", "DS", "SS"] allchs = list(grn.keys()) - # 为张位错计算 Vp/Vs + # 为张裂计算 Vp/Vs src_va = grn['_src_va'] src_vb = grn['_src_vb'] VpVs_ratio = src_va / src_vb From 709fb83d3ad7955b521e929087f1ea8b2ea13d22 Mon Sep 17 00:00:00 2001 From: Dengda98 Date: Fri, 2 Jan 2026 01:49:57 +0800 Subject: [PATCH 12/12] update --- docs/source/Tutorial/static/run/run.sh | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/source/Tutorial/static/run/run.sh b/docs/source/Tutorial/static/run/run.sh index 590d0937..c9bb59c4 100755 --- a/docs/source/Tutorial/static/run/run.sh +++ b/docs/source/Tutorial/static/run/run.sh @@ -178,4 +178,4 @@ for pdfname in $(ls *.pdf); do rm -r $pdfname done -# rm -rf stgrn* stsyn* \ No newline at end of file +rm -rf stgrn* stsyn* \ No newline at end of file