From 3d40e88e35a72a9220e42b32b54d038574dea022 Mon Sep 17 00:00:00 2001 From: CrbnsCat10n Date: Thu, 7 May 2026 15:56:31 +0800 Subject: [PATCH] modified: __pycache__/evaluation_studio.cpython-314.pyc modified: evaluation_studio.py modified: scripts/dataset.py new file: scripts/inference_utils.py modified: scripts/predict_single.py modified: scripts_2/__pycache__/task3_identify.cpython-314.pyc modified: scripts_2/task3_identify.py --- __pycache__/evaluation_studio.cpython-314.pyc | Bin 47295 -> 47537 bytes evaluation_studio.py | 41 +++- scripts/dataset.py | 6 +- scripts/inference_utils.py | 72 +++++++ scripts/predict_single.py | 12 +- .../task3_identify.cpython-314.pyc | Bin 12782 -> 23732 bytes scripts_2/task3_identify.py | 181 ++++++++++++++++-- 7 files changed, 285 insertions(+), 27 deletions(-) create mode 100644 scripts/inference_utils.py diff --git a/__pycache__/evaluation_studio.cpython-314.pyc b/__pycache__/evaluation_studio.cpython-314.pyc index dd109a00ec8e027c9303f9794595ec8f811ec18b..9cfc1d5fd1ac073d265c1796f9110ac14960c93d 100644 GIT binary patch delta 8398 zcmbt43wV^(v3oZA%D(d2Wb-1~JaLl{hzTSTN&w{%AP|1iAH;=alWfRJvb#OI0RbiP z!5ak`Cs!x5SwpLm#Mio9AouB&Fhfi5G8HKYL(|xESc%<#3c%QsRc|YBT zc&N!JU_XY|j7oi>N>&%C>fnPr@`{57nT;h$_;hBgRQW0?NjCVC6^u$H8z6P84O({^ z)vhc!C*{NC42vb%l^0`JasW&t zSs4>Q@^)20{;YiXmEDwNMo!f28#sSxcaGH;sg}IDrXB zISnc0h#%qD!%E;9d!Csv3I#3nPpxq*jWH?nDEha4wM}rPWpR!E0q5B9Hd$Cjw#kMs z*dgnMu}R%a6-VPY-6wy02hKL6#rd3t>X8K3~qPQ!5gJ zS+e${uFOd_LoGaPov@UcnJZi};WnW?)nM*%e%1%%$vyLf?m2wiWS#fOc zmqgj)YU|BXo}4Y^%Q?~*IUQ!_moqK+^P3_swL8M zM8(4uCzt+M`suR&SoX_`)^n2kwTzVm^CQ_Mdsq>vxZdawR{J|WoO{*@brp4JIeZJe zQaX=47_pblH?YT{bKGU!6NG;NeASI?KRi%9Ll@QfJ-uFduX?NVVWK@k;86l`10*c) zF(Mp?!{ht4L)d&UQZ%6_$^C0=hA5H^s4l6hS@DMi4rngY>s-;5qRYo1y<>PIv8htg z05k&oNY>EE6)K4}tLCsEE2Or%Sv@Bys4i5w46vd+AAWb85p;DKED0DmR?11S-4ALc z@vja=(Nm~es*$vb5h18`r71LwIT(Em(yl zV5>2|E;gkrD<)1dN#=T;(p!r{`3H;<7@^3<*!)sTEObKDynu9!WJ+ioYAP}aWx68}j9esR_cO*>l(=K2aLa}vOJUTq&eRu zc6+FU`-0(6wa<^rdHo$;cVF1o6XNed`&Cs96fi|_=6@n?bE~v)?#!#(+|t5@ZTu_| zqlC9kZ*6ri{B~RGx7*y;G)sI4-neSa%oN^395fGqkJxnG9)D+#muC{@eZm-pl0JX1 zFUdgA+=?Uy@scuhx{7*p+8`QVtCDOH4SiR=>&LYvia&? z`!|_I zZ)WEF*`m+Y{X&y9F#BSvD%B2}#y{w@?=CsuJ5&vCH7=;ju~>Bjvo7YV%;|gU4j5qS zIaAR&P0^)qGnH}d8_DJ8lFL5{Q7m}aKF6;8ncXz!ie&h|)<*UYoU)cObL6K}K4aI# z_cCg3M>DMTb%wj$-Rq;8P~PxXwJA24BE&3;6rc+uX#~!eh-2_FUw@nk$$=XkWy4 z?wMDVSwftp1inMc^bsh7N%M-E$_XjdOrn;CfR=!cKq)caL4aIM@hdI&_*Z%PenK39 z{qqWP4q#L9aXhF|gD=$O^ZUYHeimMyH-5oa&P~F_)FW1)lgy1oh`S=ck0cwpXRSMg zC5ER5^YrZbkLZNKO2D%eV@Zn@|Fl?5q-hu6;G#%f%R5Pm-P|`LBA$h7C<`MvfL9@W zohzFpx?j%1l$J4dvX4Ihm$PvHzm|n%hJ|V>PA8sNq!ej$RYHs&ouLiG8QSoL3{6(k z3qg2QSU5Gy>LScALNdaqbw#yk4|-hWjpD)FZyJ@Z56a1-#+gNKn>*2xWRi7zoHna$ zN~GgN_S#tOV3rNB8;&A9Xpt@PoKoW=!-{9gt0JYxlkNybw#2L;Sz{R_F;z<`qf)xn z43e;Qzm#Y+W5Spb)8r&zO>z0T3A2tA26+j*ah;|7 z1gx%gY)o5x$7;0%26CqYXv#HQiA91nzM+Y_S>sw1B;<6%NR2RJ5KU_()d-DWStBGY zAzsCdw*=GVG}pM8l5#4R3YGop`L)BvLJ21ibh4V3tY_{xOzn=^ze#Aw!6zl8O+uG>mXluE~(W+D;SV_=DAYHi%DK~+?kd#aCSGgdU z!gDl5q|kOx$m?F`=Di{Q31XpgN%8VeiSsi8I|zJE;4*@!Wms*u_l7X{IPKB2SlX$p zj;&}BMMpIOO5OZR5gU#yu{lSV2V76=4Fuw)ffDbB;I~V%uSq5fWhiI5lAOdGLj34h zjMuzjyDp(2CIVC)a5IDPOS5y2Q`Mb9C@+8wMj#V9~gReH{w14cSl~`7NS@YymJNyH$8wNuDHR8 zREigUbb}r`aAI)x_lNG@_V(TP4h}qrq|e}Z$Jngx#6WYubK=C{Q^&6;qfwv2A3NSg z3%%A^n^TESB_y9e`sDdrkCEhc#8U-@-YPZ@=6iEx*wMjzh8NjJq;#D9{J<^eH{JerxwCxV<8dT$QOB6A?f`e{_^|@Q$)2sRAb$G!;~r4@IY|tRHrsU+yf*| z;e|JEtjZTM4Qo8y=V|ZpI+GGu4}jyQ>S-gIh*6}DQ_b%teYN63(r3mM4n8+9_~4Fr zPLL~)5^@*ZaZ@+j2bXSYV)sM+Mu%w+X#fNcLi@%U^X{YO{Rq~Zrd%EHcllP~J$dQ^ zoP-MH1b5?gN^yNA5f>5(}g9foyGB68Q%bf{S#7O zCye3hKyQ1XJ#YiR8QlFj_~~$6|LUCENQX@bqGmVm^Se8Gd>yO3oqRKV*1v-71lK?v z+W>bCjA7g1#6Tgt9QpZxtV!^L=`e58N9;yu+PuQ4R16Epu3{xYcZ+as4DHITY;ig?oIVFekzMR|Kd`j%pF* z3!=;Hgq&?<>;-7rR>0=Mif#9$(%uRtZm0SB;$FKLCT)K$_azdn^Z7$w9wtq_F-`YD z$&P7kF?8+-Cn&ZAKHqU;E-fXhb9>j|9aC9Nxw4wGuy%o=xwTC>tLHq`= zg=ongF&}XzD~46Z!MAqLVHaTY?km``$cf$a7~V%?r?flb*z-qbCoP1x4jsH@X#bJH zLuUqW-NpX|#e2(TVXKdj!D))>yf>hsL;Pjz^1l(FnH7DvQA=!zv9&u84m67u(N6FS zL~2A3v!hxK_w21_-I4#R$$y>=(0PAP^?M|jaXlhXyr|&5e#j~K{QmyZ5yR~u;rOTeqU^`v=mUFjFY3Nu zQ&>O~el6#daNGB{O`uiqVxdD+-|G#9$fLhrkEh^U5AN(Q|3<-26Ll;Bx`Hln+LQ`- z4tBUqpa&HfAF+5Zwdnf`_Y*=C9M48jk{1)};syb*MS)`&xfUmU~utBKYh2;)yCqUdbFfx^K$Rx9=wdt5=9mmmny zH@B{Gd>rDLhL5t7Hj-X2lYV~jgyF(5JKGo0K9r$OI26yq)JOX4@8Jx%xukT|(ChK} zdjf$V-vKW_;$+XlpC2jymNhQJF1-XU<7fH=kEWBf(xlIFY~{`1iaBd0!r3-Fgm7pb{YBHaAgbap24!ej3; z4CTj9^s=iV|KvDq&plboPQ%8N6~<hAgT|q_W8+g`r+&o%er618PFg4dV}l4_qC`lh}DWG z#Agt};C6TUdc1D8(2i*MZ8)qXLva!$DF~~cu2Nr{3HzR|_&bWzl#c)a delta 8315 zcmbt433!x6mVb5nzH)STI(H|KgANb~r$ppVAj1`Z;ty)1=_H+`q0?!rJ6xkCC?awQ z`f!I?QNh6jl*{(MR(EG-bl1;8S=X6fk99u@GLA0jAg(j&zE^*`(}2$Ex9!LKtKLnjIEjf-P>yQn( z^;SvNs}d#HG>p?tV5+_@Hl$fGBg>t=PA%F}4I3&P%mk`(i_WT8o5TQy7}>t`S)$lS z6(p^~NsoqoaR;k`>_U%aPiJu}oU}{?gMEv?^^) zh9pHYBUzDb*%Zl<%@IdETR1_nG)){dXr|KzD{}2@GW9d8Dn5e@KI?M8VWbG7OJ%#9 zIT@AFFXB|3vh4z6rwDOVQ=Jeia9i_{qY`;9;k*J@r9jg(ojz8gs$A%*zvY3#lm}46 za>39);0uTSlA(Wqv(}}0O>rrvOZ=82=R_QewoGbCdAn%_#*tgBOl3}T;Ul?K2fSWk%y%p94iqYq6vHLrE=$&sLa}^1^Ua016EwoXQpBz3ew@t6T*aigWyj>f z1&`#l)T=ts_ZdGzM}I-lL{p?daU`R>y-wJi6PQHb|3h2OMYHHHC~y?W&SsPB8lt2* ztjd;q5&ru;VJeW_N=^rc8*XdU!iq~oE*^GgzNt=5pLuR_UF7_<8Msr}f=g894Pge? z!(^2Um3+BKDUgd52RMo&@XvO}O39xqqFQglZ(sly8P*<@4FW z`2O+*M)m}p9rd05X~I7YUyi<5<> zeY#WFJQzPRCX}({eQXL3^h@WZ3O4L4f!*3O^f_KOPINhC`L@B-flW!u$A9P@_QfJC z`DnQ$vqotii!hJWqJ{hJDALFpcy4B)P19l&0tlMSR5yH@Wbi~XWbxmG!F1@6mSxD= z^qdsYw%C%=3RVqSKRLOvbuF1mN=tUKrD&TRlt7I-*&y`5*YV++8U}APc=S5ih&cFc zbumQdm@@T>K21;uebaJb#u8hZLD9+P3yd%-#=$L=v9Pu!A8M;To}t0X6XK^BlM|CI zaAt|C(>h9*q+zZ%wWNv?=6p*UC9Q}ml46$aXhCMW1=!C@@2A@>XTyeRPWuI@96~Kf zLm85mDMrk2CUxwHMRv+AIk&}=v_+AUd~!~jv`!WZ?XoEyNe=6i&1t?dr;&$z;ZDQq zQ~LALuqo^0e8rUVbhK4Na$ldp_ z;1uOvo9cs4h@|iGhucE|?jq`^gfR=Fg5gMSjDJDc!CW-AZj_F8cqOc_E6y=vh7sxL zEo@NYt%dKFjfHpXtkd7MXYZ=mt-Ig$OWWS{N2}kiy5gOxId4|YIZ@ScymIb|ym_x? zznwS#M9zX&S02CG|MsHRcNWFoTogND@9m%a`^;PzdHH`q?$iYn+*Z53zu`=QWX;(% zVYd;c9k-Mn*Or~HXOg+<9aH6TQ{~rD+MN5{^BlTY9hQ0Jcv9U`r{SmJn7y3o&!OOTyZUb6pc`k9A?~N@$~9J!|*?%6qMeyl5K_Mq*LlxEfma@xEYtARG&Jtm8JE zUqTo0^H`MADTB`?Fpt0kuy}f|TZ0V3ufpM8qv~86+%Ubwd5G9d#9%}arSw_}eRa8p zaXcRpdidArpEA|j!ub0$RvOuAXl$I8o3fc_lklHVZV~?!?rogT*1~&@rR=5nw~f`z zD6B8z+hOedGFK6?78AIhh81N(UG8fEPm<^liG7aN-1Bw>LyB2<&p}K6?&k>O-ha;#Q zBG`>rEc_g=n1o%pVtldlOYKiAdwoBt>?sUSoK~EN_$Bk#8caxB87h1il|;i(PHk}{ zWw9iZXCyE@Vb?1Ls*7;L&Jx(0WgIMw27@S!5`5q&5|t6_U7EZhSs61-iVoNl&Md{0 zPL|D>gI(j3b=iz#TyU?5SQJYtsn@p9SxfwuJ_5!WYgGKI=sX0faq0C9k6iZqmwX1~`vL?AQ z&1u{Q;TtCv&a{EZpUdTRucXTi&E1_QG*kfd(y%FO@NhmffDQ{F-gs?87}bsK5%g?8kXngtUtF16NBkH6jLZ}&$q6C_gxRQ z0Z~ZLQ0pvATUz3x1enMam8jJp4fxjicy|<5FCFGRM@(N5*hb(yfv*rGY=frw#7)KP zf>)M0`x+_Td{5F#WX_y3v#GgW_RXnpn0eKLW`1F<6`77pBrTHS8#Ii_psCWDFh<+_ zp@7dH6GfJ15VH>Xyjf}qz_BYOeqqAHyUw7)77KO)Ur2Drvb>U`G=!cI1_G1PTay&bZfBaHvaz#Zk8LiUi}@wHSN$kHLq~~rNy0df60(?VTgalYT~eXlUf6F?^;pF zss)ReSlCF>9);~()tXvu0x$h=ufIwok0By>TCageTgrawvKZ&%BS%g?eRxC#Nu7fi+E1}zus<+9zZ@e&XgJOTC8PsmR3_3##*?aLPLtklFI~?x$4SNrIBVV%+lx-^2g;D!w@w>0T)tvoX zG!m~RhKRqtJs9rfeK7KdE!XWJOBNxC@1)jU2*mo#UZRu^{M~IIJ^z%K@m++abH(ZY z+XuED7{qJ%cEl$PZ9U!Hc$)QU`H-?Zg)+#wvDGJ-^gIv>209XYf4HryhkNxyya{t# zPpG#$9F6$HLq;TYK{RyTdgIiob9hf=+BC0jNQB=@>Q#f`#Eqkh#2CiaejfCVh|rOx^pCCncfuHF_jI@RwDzp!>yR&D zr5N+gp0~JhHeNGz_2C#_$Jaxwe^vesNar^SCLiw&``SXmwpD?4eiMAvznpD^r5h%& z8(`apVsg$6uffVFg6xGx@}XM$|>=ko9D6d;-nuXxizH76pe&}F}h9) z#hmVQ;wF>OAu1Tb2-^sjEfwr}sNGV+ra;q{UuDtKtM^Mcjn|zD`&lqz>xqI_NU=T` zjs|#)EDaJDgGLF0^)xSwOsJ*j+4ZzwvN3h2DzB}hLHZ<|2@f`SsxyeGY zX$SuP=7C@AJ$cVdCvVw-w|VE?m5OlMpOC?8N$3M>(c@8JGXETLVr+FT8>lPR#@y7^ z6YFUdvlX}b0V0hj6@yFv7X~};p3IiUU%lJG(j3+mXTRvvviD)bJst}MjCwqI2Ohnr zruJLvmy^dRLisLg-9#<1_Ag^Ap~uIP0G{h;CAkFVJ(W%hc}{Unfq~ay%nci1&Yp?v zMT|B#MD}ZqgHhfUzjx1DtX$+*)^YxAz5L9Rv6OIz@x z=bt~hWgjN>>tWBYZZCPCh^2@>py6aC{{R7cLbNx+?EUlEeXwKyw###I6yy=ct2_xy zPj3ugG<-2S9n!858>YWs-7tA%Rz8CJL3uD#bCUEjE+GQdi<8(7i|mCn5B8N0`R^7| zj{l^Kb?*as=q_B6%7<$`ZW{0hIUj`858pnPX2A=E4GBYcAR48h|6w~Gf-#S5>+}3g z#Sas;lmNXFFY-FoGwnR=aD_lS8qVo$hX)AIBL(ju&_JM*03AZrm=P~Z)g(ssl+&ay zPw4y|bRVD~O_=CU27y+8+p2`QgPy8lcsOtKW~$pNM3N#%M0{jH?*Y#oalXVC9yQ9^ zR;k_>D?FU+E>ZKKDZ8kl)=w?U6>5ahyw&|9 zw#bc8e`xCXTM2nHfu{&OPvAoWza}8|F@+a@k-B6y9}7<$sy45rK7n`P<3m?!R;tI1 z=;4{{Q2fcme`Kr<_8#qKQ^9_0G`6Q68_(hpIaY0+MLGx^g8j!{Q|Ao-Y1Fdk;_ znm$1n+Ru{Tg{fk{e?udRX^Oo{HyCatWhMe*X;RgUPa*19B;dZ#;U|&dk3_^)DlJu) z;Ganyf#;AMeH{7wrFShZ*zx3_FP}~tMFIYA;-V)gIvT~gLcvz`3C&4(3$ay_c0+7c zFihdJ23~sVny%lGfG*tI9a$%yrV@HQzO?(Jd=?RmK3_)=FJwMpIng-2;qf;~rMLoK zNiA`Sn~VHvtOp1uB3uM_oVw3I|I-t@D1y0&Zm}~~lE^^<;sz-8h}a%tXNXN9c7brQ zaJKNoIGhHil6m-{IFZyqO&F)*dvtFoFpV#Z7eAB3s_Fj{IJZnw%iQO>82&luTocP- iMdzxoV+dzt;EUEdh2gF0d dict: calculate_rms, estimate_sampling_rate, extract_core_features, - extract_frequency_hz, get_dominant_frequency, get_middle_segment, ) @@ -143,9 +142,7 @@ def build_prediction_input(csv_path: Path, for_tmd: bool) -> dict: has_truth = True sampling_rate = estimate_sampling_rate(np.asarray(time_middle, dtype=np.float64)) - frequency_hz = extract_frequency_hz(csv_path.name) - if frequency_hz is None: - frequency_hz = get_dominant_frequency(np.asarray(x_middle, dtype=np.float32), sampling_rate) + frequency_hz = get_dominant_frequency(np.asarray(x_middle, dtype=np.float32), sampling_rate) features = extract_core_features(np.asarray(x_middle, dtype=np.float32), sampling_rate, config, known_frequency_hz=frequency_hz) x_rms = float(calculate_rms(np.asarray(x_middle, dtype=np.float32))) @@ -221,12 +218,20 @@ def save_prediction_overview_figure( def run_task1(csv_path: Path) -> tuple[dict, Path]: from scripts.config import CORE_FEATURE_NAMES + from scripts.inference_utils import predict_task1_tr with CHECKPOINT_DEFAULT.open("rb") as handle: payload = pickle.load(handle) model = payload["model"] pred_input = build_prediction_input(csv_path=csv_path, for_tmd=False) - pred_tr = max(float(model.predict([pred_input["features"].tolist()])[0]), pred_input["config"].normalization_eps) + pred_tr, tr_source = predict_task1_tr( + model, + feature_vector=pred_input["features"], + frequency_hz=float(pred_input["frequency_hz"]), + normalization_eps=float(pred_input["config"].normalization_eps), + project_root=PROJECT_ROOT, + prefer_curve=True, + ) pred_rms = pred_tr * float(pred_input["x_rms"]) true_rms = pred_input["true_y_rms"] rel_err_pct = None if true_rms is None else abs(pred_rms - true_rms) / max(abs(true_rms), 1e-12) * 100.0 @@ -250,6 +255,7 @@ def run_task1(csv_path: Path) -> tuple[dict, Path]: "x_rms": float(pred_input["x_rms"]), "true_y_rms": None if true_rms is None else float(true_rms), "pred_tr": float(pred_tr), + "pred_tr_source": str(tr_source), "pred_y_rms": float(pred_rms), "relative_error_percent": None if rel_err_pct is None else float(rel_err_pct), "truth_available": bool(pred_input["has_truth"]), @@ -378,13 +384,21 @@ def run_task3(csv_path: Path) -> tuple[dict, Path]: def run_task4_tmd(csv_path: Path) -> tuple[dict, Path]: from scripts_4.adapter import load_adapter + from scripts.inference_utils import predict_task1_tr with CHECKPOINT_DEFAULT.open("rb") as handle: payload = pickle.load(handle) model = payload["model"] adapter, adapter_extra = load_adapter(ADAPTER_DEFAULT) pred_input = build_prediction_input(csv_path=csv_path, for_tmd=True) - pred_tr = max(float(model.predict([pred_input["features"].tolist()])[0]), pred_input["config"].normalization_eps) + pred_tr, tr_source = predict_task1_tr( + model, + feature_vector=pred_input["features"], + frequency_hz=float(pred_input["frequency_hz"]), + normalization_eps=float(pred_input["config"].normalization_eps), + project_root=PROJECT_ROOT, + prefer_curve=True, + ) pred_base = pred_tr * float(pred_input["x_rms"]) pred_tmd = adapter.predict(pred_base, float(pred_input["frequency_hz"])) true_rms = pred_input["true_y_rms"] @@ -409,6 +423,7 @@ def run_task4_tmd(csv_path: Path) -> tuple[dict, Path]: "x_rms": float(pred_input["x_rms"]), "true_y_rms": None if true_rms is None else float(true_rms), "pred_base_y_rms": float(pred_base), + "pred_tr_source": str(tr_source), "pred_tmd_y_rms": float(pred_tmd), "adapter_scale_k": float(adapter.scale_at(float(pred_input["frequency_hz"]))), "relative_error_percent": None if rel_err_pct is None else float(rel_err_pct), @@ -424,6 +439,7 @@ def build_key_metrics(task_label: str, result: dict) -> str: lines = [ f"预测RMS: {result.get('pred_y_rms')}", f"预测TR: {result.get('pred_tr')}", + f"TR来源: {result.get('pred_tr_source')}", f"激振频率(Hz): {result.get('frequency_hz')}", f"输入RMS(x): {result.get('x_rms')}", f"是否有真值: {'是' if result.get('truth_available') else '否'}", @@ -435,9 +451,15 @@ def build_key_metrics(task_label: str, result: dict) -> str: if task_label == TASK2_LABEL: lines = [ f"自振频率f_n(Hz): {result.get('natural_frequency_hz')}", - f"阻尼比zeta: {result.get('damping_ratio')}", - f"精确阻尼比: {result.get('damping_ratio_exact')}", - f"包络R²: {result.get('envelope_r2')}", + f"局部阻尼比zeta: {result.get('damping_ratio')}", + f"全局阻尼比zeta: {result.get('global_damping_ratio')}", + f"最高峰后衰减阻尼比zeta: {result.get('max_tail_damping_ratio')}", + f"局部精确阻尼比: {result.get('damping_ratio_exact')}", + f"全局精确阻尼比: {result.get('global_damping_ratio_exact')}", + f"最高峰后衰减精确阻尼比: {result.get('max_tail_damping_ratio_exact')}", + f"局部包络R²: {result.get('envelope_r2')}", + f"全局包络R²: {result.get('global_envelope_r2')}", + f"最高峰后衰减包络R²: {result.get('max_tail_envelope_r2')}", f"采样率(Hz): {result.get('sampling_rate_hz')}", ] return "\n".join(lines) @@ -462,6 +484,7 @@ def build_key_metrics(task_label: str, result: dict) -> str: lines = [ f"挑战任务预测RMS: {result.get('pred_tmd_y_rms')}", f"基础模型RMS: {result.get('pred_base_y_rms')}", + f"TR来源: {result.get('pred_tr_source')}", f"适配系数k(f): {result.get('adapter_scale_k')}", f"激振频率(Hz): {result.get('frequency_hz')}", f"是否有真值: {'是' if result.get('truth_available') else '否'}", diff --git a/scripts/dataset.py b/scripts/dataset.py index 31b2567..f74a676 100644 --- a/scripts/dataset.py +++ b/scripts/dataset.py @@ -221,16 +221,14 @@ class LoadReport: def list_harmonic_files(config: DataConfig) -> list[Path]: files = sorted(config.data_root.rglob(config.harmonic_pattern)) - return [file_path for file_path in files if extract_frequency_hz(file_path.name) is not None] + return files def build_record_from_file(file_path: Path, config: DataConfig) -> SampleRecord: time_values, x_values, y_values, interpolation_count = load_aligned_signals(file_path, config) time_middle, x_middle, y_middle = get_middle_segment(time_values, x_values, y_values, config) sampling_rate = estimate_sampling_rate(time_middle) - frequency_hz = extract_frequency_hz(file_path.name) - if frequency_hz is None: - frequency_hz = get_dominant_frequency(x_middle, sampling_rate) + frequency_hz = get_dominant_frequency(x_middle, sampling_rate) x_rms = calculate_rms(x_middle) y_rms = calculate_rms(y_middle) diff --git a/scripts/inference_utils.py b/scripts/inference_utils.py new file mode 100644 index 0000000..ac59cb4 --- /dev/null +++ b/scripts/inference_utils.py @@ -0,0 +1,72 @@ +from __future__ import annotations + +from functools import lru_cache +from pathlib import Path + +import numpy as np +import pandas as pd + + +def _task1_curve_candidates(project_root: Path) -> tuple[Path, ...]: + base_dir = project_root / "evaluation_outputs" / "task1_final" + return ( + base_dir / "task1_final_dense_curve.csv", + base_dir / "evaluation_dense_curve.csv", + ) + + +@lru_cache(maxsize=8) +def load_task1_dense_curve(project_root_str: str) -> tuple[np.ndarray, np.ndarray, np.ndarray] | None: + project_root = Path(project_root_str) + for path in _task1_curve_candidates(project_root): + if not path.exists(): + continue + try: + curve_df = pd.read_csv(path) + except Exception: + continue + required = {"frequency_hz", "pred_tr", "x_rms_interp"} + if not required.issubset(set(curve_df.columns)): + continue + cleaned = curve_df[["frequency_hz", "pred_tr", "x_rms_interp"]].copy() + cleaned["frequency_hz"] = pd.to_numeric(cleaned["frequency_hz"], errors="coerce") + cleaned["pred_tr"] = pd.to_numeric(cleaned["pred_tr"], errors="coerce") + cleaned["x_rms_interp"] = pd.to_numeric(cleaned["x_rms_interp"], errors="coerce") + cleaned = cleaned.dropna() + if cleaned.empty: + continue + cleaned = cleaned.sort_values("frequency_hz").drop_duplicates(subset="frequency_hz", keep="first") + return ( + cleaned["frequency_hz"].to_numpy(dtype=np.float64), + cleaned["pred_tr"].to_numpy(dtype=np.float64), + cleaned["x_rms_interp"].to_numpy(dtype=np.float64), + ) + return None + + +def predict_task1_tr( + model, + *, + feature_vector: np.ndarray, + frequency_hz: float, + normalization_eps: float, + project_root: Path, + prefer_curve: bool = True, +) -> tuple[float, str]: + features = np.asarray(feature_vector, dtype=np.float64).reshape(1, -1) + tr_value = max(float(model.predict(features)[0]), normalization_eps) + source = "model_direct" + + if prefer_curve: + curve = load_task1_dense_curve(str(project_root.resolve())) + if curve is not None and features.shape[1] >= 3: + freq_grid, _tr_grid, x_ref_grid = curve + query_freq = float(np.clip(frequency_hz, freq_grid.min(), freq_grid.max())) + x_ref = float(np.interp(query_freq, freq_grid, x_ref_grid)) + x_obs = float(features[0, 2]) + ratio = x_obs / max(x_ref, normalization_eps) + if ratio > 1.15 or ratio < 0.85: + tr_value = max(tr_value / ratio, normalization_eps) + source = "model_amp_adapt" + + return tr_value, source diff --git a/scripts/predict_single.py b/scripts/predict_single.py index 028bfb2..146cd4f 100644 --- a/scripts/predict_single.py +++ b/scripts/predict_single.py @@ -10,9 +10,11 @@ import numpy as np try: from .config import CORE_FEATURE_NAMES, evaluation_dir, make_experiment_config from .dataset import build_record_from_file + from .inference_utils import predict_task1_tr except ImportError: from config import CORE_FEATURE_NAMES, evaluation_dir, make_experiment_config from dataset import build_record_from_file + from inference_utils import predict_task1_tr def parse_args() -> argparse.Namespace: @@ -87,7 +89,14 @@ def main() -> None: file_path = Path(args.file).resolve() record = build_record_from_file(file_path, config.data) - pred_tr = max(float(model.predict([record.features.tolist()])[0]), config.data.normalization_eps) + pred_tr, tr_source = predict_task1_tr( + model, + feature_vector=record.features, + frequency_hz=float(record.frequency_hz), + normalization_eps=float(config.data.normalization_eps), + project_root=config.data.project_root, + prefer_curve=True, + ) pred_rms = pred_tr * record.x_rms out_dir = evaluation_dir(config) @@ -99,6 +108,7 @@ def main() -> None: for feature_name, feature_value in zip(CORE_FEATURE_NAMES, record.features.tolist()): print(f"{feature_name}: {feature_value:.6f}") print(f"Predicted TR: {pred_tr:.6f}") + print(f"TR Source: {tr_source}") print(f"Predicted RMS: {pred_rms:.6f}") print(f"True RMS: {record.y_rms:.6f}") print(f"Relative RMS Error (%): {abs(pred_rms - record.y_rms) / max(abs(record.y_rms), 1e-12) * 100.0:.4f}") diff --git a/scripts_2/__pycache__/task3_identify.cpython-314.pyc b/scripts_2/__pycache__/task3_identify.cpython-314.pyc index 53c676f726a99d2260bf95045668d85c5a203ab9..876713108da144d00c0d6068765ad0746de85184 100644 GIT binary patch literal 23732 zcmdUX3vgT4ncf9)@diQSN$?H635gGhq$r7c(0WrZin^pgTcQMm#3cm_0^keKGGWK2 zGuw@uq#MO{D_X0r*wj0rJKhPCaT{gEooIK{iFPKPTn0eM+%Q{r<2Lnn+CsHGZqv-P z|9>v-#RUyYv^U#M?*QkX$A2F8oOA#4o&W#O-EJv0QgHpR?%##}@GwRF8Uxa$6FsVr z-d9o7JavU)R41t+RYFyd{*sd=i4wT0PpXIL1f7f1Bs3&Vo6ug?j##9k^7jUcq5Ggb zymEQ;LxzMg=QkzHmrGeaMb%K3OBv0cCn+j;f`l|9HDWHMJ5sL) z3SnC=+{D-+Re~uHi9{2DL?{}G4?~z92_&ux>O(U#0_8Ob*0FOVe&*Qd$hqOsWBx-= zo*WhI1yQ5NhDXmaLfK>>5{w1nasPBEoM1V@7Gx9bRKg!*rvh{S7#p}27pj8M>!C;> zlJHM+?6b3MWNOYIiuf_BU}YomS&sF`0vrnj!4N07*mxpzJ&<7i@xb+1I24)ja}W}2 z92*Gw%M_nn|^9L~J zAg+cwuez#YR7^>?s)_QdKSWwV%dv@BE|RnbxEW|!ob4J8Txa94z!dv7_MUb+%0Ua` zSV5yO)O{h&a`EmdZZgu>J;FtwfhNVf56_0eL2Uh4Abw2_b;qZ;P%IJm_jD&9zSkcL z!h?sV=elBZg4yrSjSK7db1o>o4*qcmIP=t}RF$g!Av(OHmE1pqcxnnA>4JAkeA2&x zrp;3mQV&oQ_*~?YX84Fkipk%l88f2I$4S1z_a#P6U|r-IA>2g@iFe=t(5Jg|^~Y2` z)ktX}1%?>)gv>oA%A|B3aE}g?6ap$uK}cH=(kVmIvoU(cu$vZ3sHZlm+qg;HCS{3? zSs9YcFWo}XX@ggm^i+-gmjS|dP=ZDzC$ulQ&lgN{} z%6g_lTs+~41YimW!XEjZ^T7M&3BoG{bJ7EER458*(d$USc9QoMQWJqLb; zG^Eam=g6qf6NmQ)ee~sT-kS;_7dD$(0&RK0Oo&Pv;Xn?L$egCpt}@ z=J0F+*x(H)>#7H{b|Q0x5kcnY32@O_7y}%m+*==l8b32oUD)OEtI->RF%E5J{RwWC z#b=1e!ht!GeTI#&$P}T|LNmZWp+h}^U@(MFgPk4^JcpwkN_2Tk1Z`jn<2esZ6|a%Q z_d?M7kz)gw5p*0IkA|OvUKV+9T+mBgPp}-~W@fJge;Pq1$O&es!w;z>PEf}aoS>P8 zDUlG2giT^v9BMSjL?-CxqS1t)C23q7S&+xWIpM~2Ne*Ff{}n>>)VjUmj{fz7t9I{- z-J7y+U$ytH*n5``W$gX)XVy!dtEJ69D{ba?zWgyoSq_!VA77)5t8~Q*U9qMwUoS0R zFSp-3v$!wg=vuDKI0o)r%s5W)r@+0KD!;T=dOmBQH13ZnjYgX_QyP<$)g@)+n|5a$ zyZM1r8ONzq`RTRNW1Hk~u9w=@%I#~V(M<@-*UG9fyN*;yLxo1lY+E(ete9$2raGRk z<0{~h3l=FH7I^VhA9+Y;6kntmHA8p9V(!&&Rv;4$WDsNUs*z#kY#2$y+sSgsIJ_zj z*A8)=5<4Cj6|%DI!vsV?xU$eH<3n_of?2#RCrBvYNFXB@t_!4*cfd1GF|d}7h$#hb zIdwrbTe5?iEjbP?de%CU<7i4EiHlt;ISp+zz)Ck;5>xF3eip36b$LQwpP-H+>ya5& zn4S`eBf}aF-7;HJPf7{FO3u|nZOI8wDDDBIBXOL7XLEkzIUePv*kEU#RV4zUFfa;u zzkqE_hJwIl2?Ig);_2-g_R@j|CU+t(Xu{DO@PZh9z>C;p92eyTIs%x%Tf+4~I|Q1L z;TzB@z5zlVu&odXvr7>GKQ{R(^;4Iet+9j*{^3i!rH-#l^5z#lQkRq&vs8(}_)tyhkE;IX z&nl0^u}^Qh4;i!tQ%^<>?WC$rr(Ov|GbN0M(RR}?1$E*SGyzzH=Hqk6e`}1yx+1{< z#|7qqB>?*KV=6G@=}^K?&|G2`))q2PleLN2SeU((&(h(^`!7Fs0w0FP#&?hwe@f|9 z4G+;-vsP!-aDTtNDyxBeR!ceSvpV$YDQ(k31Ng*gt(=l_2P}ylpei=00W`v77AL1q z0jNgs&E0S{z>6VG&0d(8dH7XqDmPH$s`17vlyXDjD`85MVMet(_d>+f<7zAura1Bd zV7GLbolXcEE;MsB@!|i53+5#q4<%WlBoY&hNtTPo{o&9xR?vnbun@((G`AnJatF~F z1V_-bk>}WOG{(m9X;0na=I+)Z!$%FD03ykJf1 zs@BS??i{{*_V(GO?Wvl6zO0|u_pg^(zcKgKxrOTGF5WkuIq@X_bcBySn|Yc`m&NCY z)=VzmT+;DZsYYs&fZ40d| zO;HxEsv9oY0ez0zu9S2r#>uBH40kBXmY%F=`!wTyK25J0x%1m9K@&vyJR@aJo;sB9 zgxLUW4iYz_o>(*vc=0*bgD;CG|FV$v4maw6nE_t{UY0fuxs5qs$*$o64pwPucvvd^ z3GN7ly=Lwh+yq@9E-nH(@`m;I5tw;35M#Ncn63^)=D6L&7mm(wo#>;lvw?_E5(@$E zdX`IgwcKfpHs)IWK))DUO&C=thnWS&9`uBF}_r?aBk_D zd$vD|WNP-OoCoHI*DLGpaCdLqzOndHs=hx}xtph5>yApks$*IG?&SB`@3Q~Wl&RXE zavYdHyIxhlrmwog+#SC?zBu*H;qM%K``DjeTMlGe0LPwr5XcOi%e0(J)sAHJBLE1h zwLDz`OQN}s?>fYr4)OG%Jg<^CE^-0FadFmdCM|*Y7tsJ<0Y-xL>Xm3hn>kuMcxGetPRPfg&?MlP?CT60yUXiZTr zKVDSA8lxu+(qjO!FH$oA9|a$l>WgW{M~YFCRswwuHLhXQD(b9ijEZUVe0oBbpeWi~ zg!eZnd$;KG_~dm}YLj#aX?Jctj_G`gHGA@S(#g^ z17(*n{kU7Xa>fkfLj_@@&rqzla-&N6YfKe0jL|Wq#NFtIOVL}(xu2tw+UF}#u0t`C zBA0tZfexE-5?6 zJq<(zl!$vK1E?XModxMYkb4SKaIN7kp>qTrubn_4?i@zoev&(n&J;RF(GlTKU@|VK zuLqtNmq6`xkcI*Z(gH{W=*Eo#p`e50(Dm8tAfbzf=U^udlDFCG+!z*u^H)%V^wR6% z{21{xIstSh(HRFvFrpm8A5Va!LNG<5p*W~k!mtcW0l;V=uQ1KV+3DG^|2haa1WPU; zt=HCESOLbwjlZaF5p;;4AW_i7aTh5V$VPLDjU@yf?)@QI(1Mx>gpK5$h$gse*ucV` zizASR-KPl7JjzQ%YPGU*NF3grqDew*uRwOP2;g)4Hh=YCvTs8Yx3LE-!lxp3Zsu@VR_T3vx z*$&U2UbEHQi7)mpRcD&|mV@tIxIggYgMWT7b!0SkU@X%=ma4z7V!Hrj&W4P$Ykqjm z-I#HA^R#nK=UmlQt>~&2o=WRl*D70HwXD^&^7_j6ZJt$I+lsC2os#bu-!`Ufz021# zwjsw?U+BdX0c~1 z4T~o;mQF}7wY|}EtN*qBh25Fzj%916swZXdO_lb|A6F#3#e94GTjNXPsg_-t>RtC5 zGFAIh_Wh~S10cKwY3(#TX8+QzRQ=A3YiFi>=lsd_0}MYJOdViX=mEgJ3#EM90KhXt z#e!u`U&&W@tm)l%8t=B=Ze6VSSKM(fd?{7g zxu&nZGt0N^N!9Pomg&vfETz}v;2eV8zcb+FuHX!QoS^j&fJla*od`j}ArKVq1cJgn z2SL#%LQwRH5EOimfS|wozX3sICHwCcg62#8e+WSp003eAh6@%o*^oh6gaC*-nGzL% z%<*I7uJdUh0kt*)Xth)h0O-d7VOj}dq270ED$M^O!D%lql5La4C9&-jD=AsRuRyVT&SHQJ||xula1LFz-qkz0{Xd9ROTq z(vjkD^}iJY;aGXh>MIBMtJqNatOa1$mLCa4=^nWv*@FPjXH!DQ$G~Cx=P6lHccWaT z>`{Wt&(Q)0fLD8d7Wy3dKHU@!J11JCYJH_r&l5Q8idDqi0ENwr2}VZ=;9I2}`8%{& z34G;uFAruJOTK@6E~zG_Oc7fl$CfK%-BK)p!!{8PA9@@d1`yeaBA8ADgT+wsaB@Q# zfoW7)qqaH+Kt)L*N(476gJYowxM_4yFvwj+Cxp&3=sb-MO7#eUE(QQa=3gr+G;{W>q#Szv!K}SdV{+)_uAZIRjQ&RW$y$?ZLeP} zxVkOl?%NDkZ~wjF>R0yyvc7q(D7Ic4|Nk9bx5MlMVS?%-M5MUR`taHYOQncdrAB0B}963WCg_Yygq`0lO zhoe)0u)jbRA=q;}*n&D-%E){Tf-@Iis7DtZ(k>WCbB2_m z-MnBUkp*ofqWH&dbl!zti2p}2%AZn)RH}-H#9b@1XVn;dzub}4pa(WPAbHmjubu*p zy@7a*l*^Mf5wDrD*JewJ*Fx2_X3L1ToT~R`t;B0XwZ5Hr9aOa^>m*(mRaKv@AYM0B z)ts#)-YTlLEn7{zHCR|J@z!CsdgAp^?#65b@itP`E!igGZN?H>h_{uhtjo3$Z#$(m zKJ=nj+>|MyE5dm))!zdJlhXRqCriy^s`1tWcxs=39OC6G@!5TPr5qFGZxb@;9D!0d zUJcY0^%WJ0a)?5M&u1fJ|CklE_F9kd21T5Jgg`*M<`~!Io>2`fF#p| zBKc&iPDKm!dCN$!;xSCfx?Bhudo&WfDN>pYqZFl?9zkD3t}U_>Nq=6Xjv`xtzFc*L zxeYBX2v4 zJj3IAwXp45>(#Pih$y`ldKbCycu zkiYvJJ}21b*D`fX;aF}drfr+9^bAU*O(jLw8GMgJ%JlG9ZfI;Dwt-b zC12KO@)mmrH*8j#pdGDDTR}Uz3))fnnX-ztqtaK-w3BvtKX*GSK^xHy?bya_SGMCb z^%k4Q9iV^cgn8WMvygkY&xZFNStD0uCe-G~QSvc63V20fzd!yhD5T%oTTqzik#B*5 z(-+;`=ef{$`G`k)ru-ea=gyu&eO#5I%;MYXsw?j@Jxs5{=DZ5v)2G%ai470owGqlV{WuUe$kWECjyfLGhhoyc7MroCU66c zlVcuG0>t3l2@3|^nyF|QbPonN=D`v(Yp|pYOtRsmWef&%L7SAb+n)AEJO@0#Rd@AGLvx0bGynD{czrOKER#t^iwCf; z!=ku+qu%P7w<9!z1I&ZbOefcUDz8{^hp(AfoxUXaIAEAT$5AKi9`35>S(fKAiuc7lg zI&Yx!kJ0%j=-fi*AUHy)G>pVlmm9)3HB>Lqf#*S4XU**}$ZkS3<`tzMiU^Z<2muap zC9cMC*+8kqhD%cAiE;}oK4WgI#LV1%$oDS%<6XddKy|l`vNbIZq^w<_O0(8vti7OT ztLs_bci;G6I8}FPg|41IdvnhR<@Q_V*USs&7gbAp?j1{)AAaTJ{E=53YqahwXI?qO z+xOpl>W5c;aE0H0hPR$g)5Gu6`d6+0-~=Y#tX?+WpW@F?&7Vor!8I{qm3FSs&IL=l zZGdmxcds`+IKm%1&!1;_*Jzp^lc`)QwCm!)e-4SmnhB zy}a)kB)cY&oqOfn&C?l2=hDlW?lb(@(|mB6kGwR0E=|8I*N3g*-Ccb5k^4vZuAv8N z-UUURmDB%j?YYQL#<58WWs^?vo;~*s4;X$thHZIPO2e+JUDPj|?#Fp1jM=UiCb<_4 z%hk}r@k{*FmoWXMHFs0Sy$3}px$cKY<{d43+kyKN{DpDeKgD00!%C780S=zEk#F96 zzlMKmg1;KhQq-*KRLT6Ak10xXx@67Oka6{a{3u`G&2xOq-uu=32*Y2z!cV_|<-I7S zW25W&#$ER=@nhrs)HDaWQCV7FF4_l{Rn3H z7?fXhKYFKw@rCVJ!13;d>V`Zzy;$yaPd6YWdO(I?cV~6PtEWt5Sp)Hc%GCSNgx=5UO!2{hx2Qh)9`Nge?YXSQdDKA_ z(0$5Ul)QC+JvD)hpIl(JjVl4UCU8w57Z~Cpn1YL#CLNoQL3&=d2nK61?Dc7UTHMZG zDQDEucJBsk?YZZ<+X~Yjpz`T_C0A&MmZ-q|N)|_q|0~$^Xr!9Kgi40DF}6gAm`?SoFaRhW?VS&UxVM#n5$D8*1vO1{<_v3HcNqU5qIlx+M$C6`asOLgW= zmu|pV#jG(~%s$Z|r<#0b=rfbtXAXtQ)W$T0{)p-8hU4EpYEs3N`c&Xrd^EUaJ`K3# zJ`=b)pBY@M&jPN^R|c-pR}L<0S;2MqY~VT{+XIfrj}VmO+Cn*BXpErjiY=7=g+>X= zc2Bg)BjI2wm18Y!=YXo}u`O_Jsh6s^P|g?arJ602{Y87JR?$n&&3cIj!Jtnqb!QOX zK{#mFF7@unHp!R2D?cNkpzO#~NIp35P$!K?J%uxWqFYMK^Gg`1`aG9}kSEV6A=IFF z3VCx6@peVLPO1gs8}sqwU)^M8G%09qZ$i^NA&WaEcFL`mO)wQTwiKnbNEFf?+Scl` zf?G-6;WnQg{O!f^K6XxQ;q8`cmF`fNH~-Y2_bw1RE2y^iZfz?lt8J4$++IM{u?bb@ z7wUs8yt~p{<9$>l)Xo~U{W=8={(Whs27SI%ZcBV%`Q_2*l=~cWlvQ0p@ z?uZJ7GO&d;<-UJ06%De2nVh=xg9v5zdSpg4T>Zh|j5vKAPr$}M0DeJ>=PIry;SQQ% zxZCkSAZ*s)yf}QFVj2#_L^+ApFa-x|SPoI%UtnHMl!K%AaEg_bWQL=vxes_;O~T<~ zjO3n1=PWw1W*cRA4JK6ju_R;%SwCXIjiG`y?;V}1ONW~3nC@z#2oGiiT zEt{TBh?=}fu_J79#C16AAHpN4a84*E3=%1);G}*+uoo%A6%I%GdPzjZ} zv?3>zq-+H$Ht}SVV9RBbH|1Opu?-MvaPWv`u2^l=a%z$KTykPTjfx@*;&Ux=$hj14 z5lo5z_d^^LM4Vtp0RWZaaA6?_mpQZ>%Zo)gguI*y9n?s1h>bY}yb3#yDGZLFGl9-4 z=)8)~eRMj}`7?C>5*!d4hNJL70$hazb(oEKD-^O4oVbKcAvGbB5N8;CumI+8TrkQ* z4}K!c0OJn!`TV>UOe9vcZ{8?tc?zo_Q<6*si#*um;8@ODM6$oYZi=I{7)2tI97J&! zMIx>Me1~E#x$!tS_n#mm*%kSqi-}3I)p_%&*UMJR+g8fkQsrJSnxz{T_onI259)R- zSEh$9rJs2rUH2kNNKUUiwk_>TId;QQv4OEQy=~3av{q)rUrpV|Chb1$3JvGW+>2Eo zTd0cG9D%JBOpDhnWk67&oj(Ujj@HE}ADB(spMx(k*tRVV@=vm9+w}bD4{GY}p1OT% z@rC7snYP0}X-d}&^K|8UeZ$@4?c~ySzW+pK+faIhSsl5!GIB9pe~G7S@bnv)8!wvp zj{OXV03u4IAWY`dCN1wDZFs8!5+bWu!w3 zVgb+@H~)9aVb}h`24}$G{_jDDtbzMOI0Np<>d>dB9CcX(`ixX%d)9 z)%Rt~&{s}14rHz9vr+EWtQ~z0s;(#NM4yYQ-IlEY-^XsMww^ryOUsc=+p(XFW@^r7 zD>0>t@-&h@=;HUD%4|EG9=)(SdU<8^a;E-Dwi+pFDCnhUZa=eB$L|@+cu%I!kFK7- zv~vDZrfwozixhQKUDM*hf`0JyUfy+kh#J zRGsJU#oHH`%>1t7nf4Pu;eML@t7NA3VzvpBo290{$oCz^wlY7x_*WM*HBV++Fu9eg zX;`dTie?%PrmF|DZ4ftq`qOQc-n_0ieMA3M{nw2sd^q{#ldrNM)5_A7W-Xduk1kZc zekH4h(EEC4x?)>a176si8_Lq<4Otxq^pv6c&auVLRLzd80V9l*q5drm-#U=;?8}-k z#!MNi?(AB0W*T}@)xFtLjI&UN+B@vxsZ8@gs%~Gl4CBfvL*1SDTfK`frkeUv^}DiG zjJ4%r`L_M3rUR+^gIPO9JD>o*Zb#OM0T*R(FVto$(C-F%zP=}0iGeE0(DK&QJ5BFQ ze&?Bg@yxxt2bRp<@w9(3-8Gfk7EHCW*=i)Kp$yG$U06EE4-99vpG&ojWNR_Hjy%d@ zWZ9o-J(ltu&(>qC2g_uahWULXnU3?RRwmnk(T(C04KMG@G!Lfgp2#*~Y%^tOdq@4w zk#|eJZ~UI|-hODnKGJ|4Q>l(%s-0cax4bpGG=VjDo=dfjtm!>(wJe#T%JyTa#^Y=H z)^{{ZVXUj`e5##U)7QPFSqwuZt%E7g6Ki_s8?y_S;4ne$u2e-o2IC9+QWb3vcVg{S zWigX@r3!R$h4u~#O0bV8Je3Z5H4VgR zm4Fs<%ZA|v_2__7ucTxQptm$C_3w(hOh9cepe}s`wFRilHc~6CHWVe71Fdz7v^Jo% zZ;{pkw9YNkx`4J~i?nW_Er8z)^u#95Rt2=xTP(K*Xlu7fTL-lDTcq^>ZG)0_)A`g0 zlucWrYzE2#8<0)QYz4}KjmV~y?Lg^O*8Yfo-3GJ;)(MZG?Eu=&EtcB_v<21%g{@_F zz~*LWa_4~q2Rx#3dB8IY-@KXuvC1ekQ9 z-dIvFd&J!mq~QjEyNeDFPO@y2Siy-OLc_c(nDiDpe~Rf9!=L$j8fbrk4k08cCi|z8 zo-=6E=Rv#;TFV(Ydzkwy0rt=$s~d4U|^WU{&# zIr(S-P$1S!Hk=T(w2!EG73sV=x-HcFI?~q^rI)H+Lz1SmlHT>K91WAa}&3bbpO>*5l%r z-%yhZ-;|3&Ex6IZk9lfh9B5u(fBv@+K|cP~D(*4hyKImEhcNMDciex0RI>P1l;9G& zP8)dRfi!&(Jg=Tva4niLl^wjhlQ(vy>26FlH!QX;)n*!Z^9_6MsqY=X-~X=x;+*D< zsENc}=N69hO}jER{mXH_8kDEJaWG9kL9#Y29_P0o$}}Io-^(|h0L`KCbecXxk~)@r zeBV%}^W=jjzGIj-g4%9`B=jx&?w!l*I`d$Z?>o;MnKV60;@g(AeD_g6KcKyBJEbV8 zb7`92g)oR`F7TZ`-Z-A7pCnm!El=|YM&a-zf8j}f*HgUlVw%20QreeX{Ejmj@7V`a zd^^J%QOipbx|g;5?y<~{3%u_szWX9?yp*OVNXqVeT7Gamv*$_FOzytS8?U5kKS|lP z)X(>fWVVwAY*RL%X2A!_%uMY7U$c)l?oZPP@|B#)>^RGJE1!JFvWwp{mf3lM-=Qp| zXStsrcrw%b6yKw457fHcnc4oteLKJH7;ikDrcaQCmPRlXV|H?bj*CSRy(BM4A5InHCBz8ncC>>ypdl&ZU8cC?bv zuz(J0*l@{oD9j)iSksy%fBr7h<<)xBNO?Y8@=0`{B&5Z-m@2Qt%i9gZH+)sFN|u0} zE-y-x1s%hnw(WxKL|7o&yVJNRiV`=Jod^b5gUkI0DgP23!nTWHrLI>t+{W`z1|I9Zd{=6mc3b)vbV20s}}n1 z9(?m)%DHX5)Oxe|miL?9ulYYVXe`?Kr*Lr`IWIW)mCD4`Fc_U-cfmpUIU>qQphXVE3Til~E?D8)ytxSY zsxC3s(TaB>LJ~BQ+3T@6B6|`{evpJtT?@?Mw=m%=eee}zZ6Fo{iyb0e;##q0BG+oi z3SIK&d%HwaFaK-;z7a@_e?)m#Ucg1GIy=qq5E@y#P|(hE=iv1QfW zwc_qdxp%C(_pG@0+_PrfAWK5t^(k#smEbqfNdzHJMgd{({#3;c@sXWAEwsP5=B}nW&d-*gmhna zB>&F;b^rH&&iy}p@3+i1Q<^=RY9)fz^xc1w-+HfWYM6VUG^CgT3#+vvR%b;iU5d%p zv3lSPl*{^9BhV(1u0OepHSZrqC;?6!=w!VYwNP)c+=!L<#9zRAS1tKOt*ilQ1wzcf z`^bvIMeYQ5%Oqj74+V7ONwCRX)J`W;SswF=$xL!49Zwb2^L+eAlgV@NnH)e4Jw+kL`iMe7Yb;g3{!nKsRl;7XLIy`kMkzVc*OnBNQX-3|q@t7> z8J(qSN@B&6c9)ElGND>e$xJB=GTBO2pc+aws7m?B1|+sZ6jLgI z@=^dH9`afx$uw)%xk^RI#?&` zVrBoUv)vCM(7IU<>t*ZOs$M)x^36W8k5jVs)u>t0!tBuev)QNJwNnD&Y{n@e5GRCm zAu^};H%ngGal`XE@?FGh)lJ-E8`#}#TlWc3>L8PXq@^XQnAv(KE+RE`GhCIPsJBK8GdG# zPiH+BJc;=1d@?=b!SQTz&a=d4J51sJ?^$qpSp**v^|I?mzUDuT(x$e7<_xCWQ zqutiOqFT)eY$|tbT{>Hk%73hX_l>K^t{l7O%qv!EKd%r*#)YYiRC{Sl`?b{krPM6c z^sZhMj*kkxV}j+xx^x`0WahxvOh4)DnC-g+ab{izpk><2FjbF46bqHIVx>|7O4V3X zl8K}osk9}9NGefP;E}4AA**ha$_0MJy6UUEYdk2EFc#vnxzwDNRk3aca+Fn*pX>2H^jW)M_inNt&Um9QM-c#^PTRA&#<|-Ptj1v~xC1XL--{m5m zT`E}4=S;?<@fvf6T}&^~psrm~ttWjU*l$;(pK5EFYh<=LGFsM@wR?KqZfTWW?& z`dgPv!&+XbM}9_z4x2QrHC9{pAB&KtZo7u9c}Xq)#-IuSxi;IdWs4pH*o?gOBl&_q_ z2_6R|Sk$VgF@NSPpPpI-g8(B6#exbj#mYo>ah@+~#i=1W&z?>28L_~k?brt1qS3A1 zfDQ#z->=Z&|EEN0NEPu^YV+aCl$FNMCo@IeOXs2C)0!sFd*5V>hFLzIo=nZnOhOmI zJnbNurH>l7x(L04Bl#T6(YeZ3-_s6TMwL zM3hH-z1vNys1mPo0VeR=!uL8>QS;JqGWa@mKo5m)P*@<#j(T$ENR8@4U^iQVgnV)2 z?{2g81MUVVI9gVQS566zSA~wZgySO{9itnL(KYiJ`IYfh_^})r>kEe3g3(zpxhd8( z6wKa&D^l>jTCnaaIGYL<-(#!VSe5I4Y9}4XY?rH$!u0;|O&ggg$j$54U_q`~*Le$a z`>mdQ{Qd7gCVqUddWUOCWh@Y`dOCh4wL7@ zbu#*OmduCA%3J0x`jM7&hxzR0bfX_S)9%4PrDEsy&=bX;={f8v=Oox%qc$}?S$6f64=g$eDu8oH7 z)o#JxCn%4sOOHN~Dg tuple[float, float, float]: + amplitudes = np.asarray(peak_amplitudes, dtype=np.float64).reshape(-1) + if amplitudes.size < 2 or np.any(amplitudes <= 0.0): + raise ValueError("At least two positive peak amplitudes are required for damping estimation.") + log_decrements = np.log(amplitudes[:-1] / amplitudes[1:]) + mean_log_decrement = float(np.mean(log_decrements)) + damping_ratio = float(mean_log_decrement / (2.0 * np.pi)) + damping_ratio_exact = float( + mean_log_decrement / np.sqrt((2.0 * np.pi) ** 2 + mean_log_decrement**2) + ) + return mean_log_decrement, damping_ratio, damping_ratio_exact + + +def estimate_global_decay_metrics( + filtered_signal: np.ndarray, + peak_indices: np.ndarray, + sampling_rate: float, + min_peaks: int = 5, +) -> dict[str, np.ndarray | float]: + filtered_signal = np.asarray(filtered_signal, dtype=np.float64).reshape(-1) + peak_indices = np.asarray(peak_indices, dtype=int).reshape(-1) + if peak_indices.size < min_peaks: + raise ValueError("Too few peaks were detected for global damping estimation.") + + amplitudes = np.abs(filtered_signal[peak_indices]) + tail_start = int(filtered_signal.size * 0.85) + noise_slice = filtered_signal[tail_start:] if tail_start < filtered_signal.size else filtered_signal + noise_floor = max(float(np.median(np.abs(noise_slice))), 1e-8) + useful_mask = amplitudes > max(3.0 * noise_floor, 0.05 * float(amplitudes.max())) + useful_indices = peak_indices[useful_mask] + useful_amplitudes = amplitudes[useful_mask] + if useful_indices.size < min_peaks: + useful_indices = peak_indices + useful_amplitudes = amplitudes + + time_window = useful_indices.astype(np.float64) / sampling_rate + log_amp = np.log(np.maximum(useful_amplitudes, 1e-12)) + slope, intercept = np.polyfit(time_window, log_amp, 1) + fit_log = slope * time_window + intercept + ss_res = float(np.sum((log_amp - fit_log) ** 2)) + ss_tot = float(np.sum((log_amp - np.mean(log_amp)) ** 2)) + r_squared = 1.0 - ss_res / max(ss_tot, 1e-12) + mean_log_decrement, damping_ratio, damping_ratio_exact = log_decrement_metrics(useful_amplitudes) + + return { + "peak_indices": useful_indices, + "peak_amplitudes": useful_amplitudes, + "fit_amplitudes": np.exp(fit_log), + "r_squared": float(r_squared), + "mean_log_decrement": float(mean_log_decrement), + "damping_ratio": float(damping_ratio), + "damping_ratio_exact": float(damping_ratio_exact), + } + + +def estimate_tail_from_max_peak_metrics( + filtered_signal: np.ndarray, + peak_indices: np.ndarray, + sampling_rate: float, + min_peaks: int = 5, +) -> dict[str, np.ndarray | float]: + filtered_signal = np.asarray(filtered_signal, dtype=np.float64).reshape(-1) + peak_indices = np.asarray(peak_indices, dtype=int).reshape(-1) + if peak_indices.size < min_peaks: + raise ValueError("Too few peaks were detected for max-peak-tail damping estimation.") + + amplitudes = np.abs(filtered_signal[peak_indices]) + max_peak_pos = int(np.argmax(amplitudes)) + tail_indices = peak_indices[max_peak_pos:] + tail_amplitudes = amplitudes[max_peak_pos:] + if tail_indices.size < min_peaks: + raise ValueError("Too few peaks after the maximum peak for tail damping estimation.") + + tail_start = int(filtered_signal.size * 0.85) + noise_slice = filtered_signal[tail_start:] if tail_start < filtered_signal.size else filtered_signal + noise_floor = max(float(np.median(np.abs(noise_slice))), 1e-8) + useful_mask = tail_amplitudes > max(3.0 * noise_floor, 0.05 * float(tail_amplitudes.max())) + useful_indices = tail_indices[useful_mask] + useful_amplitudes = tail_amplitudes[useful_mask] + if useful_indices.size < min_peaks: + useful_indices = tail_indices + useful_amplitudes = tail_amplitudes + + time_window = useful_indices.astype(np.float64) / sampling_rate + log_amp = np.log(np.maximum(useful_amplitudes, 1e-12)) + slope, intercept = np.polyfit(time_window, log_amp, 1) + fit_log = slope * time_window + intercept + ss_res = float(np.sum((log_amp - fit_log) ** 2)) + ss_tot = float(np.sum((log_amp - np.mean(log_amp)) ** 2)) + r_squared = 1.0 - ss_res / max(ss_tot, 1e-12) + mean_log_decrement, damping_ratio, damping_ratio_exact = log_decrement_metrics(useful_amplitudes) + + return { + "peak_indices": useful_indices, + "peak_amplitudes": useful_amplitudes, + "fit_amplitudes": np.exp(fit_log), + "r_squared": float(r_squared), + "mean_log_decrement": float(mean_log_decrement), + "damping_ratio": float(damping_ratio), + "damping_ratio_exact": float(damping_ratio_exact), + } + + def save_task3_figure( file_path: Path, output_dir: Path, @@ -85,7 +188,15 @@ def save_task3_figure( selected_peak_amplitudes: np.ndarray, fit_amplitudes: np.ndarray, natural_frequency_hz: float, - damping_ratio: float, + local_damping_ratio: float, + global_peak_indices: np.ndarray, + global_peak_amplitudes: np.ndarray, + global_fit_amplitudes: np.ndarray, + global_damping_ratio: float, + tail_peak_indices: np.ndarray, + tail_peak_amplitudes: np.ndarray, + tail_fit_amplitudes: np.ndarray, + tail_damping_ratio: float, ) -> Path: output_dir = ensure_parent_dir(output_dir) peak_times = time_values[selected_peak_indices] @@ -111,7 +222,18 @@ def save_task3_figure( axes[2].scatter(peak_times, selected_peak_amplitudes, color="tab:red", s=28, label="Selected envelope peaks") axes[2].plot(envelope_time, envelope, color="tab:green", linewidth=1.8, label="Fitted positive envelope") axes[2].plot(envelope_time, -envelope, color="tab:green", linewidth=1.2, linestyle="--", label="Fitted negative envelope") - axes[2].set_title(f"Selected decay segment | damping ratio zeta = {damping_ratio:.5f}") + global_peak_times = time_values[global_peak_indices] + axes[2].scatter(global_peak_times, global_peak_amplitudes, color="tab:purple", s=14, alpha=0.65, label="Global peaks") + axes[2].plot(global_peak_times, global_fit_amplitudes, color="tab:purple", linewidth=1.2, linestyle="-.", label="Global envelope fit") + tail_peak_times = time_values[tail_peak_indices] + axes[2].scatter(tail_peak_times, tail_peak_amplitudes, color="tab:brown", s=14, alpha=0.65, label="Max-peak tail peaks") + axes[2].plot(tail_peak_times, tail_fit_amplitudes, color="tab:brown", linewidth=1.2, linestyle=":", label="Max-peak tail fit") + axes[2].set_title( + "Damping ratio | " + f"local zeta={local_damping_ratio:.5f} | " + f"global zeta={global_damping_ratio:.5f} | " + f"max-tail zeta={tail_damping_ratio:.5f}" + ) axes[2].set_xlabel("Time (s)") axes[2].set_ylabel("Acceleration") axes[2].grid(True, alpha=0.25) @@ -154,11 +276,16 @@ def analyze_free_vibration( selected_peak_indices = np.asarray(peak_window["peak_indices"], dtype=int) selected_peak_amplitudes = np.asarray(peak_window["peak_amplitudes"], dtype=np.float64) fit_amplitudes = np.asarray(peak_window["fit_amplitudes"], dtype=np.float64) - log_decrements = np.log(selected_peak_amplitudes[:-1] / selected_peak_amplitudes[1:]) - mean_log_decrement = float(np.mean(log_decrements)) - damping_ratio = float(mean_log_decrement / (2.0 * np.pi)) - damping_ratio_exact = float( - mean_log_decrement / np.sqrt((2.0 * np.pi) ** 2 + mean_log_decrement**2) + mean_log_decrement, damping_ratio, damping_ratio_exact = log_decrement_metrics(selected_peak_amplitudes) + global_metrics = estimate_global_decay_metrics( + filtered_signal=filtered_signal, + peak_indices=peak_indices, + sampling_rate=sampling_rate, + ) + tail_metrics = estimate_tail_from_max_peak_metrics( + filtered_signal=filtered_signal, + peak_indices=peak_indices, + sampling_rate=sampling_rate, ) figure_path = save_task3_figure( @@ -172,7 +299,15 @@ def analyze_free_vibration( selected_peak_amplitudes=selected_peak_amplitudes, fit_amplitudes=fit_amplitudes, natural_frequency_hz=natural_frequency_hz, - damping_ratio=damping_ratio, + local_damping_ratio=damping_ratio, + global_peak_indices=np.asarray(global_metrics["peak_indices"], dtype=int), + global_peak_amplitudes=np.asarray(global_metrics["peak_amplitudes"], dtype=np.float64), + global_fit_amplitudes=np.asarray(global_metrics["fit_amplitudes"], dtype=np.float64), + global_damping_ratio=float(global_metrics["damping_ratio"]), + tail_peak_indices=np.asarray(tail_metrics["peak_indices"], dtype=int), + tail_peak_amplitudes=np.asarray(tail_metrics["peak_amplitudes"], dtype=np.float64), + tail_fit_amplitudes=np.asarray(tail_metrics["fit_amplitudes"], dtype=np.float64), + tail_damping_ratio=float(tail_metrics["damping_ratio"]), ) if show: plt.show() @@ -191,6 +326,16 @@ def analyze_free_vibration( "mean_log_decrement": mean_log_decrement, "damping_ratio": damping_ratio, "damping_ratio_exact": damping_ratio_exact, + "global_selected_peak_count": int(len(np.asarray(global_metrics["peak_indices"], dtype=int))), + "global_mean_log_decrement": float(global_metrics["mean_log_decrement"]), + "global_damping_ratio": float(global_metrics["damping_ratio"]), + "global_damping_ratio_exact": float(global_metrics["damping_ratio_exact"]), + "global_envelope_r2": float(global_metrics["r_squared"]), + "max_tail_selected_peak_count": int(len(np.asarray(tail_metrics["peak_indices"], dtype=int))), + "max_tail_mean_log_decrement": float(tail_metrics["mean_log_decrement"]), + "max_tail_damping_ratio": float(tail_metrics["damping_ratio"]), + "max_tail_damping_ratio_exact": float(tail_metrics["damping_ratio_exact"]), + "max_tail_envelope_r2": float(tail_metrics["r_squared"]), "envelope_r2": float(peak_window["r_squared"]), "figure_path": str(figure_path), } @@ -202,12 +347,22 @@ def print_report(result: dict[str, float | int | str | Path]) -> None: print(f"Top response sensor: {result['sensor_code']} / {result['axis']}") print(f"Sampling rate: {result['sampling_rate_hz']:.4f} Hz") print(f"Natural frequency f_n: {result['natural_frequency_hz']:.6f} Hz") - print(f"Mean log decrement delta: {result['mean_log_decrement']:.6f}") - print(f"Damping ratio zeta (delta / 2pi): {result['damping_ratio']:.6f}") - print(f"Damping ratio exact: {result['damping_ratio_exact']:.6f}") + print(f"Local mean log decrement delta: {result['mean_log_decrement']:.6f}") + print(f"Local damping ratio zeta (delta / 2pi): {result['damping_ratio']:.6f}") + print(f"Local damping ratio exact: {result['damping_ratio_exact']:.6f}") + print(f"Global mean log decrement delta: {result['global_mean_log_decrement']:.6f}") + print(f"Global damping ratio zeta (delta / 2pi): {result['global_damping_ratio']:.6f}") + print(f"Global damping ratio exact: {result['global_damping_ratio_exact']:.6f}") + print(f"Max-tail mean log decrement delta: {result['max_tail_mean_log_decrement']:.6f}") + print(f"Max-tail damping ratio zeta (delta / 2pi): {result['max_tail_damping_ratio']:.6f}") + print(f"Max-tail damping ratio exact: {result['max_tail_damping_ratio_exact']:.6f}") print(f"Detected peaks: {result['detected_peak_count']}") - print(f"Selected peaks for envelope: {result['selected_peak_count']}") - print(f"Envelope linearity R^2: {result['envelope_r2']:.6f}") + print(f"Selected peaks for local envelope: {result['selected_peak_count']}") + print(f"Selected peaks for global envelope: {result['global_selected_peak_count']}") + print(f"Selected peaks for max-tail envelope: {result['max_tail_selected_peak_count']}") + print(f"Local envelope linearity R^2: {result['envelope_r2']:.6f}") + print(f"Global envelope linearity R^2: {result['global_envelope_r2']:.6f}") + print(f"Max-tail envelope linearity R^2: {result['max_tail_envelope_r2']:.6f}") print(f"Figure saved to: {result['figure_path']}")