RDown excitation time

Quit[]

Setup

Load package and set root paths

rootpath="~/git/rdstackingproject/";

SetDirectory[rootpath<>"code"];
<<RDown.m
<<mcmc.m

mysxscase="/Users/xisco/SXS/BBH_SKS_d14.3_q1.22_sA_0_0_0.330_sB_0_0_-0.440";
tshift=0;
npoints=10;
ntones=1;
ωfact=0.05;

Setup paths and other notebook folders

(* Change them accordingly *)
(*Notebook directory*)
notdir=rootpath<>"code/";
(*Directory of the l's Ringdown data *)
modedir=rootpath;
(* True if you want to export the plots *)
ExportQ=False
ExportDir=notdir

Plot Formatting

(* PlotStyle *)
markers={"●","■","◆","▲","▼","○","□","◇","△","▽"}
colors=ColorData[97,"ColorList"]
l1=.025;
l2=.001;
s=l1/2;

Other functions

Options[pdf]={"SNR"->1};
pdf[ansatz_,data_,tend_?NumericQ,x_?(VectorQ[#,NumericQ]&),OptionsPattern[]]:=Block[{ansl,cfit,cfitd,h1red,h2red,norm1,norm2,myTable,snr,tmax},

snr=OptionValue["SNR"];

tmax=data[[1,1]];
ansl=ansatz;
cfitd=Transpose[{data[[All,1]],ansl/.t->data[[All,1]]}];

h1red=Select[data,tmax<=#[[1]]<=tmax+tend&][[All,2]];
h2red=Select[cfitd,tmax<=#[[1]]<=tmax+tend&][[All,2]];
norm1=Total[(Abs@h1red)^2 ];
norm2=Total[(Abs@h2red)^2 ];

myTable=h1red Conjugate@h2red;
-snr^2(1-(Re@Total[myTable]/Sqrt[norm1 norm2]))]

Options[pdfFit]={"SNR"->1};
pdfFit[ansatz_,data_,tend_?NumericQ,x_?(VectorQ[#,NumericQ]&),OptionsPattern[]]:=Block[{ansl,cfit,cfitd,h1red,h2red,norm1,norm2,myTable,snr,tmax},

snr=OptionValue["SNR"];

tmax=data[[1,1]];
ansl=ansatz;
cfitd=Transpose[{data[[All,1]],ansl/.t->data[[All,1]]}];
h1red=Select[data,tmax<=#[[1]]<=tmax+tend&][[All,2]];
h2red=Select[cfitd,tmax<=#[[1]]<=tmax+tend&][[All,2]];

Total[(Re[h1red-h2red]^2+Im[h1red-h2red]^2)]
]

Imports

RDown data

(* The QNM spectrum is obtained from https://pages.jh.edu/~eberti2/ringdown/. Up to now, we can load any lmn combination up to l=4 and n=7. To extend the l number you need to download the l>4 files. *)

(* Import all the modes you want to use *)
modelist={{2,2,0},{2,2,1},{2,2,2},{2,2,3},{2,2,4},{2,2,5},{2,2,6},{2,2,7}}

(* Files for the modes. Note that in this file notation the fundamental tone corresponds to n=1. However, for the code notation we will always use n=0.  *)
Modefiles=Table[modedir<>"l"<>ToString[modelist[[j,1]]]<>"/n"<>ToString[modelist[[j,3]]+1]<>"l"<>ToString[modelist[[j,1]]]<>"m"<>ToString[modelist[[j,2]]]<>".dat",{j,Length@modelist}]

(* Import modes data *)
Modedata=Import/@Modefiles;

NR data + EOBFits

modes={{2,2}}

mysxscaserh=Reverse[Select[FileNames["*Lev6/rh*",mysxscase,2],Not@StringMatchQ[#,"*raw*"]&]]
mysxscasemetafile=FileNames["metadata.txt",mysxscase,2][[1]]
metadata=SXSMetaFilesToRules[mysxscasemetafile];

Print["mass1 = ",mass1=("initial-mass1"/.metadata)[[1]]]
Print["mass2 = ",mass2=("initial-mass2"/.metadata)[[1]]]
Print["χ1 = ",χ1=Chop[(("initial-dimensionless-spin1"/.metadata)[[-1]])]]
Print["χ2 = ",χ2=Chop[(("initial-dimensionless-spin2"/.metadata)[[-1]])]]
Print["m1/m2 = ",q=Max[{mass1/mass2,1}]]
Print["m1·m2 = ",η=q/(1+q)^2*1.]
Print["af = ",af=("remnant-dimensionless-spin"/.metadata)[[-1]]]
Print["mf = ",mf=("remnant-mass"/.metadata)[[1]]]
Print["Tag = ",tag=("alternative-names"/.metadata)[[2]]]

sxsrhs=Flatten[Conjugate/@GetAsymptoticMultiMode[#,2,modes,"ReSample"->True]&/@mysxscaserh,1];

Results for the paper

Vary the fundamental frequency

t0=TimeOfMaximum[sxsrhs[[1]]]+tshift
data=Select[sxsrhs[[1]],#[[1]]≥ t0&];

tab=Table[
ansatz=OvertoneModelV2[ntones,{η,χ1,χ2},t0,"Fitα"->{},"Varyω"->True,"ωval"->{-ωfact,ωfact},"Export_ωval"->True,"ModesData"->Modedata];
cfit=NonlinearModelFit[data,ansatz[[1]],{x0,x1,x2,x3,x4,x5,x6,x7},t];
cfitd=Transpose[{data[[All,1]],Normal[cfit]/.t->data[[All,1]]}];
{ansatz[[2]],{x0,x1,x2,x3,x4,x5,x6,x7}/.cfit["BestFitParameters"],1-EasyMatchT[data,cfitd,t0,t0+90]},{j,npoints}];

Export[rootpath<>"plots/freq_bfitpars_mmatch.dat",tab]