Skip to content

Instantly share code, notes, and snippets.

@pukpr
Created June 28, 2025 05:29
Show Gist options
  • Select an option

  • Save pukpr/492a57af8a7929ab2eece6994ab0a5f4 to your computer and use it in GitHub Desktop.

Select an option

Save pukpr/492a57af8a7929ab2eece6994ab0a5f4 to your computer and use it in GitHub Desktop.
Regression fit for ENSO
Doodson_Set : Doodson_List :=
(
(1, 0, 0, 0.0, 1.0), -- T
(1, 0, 0, 1.0, 1.0), -- D
-- (0, 0, 2, 2.0, 1.0), -- 3Y
(0, 0,-2, 2.0, 1.0), -- 5.7
-- (0, 0, 0, 2.0, 1.0), -- N/2
(2, 0, 0, 1.0, 1.0), -- DT
(2, 0, 0, 0.0, 1.0), -- t
(2, 0, 0, 2.0, 1.0), -- d
-- (1, 0,-1, 0.0, 1.0), -- A
-- (2, 0,-2, 0.0, 1.0), -- a
(0, 0, 0, 1.0, 1.0), -- N
(0, 0, 2, 0.0, 1.0), -- p
-- (1, 0, 1, 0.0, 1.0), -- E
(2, 0,-1, 1.0, 1.0), -- AD
(0, 0, 2, 1.0, 1.0), -- 3.56
-- (1, 0,-1, 1.0, 1.0), -- A+
-- (1, 0,-1,-1.0, 1.0), -- A-
(0, 0, 1, 1.0, 1.0), -- 6Y
-- (1, 0, 1, 1.0, 1.0), -- 26.985
-- (0, 0, 1,-1.0, 1.0),-- 16.8y
(0, 0,-1, 0.0, 1.0) -- P
--(1, 0, 1, -1.0, 1.0) -- Fake
-- (3, 0,-1, 1.0, 1.0), -- Mt'
-- (3, 0,-1, 0.0, 1.0), -- Mt
-- (3, 0,-1, 2.0, 1.0) -- Mtd
);
Doodson_RSet : Doodson_List :=
(
(1, 0, 0, 0.0, 1.0), -- T
(1, 0, 0, 1.0, 1.0), -- D
-- (0, 0, 2, 2.0, 1.0), -- 3Y
-- (0, 0,-2, 2.0, 1.0), -- 5.7
-- (0, 0, 0, 2.0, 1.0), -- N/2
(2, 0, 0, 1.0, 1.0), -- DT
(2, 0, 0, 0.0, 1.0), -- t
(2, 0, 0, 2.0, 1.0), -- d
(1, 0,-1, 0.0, 1.0), -- A
(1, 0, 1, 0.0, 1.0), -- evect
-- (0, 0, 0, 1.0, 1.0), -- N
-- (0, 0, 2, 0.0, 1.0), -- p
-- (1, 0, 1, 0.0, 1.0), -- E
(2, 0,-1, 1.0, 1.0) -- AD
-- (0, 0, 2, 1.0, 1.0), -- 3.56
-- (1, 0,-1, 1.0, 1.0), -- A+
-- (1, 0,-1,-1.0, 1.0), -- A-
-- (0, 0, 1, 1.0, 1.0), -- 6Y
-- (1, 0, 1, 1.0, 1.0), -- 26.985
-- (0, 0, 1,-1.0, 1.0),-- 16.8y
-- (0, 0,-1, 0.0, 1.0) -- P
);
Annual_Set : Doodson_List :=
(
(0, 1, 0, 0.0, 365.242), -- annual
(0, 2, 0, 0.0, 365.242) -- semi
);
with Text_IO;
with GEM.LTE.Primitives;
function GEM.Mix_Regression (File_Name : in String) return Gem.LTE.Period_Set is
use GEM.LTE, GEM.LTE.Primitives;
D : Data_Pairs := Make_Data(File_Name);
First, Last : Integer;
Singular : Boolean;
Forcing : Data_Pairs := D;
Model : Data_Pairs := D;
P_A1, P_A2, P_A3, P_A4, P_A5 : Gem.LTE.Long_Periods := Gem.LTE.LP_Set;
P_B1 : Gem.LTE.Long_Periods := Gem.LTE.LP_RSet;
P_Annual : Gem.LTE.Long_Periods := Gem.LTE.LP_Annual;
A_A1, A_A2, A_A3, A_A4, A_A5 : Gem.LTE.Long_Periods_Amp_Phase := Gem.LTE.LPAP_Set;
A_B1 : Gem.LTE.Long_Periods_Amp_Phase := Gem.LTE.LPAP_RSet;
A_Annual : Gem.LTE.Long_Periods_Amp_Phase := Gem.LTE.LPAP_Annual;
Level, K0, Trend, Accel : Long_Float := 0.0;
Last_Time : Long_Float;
Freq, Freq2 : Long_Float;
One : constant Long_Float := 1.0;
begin
First := D'First;
Last := D'Last;
Text_IO.Put_Line("records from " & First'Img & " ... " & Last'Img);
Text_IO.Put_Line("factors from " & P_A1'First'Img & " ... " & P_A1'Last'Img);
for I in Forcing'Range loop
Forcing(I).Value := Forcing(I).Date;
Last_Time := Forcing(I).Value;
end loop;
Text_IO.Put("updated forcing ");
Put(Last_Time);
Text_IO.New_Line;
for I in P_A1'Range loop -- to aliased frequency
Freq := GEM.LTE.Year_Length/P_A1(I);
Freq := abs Long_Float'Remainder(Freq,One);
P_A1(I) := Freq;
P_A2(I) := One - Freq;
P_A3(I) := One + Freq;
P_A4(I) := 2.0*One - Freq;
P_A5(I) := 2.0*One + Freq;
end loop;
for I in P_B1'Range loop -- to aliased frequency
Freq := 2.0*GEM.LTE.Year_Length/P_B1(I);
Freq2 := 4.0*GEM.LTE.Year_Length/P_B1(I);
Text_IO.Put_Line("Q=" & Freq2'Image);
--Text_IO.Put_Line("F=" & Freq'Image);
if Integer(Freq) mod 2 = 1 then
Freq := abs Long_Float'Remainder(Freq,One);
else
Freq := One - abs Long_Float'Remainder(Freq,One);
end if;
P_B1(I) := Freq*0.5;
-- P_B2(I) := (One-Freq)*0.5;
-- if Integer(Freq2) mod 2 = 1 then
-- Freq2 := abs Long_Float'Remainder(Freq2,One);
-- else
-- Freq2 := One - abs Long_Float'Remainder(Freq2,One);
-- end if;
-- P_Q1(I) := Freq2*0.25;
end loop;
declare
P : Gem.LTE.Long_Periods := P_A1 & P_A2 & P_A3 & P_A4 & P_A5 & P_B1
& P_Annual; -- & P_Q1;
A : Gem.LTE.Long_Periods_Amp_Phase := A_A1 & A_A2 & A_A3 & A_A4 & A_A5 & A_B1
& A_Annual; -- & A_Q1;
begin
Trend := 0.0;
Regression_Factors (Data_Records => D(800..1700), -- Time series
Forcing => Forcing(800..1700), -- Value @ Time
NM => P'Last, -- # modulations
DBLT => P, --D.B.LT,
DALTAP => A, --D.A.LTAP,
DALEVEL => LEVEL,
DAK0 => K0,
Secular_Trend => Trend,
Accel => Accel,
Singular => Singular
);
Text_IO.Put_Line("Singular? " & Singular'Img);
for I in P'Range loop
Put(P(I)); Text_IO.Put(" ");
-- INTEGRATE dLOD -- DBLTAP(I).Amplitude := DBLTAP(I).Amplitude * 26.736 / DBLT(I);
Put(A(I).Amplitude); Text_IO.Put(" ");
-- INTEGRATE dLOD -- DBLTAP(I).Phase := DBLTAP(I).Phase + Ada.Numerics.Pi/2.0;
Put(A(I).Phase); Text_IO.Put(" ");
Text_IO.New_Line;
end loop;
Model := Gem.LTE.Primitives.LTE(Forcing => Forcing,
Wave_Numbers => P,
Amp_Phase => A,
Offset => Level,
K0 => K0,
Trend => 0.0,
NonLin => 1.0);
Put(CC(D,Model), "=CC " );
Put(Year_Length, "=Yr", True);
GEM.LTE.Primitives.Save(Model, D, Forcing);
declare
PS : constant Gem.LTE.Period_Set := (N => P'Length,
LP => P,
AP => A);
begin
return PS;
end;
--return A;
end;
end;
@pukpr

pukpr commented Jun 30, 2025

Copy link
Copy Markdown
Author

QBO30

No biennial impulse
image

YEAR 0.01
F9 not set, use default
DECAY not set, use default
FMULT 1.008
FSTART 0.01
LOCKF not set, use default
R2M not set, use default
SYM not set, use default
YTRIM not set, use default
NONLIN not set, use default
records from  1 ...  856
factors from  1 ...  15
updated forcing  2024.25000000000
Fit from 1 ..  258 and  514 ..  856
Fit from 1 ..  258 and  514 ..  856
Singular? FALSE
   0.36859574510      0.09176074234      0.59460567384
   0.42236089010      0.32688457329     -2.35230583931
   0.11844980964      0.04146858432     -0.79565833396
   0.20904336480      0.08765088043      1.07586424669
   0.26280850979      0.09061042109     -0.84763216129
   0.15527821980      0.03000278697      1.08641949236
   0.48878860942      0.16136965430      2.54309781449
   0.05376514499      0.00348794471     -0.40172527475
   0.22598009963      0.09207088321     -2.80249182601
   0.48158579492      0.03736068208      2.96113138962
   0.32203341461      0.15972268766     -1.45987905213
   0.27974524463      0.11725149760      1.16570379376
   0.16675519481      0.05644752208      1.82696555405
   0.11299004982      0.05337525575      1.82122542732
   0.03682841016      0.02894589673      0.52946601823
   0.63140425490      0.01363139922     -1.55907366928
   0.57763910990      0.04370722457      1.15689473836
   0.88155019036      0.01705535440     -2.13231138399
   0.79095663520      0.01849567088      0.39927058578
   0.73719149021      0.02919636448     -2.94765275019
   0.84472178020      0.00863529718     -2.61180085440
   0.51121139058      0.10191570270      2.40240040688
   0.94623485501      0.02158731072     -1.58464801698
   0.77401990037      0.02335346780      0.88461412244
   0.51841420508      0.12434780154      2.91470958046
   0.67796658539      0.07068836025     -0.79230747428
   0.72025475537      0.01543273284     -1.37304516848
   0.83324480519      0.02309086967     -0.86750995678
   0.88700995018      0.02211781964      2.30615161610
   0.96317158984      0.06046584490     -0.76933823746
   1.36859574510      0.01748735885      3.07118226465
   1.42236089010      0.03568645917      2.63079489203
   1.11844980964      0.05037841260      0.07481691527
   1.20904336480      0.03197060778      0.45660069537
   1.26280850979      0.03922421733     -0.50303774794
   1.15527821980      0.02602688112     -3.13237305912
   1.48878860942      0.04163727986     -0.18642358879
   1.05376514499      0.04363486400     -0.78658794094
   1.22598009963      0.01191088670      1.55623471059
   1.48158579492      0.02816128580     -1.12242566585
   1.32203341461      0.00587704250     -1.51400449582
   1.27974524463      0.02059617864      2.52221789575
   1.16675519481      0.01676918483      1.95591770443
   1.11299004982      0.05476522123      2.33624184808
   1.03682841016      0.06126020758      2.48835585705
   1.63140425490      0.01583669221      1.68633302280
   1.57763910990      0.00696831673     -2.22978299116
   1.88155019036      0.00115465531      3.04959628475
   1.79095663520      0.01581801038      0.61267629123
   1.73719149021      0.00421107741      2.85577521474
   1.84472178020      0.01997713202      1.10052198164
   1.51121139058      0.03152768717     -2.80545097340
   1.94623485501      0.00980297825     -2.37195910271
   1.77401990037      0.02380405486     -1.03647150618
   1.51841420508      0.02107816206     -0.89125657713
   1.67796658539      0.01049575892      2.31007820272
   1.72025475537      0.01186077613      1.78524387871
   1.83324480519      0.00863747263      2.86627526444
   1.88700995018      0.01106126067     -1.71305895022
   1.96317158984      0.02151625297     -0.26675814845
   2.36859574510      0.00088650127      2.09573936987
   2.42236089010      0.03195288755     -0.16488820558
   2.11844980964      0.00889059986     -2.67334580406
   2.20904336480      0.01211196789     -2.61203648015
   2.26280850979      0.00697962904      2.45314040851
   2.15527821980      0.00630024838     -0.07544005133
   2.48878860942      0.00766615024     -0.72265449035
   2.05376514499      0.01726645072     -0.10405328436
   2.22598009963      0.00513756351     -1.33466723729
   2.48158579492      0.01425318869     -1.42498651973
   2.32203341461      0.00616410171      0.27701111322
   2.27974524463      0.00660483303     -2.73513012879
   2.16675519481      0.00252149463      1.37171315637
   2.11299004982      0.01096546620     -2.17412998162
   2.03682841016      0.00671920614      0.61837585542
   2.63140425490      0.00496943659      1.23934645015
   2.57763910990      0.00163968616      0.57067576195
   2.88155019036      0.00683582745      1.57581296289
   2.79095663520      0.00625823573     -2.72804977197
   2.73719149021      0.00525229732     -0.72799499066
   2.84472178020      0.01031803268      2.95282533173
   2.51121139058      0.01987970332     -3.05157466083
   2.94623485501      0.00557032909      2.14995076037
   2.77401990037      0.00607393623      1.41971621276
   2.51841420508      0.01738981250     -1.81264901414
   2.67796658539      0.00544440923     -0.57197002469
   2.72025475537      0.00794737072     -1.96805326245
   2.83324480519      0.00130390661      2.89642092636
   2.88700995018      0.00541852940     -2.21936888869
   2.96317158984      0.00514681211     -1.72465462714
   0.99997194158      0.09654709024      1.10223953481
   1.99994388317      0.02951583341      0.80607667310
   2.99991582475      0.00496906360      1.35152142200
   0.63741991609=CC Fit from 1 ..  258 and  514 ..  856
Fit from 1 ..  258 and  514 ..  856
   0.75348350761=CC(training) Exclude Fit from  258 ..  514
Exclude Fit from  258 ..  514
   0.43307680115=CC(validate)  365.25224840000=Yr
1974.41666700000=start 1995.75000000000=finish

image

Chandler Wobble

image

TREND FALSE
YEAR 0.01
F9 not set, use default
DECAY not set, use default
FMULT 1.008
FSTART 0.01
LOCKF not set, use default
R2M not set, use default
SYM not set, use default
YTRIM not set, use default
NONLIN not set, use default
records from  1 ...  2381
factors from  1 ...  15
updated forcing  2019.05215400000
Fit from  1310 ..  2143
Fit from  1310 ..  2143
Singular? FALSE
   0.36859574510      0.00175027426      1.06715843180
   0.42236089010      0.00147009444     -3.09471412452
   0.11844980964      0.00486478974      0.57638151116
   0.20904336480      0.00119706206      2.81616730234
   0.26280850979      0.00142428472      0.22266109211
   0.15527821980      0.00349385901      2.78446176124
   0.48878860942      0.03145036924     -0.62589361457
   0.05376514499      0.00207444041      0.90425171576
   0.22598009963      0.00118499213     -2.80222468337
   0.48158579492      0.02350906281     -1.62505991966
   0.32203341461      0.00198616497     -2.59119567908
   0.27974524463      0.00384067096      2.80689614773
   0.16675519481      0.00386819570      0.53417725024
   0.11299004982      0.00495668886      2.81839232184
   0.03682841016      0.00172148097      1.85771547072
   0.63140425490      0.00347730698     -2.69663369455
   0.57763910990      0.00332495708     -2.46657468023
   0.88155019036      0.12966453438      0.16576789433
   0.79095663520      0.03899642747      0.66674012944
   0.73719149021      0.02768103991     -2.95894554589
   0.84472178020      0.29868971727      0.64613987502
   0.51121139058      0.03775740486     -0.88711757472
   0.94623485501      0.04490002007     -0.14242253545
   0.77401990037      0.02617756791      1.33795581380
   0.51841420508      0.02667320891      0.23858593110
   0.67796658539      0.01369371245     -0.29861973891
   0.72025475537      0.00531333220     -1.44123625761
   0.83324480519      0.08043550317      0.28033287656
   0.88700995018      0.11353150155     -2.06865441425
   0.96317158984      0.02876998843     -0.88408909477
   1.36859574510      0.00349364883      0.76649906955
   1.42236089010      0.00107945770     -1.35234681047
   1.11844980964      0.03056776685     -0.04236044476
   1.20904336480      0.01234032479     -1.18574905092
   1.26280850979      0.00973360197      0.33583680655
   1.15527821980      0.01667438357      2.51519243834
   1.48878860942      0.01080009324     -0.51079561641
   1.05376514499      0.01115291086     -2.07685695636
   1.22598009963      0.00673950480     -1.23615296943
   1.48158579492      0.00982021023     -1.23727093636
   1.32203341461      0.00705405370     -0.88742165655
   1.27974524463      0.00793632234     -1.31115752369
   1.16675519481      0.00714876232      0.29556924410
   1.11299004982      0.03969265870      2.17845509333
   1.03682841016      0.02055401761      1.02208692612
   1.63140425490      0.00275013652     -2.51756059612
   1.57763910990      0.00695336982      2.96798616953
   1.88155019036      0.02004064133     -0.38372031857
   1.79095663520      0.00560903404     -0.67907479229
   1.73719149021      0.00132146546     -0.23696216302
   1.84472178020      0.00716077223      2.77057233615
   1.51121139058      0.01457158437     -0.59038199002
   1.94623485501      0.00826734890     -0.02467973222
   1.77401990037      0.00136431845      1.10974891820
   1.51841420508      0.01207870238      0.78153426439
   1.67796658539      0.00468276562      2.17772843797
   1.72025475537      0.00328462937     -0.51953309082
   1.83324480519      0.00920347421     -1.19925661614
   1.88700995018      0.01918836146     -2.75647966947
   1.96317158984      0.01188217734     -1.80568024055
   2.36859574510      0.00658165103      0.60909315415
   2.42236089010      0.00232530045     -0.73618111129
   2.11844980964      0.00280586245      1.56887417930
   2.20904336480      0.00489998501      1.92663145412
   2.26280850979      0.00229949669     -0.69830917026
   2.15527821980      0.00262898059     -1.55211931770
   2.48878860942      0.01782331727     -1.32417466206
   2.05376514499      0.00546392469      2.66962915942
   2.22598009963      0.00220322659     -0.55200462188
   2.48158579492      0.01870346951     -2.22742382425
   2.32203341461      0.00235342899     -2.24325964525
   2.27974524463      0.00380091533     -2.18003636026
   2.16675519481      0.00502471435      3.11551371097
   2.11299004982      0.00554022874     -2.54307028084
   2.03682841016      0.00215508333      2.96813055258
   2.63140425490      0.00208338125     -2.48195804102
   2.57763910990      0.00382491360     -0.03957754903
   2.88155019036      0.00609622153      2.27713942410
   2.79095663520      0.00247816802     -0.81376430436
   2.73719149021      0.00396635044      0.25070056311
   2.84472178020      0.00338553344     -0.60452818664
   2.51121139058      0.00944359998     -1.44016207336
   2.94623485501      0.00514316260     -2.07064165373
   2.77401990037      0.00209431434      0.85152414697
   2.51841420508      0.00940320680     -0.12669477545
   2.67796658539      0.00273940914      2.36163806369
   2.72025475537      0.00463260228      1.64345973327
   2.83324480519      0.00459224921      2.63180898860
   2.88700995018      0.00456473895     -0.62609797339
   2.96317158984      0.00373732414     -3.00058993330
   0.99997194158      0.19365254533      0.57989743045
   1.99994388317      0.00752372114     -2.07322775336
   2.99991582475      0.00977458794     -2.22059278872
   0.40922556321=CC Fit from  1310 ..  2143
Fit from  1310 ..  2143
   0.98677625230=CC(training) Exclude fit from 1 ..  1310 and  2143 ..  2381
Exclude fit from 1 ..  1310 and  2143 ..  2381
   0.18896813059=CC(validate)  365.25224840000=Yr
1965.50215400000=start 2007.15215400000=finish


@pukpr

pukpr commented Jul 1, 2025

Copy link
Copy Markdown
Author

CW fingerprint. Note how the draconic CW frequency at 0.845 is replicated at 1.845 -- modulation of an annual impulse train.
Note the minor peaks at 0.5, 1.5, 2.5 -- that's a biennial modulation of the annual carrier
Also minor peaks at 1.678, 2.678 -- related to the 1/D+1/A beat of the 1/D-1/A 6 year frequency.

image

Re-evaluation of the experiment described here https://geoenergymath.com/2021/01/07/chandler-wobble-forcing/


Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment