21 #pragma GCC system_header 28 #include "TObjArray.h" 29 #include "TObjString.h" 30 #include "TPaletteAxis.h" 31 #include "TMultiLayerPerceptron.h" 32 #include "TMLPAnalyzer.h" 60 sprintf(NNi2,
"%s",NN_FILE.Data());
65 if((NNi2[p+1]==
'N')&&(NNi2[p+2]==
'\/')) {
68 while (NNi2[hh]!=
'\0'){NNi[hh-p-3]=NNi2[
hh];hh=hh+1;}
79 sprintf(Filei2,
"%s",TEST_FILE.Data());
82 if(Filei2[q]==
'n'&&Filei2[q+1]==
'n'&&Filei2[q+2]==
'T'&&Filei2[q+3]==
'R'&&Filei2[q+4]==
'E'&&(Filei2[q+5]==
'E')&&(Filei2[q+6]==
'\/')) {
85 while (Filei2[hh]!=
'\0') {Filei[hh-p-7]=Filei2[
hh];hh=hh+1;
88 for (
int h0=hh-p-7;h0<1024;h0++) Filei[h0]=
'\0';
94 cout<<NNi<<
" original String: "<<NNi2<<endl;
95 cout<<Filei<<
" original String: "<<Filei2<<endl;
100 NNi0.ReplaceAll(
".root",
"");
101 Filei0.ReplaceAll(
".root",
"");
104 TFile* fnet =TFile::Open(NN_FILE.Data());
105 TMultiLayerPerceptron *mlp = (TMultiLayerPerceptron*)fnet->Get(
"TMultiLayerPerceptron");
106 if(mlp==NULL) {cout <<
"Error getting mlp" << endl;
exit(1);}
107 TTree* infot=(TTree*)fnet->Get(
"info");
112 infot->SetBranchAddress(
"Rand_start_Sig",&NNs);
113 infot->SetBranchAddress(
"Rand_start_Bg",&NNb);
114 infot->SetBranchAddress(
"#trainSig",&NNnTS);
115 infot->SetBranchAddress(
"#trainBg",&NNnTB);
118 TFile* fTEST =TFile::Open(TEST_FILE.Data());
119 TTree* NNTree=(TTree*)fTEST->Get(
"nnTree");
120 int entries=NNTree->GetEntries();
121 cout<<
"entries: "<<entries<<endl;
131 NNTree->SetBranchAddress(
"#Entries_type",&entriesTot);
132 NNTree->SetBranchAddress(
"Matrix_dim",&ndim);
133 NNTree->SetBranchAddress(
"#inputs",&ninp);
134 NNTree->SetBranchAddress(
"amplitude_mode",&y);
137 sig_entries=entriesTot;
138 NNTree->GetEntry(entries-1);
139 bg_entries=entriesTot;
144 cout<<
"NDIM "<<NDIM<<endl;
145 cout<<
"nINP"<<nINP<<endl;
146 cout<<
"sig e "<<sig_entries<<endl;
147 cout<<
"bg e "<<bg_entries<<endl;
149 if (sig_entries>bg_entries) minevents=bg_entries;
150 else minevents=sig_entries;
152 if(b==0) b=sig_entries;
158 cout<<
"Error: Bg index<sig_entries"<<endl;
161 if((TS>sig_entries||
TB>bg_entries)&&(TS==
TB)) {TS=minevents-
s;
TB=minevents-b+sig_entries;}
162 if((TS>sig_entries||
TB>bg_entries)&&(TS!=
TB)) {TS=sig_entries-
s;
TB=bg_entries-b+sig_entries;}
166 sprintf(NOMEtot,
"NN_%s_FILE_%s_TS%i_TB%i_st%i_bt%i_dANN%1.2f_uf%i",NNi0.Data(),Filei0.Data(),TS,
TB,
s,b,
deltaANN,uf);
167 cout<<
"nome: "<<NOMEtot<<endl;
172 NNTree->SetBranchAddress(
"Files_name",&FILE_NAME);
175 NNTree->GetEntry(entries-1);
177 cout<<
"fine ifdef RHO_CC"<<endl;
179 TChain sigTree(
"waveburst");
180 sigTree.Add(SIG_FILE.Data());
183 cout <<
"sig entries2 : " << sig_entries2 << endl;
185 TChain bgTree(
"waveburst");
186 bgTree.Add(BG_FILE.Data());
189 cout <<
"bg entries2 : " << bg_entries2 << endl;
191 cout<<
"b: "<<b<<endl;
192 cout<<
"s: "<<s<<endl;
196 for (
int jj=0; jj<nINP;jj++) x[jj]=0.;
197 char ilabel[nINP][16];
200 for(
int i=0;
i<nINP;
i++) {
202 NNTree->SetBranchAddress(ilabel[i], &x[i]);
206 sprintf(ofile,
"outfile/%s.root",NOMEtot);
208 TFile*
f=
new TFile(ofile,
"RECREATE");
209 TTree* NNTree2=
new TTree(
"Parameters",
"Parameters");
210 NNTree2->SetDirectory(f);
214 NNTree2->Branch(
"ANNout",&out,
"ANNout/D");
215 NNTree2->Branch(
"cc",&NNcc,
"cc/D");
216 NNTree2->Branch(
"rho",&NNrho,
"rho/D");
217 char NNfilename[1024];
218 NNTree2->Branch(
"NNname",&NNfilename,
"NNname/C");
219 sprintf(NNfilename,
"%s",NN_FILE.Data());
221 NNTree2->Branch(
"TestFile",&Testf,
"TestFile/C");
222 sprintf(Testf,
"%s",TEST_FILE.Data());
224 NNTree2->Branch(
"#TestSig",&nTestS,
"#TestSig/I");
226 cout<<
"nTestS: "<<nTestS<<
" TS: "<<TS<<endl;
227 cout<<
"dopo def tree"<<endl;
243 for(
int n=s;
n<s+TS;
n++) {
244 if (uf!=0&&
n>=NNs&&
n<=(NNs+NNnTS))
continue;
251 for (
int i=0;
i<nINP;
i++) params[
i]=0.;
253 for(
int n=s;
n<s+TS;
n++) {
254 if (uf!=0&&
n>=NNs&&
n<=(NNs+NNnTS))
continue;
257 NNcc=(double)signal.
netcc[1];
258 NNrho=(
double)signal.
rho[0];
259 for (
int i=0;
i<nINP;
i++){
263 double output=mlp->Evaluate(0,params);
265 if (out<0.6) sig_05=sig_05+1;
266 cout<<
"rho: "<<signal.
rho[0]<<
" cc "<<signal.
netcc[1]<<
" out: "<<out<<endl;
270 cout<<
"riempito sig"<<endl;
273 for(
int n=b;
n<b+
TB;
n++) {
274 if (uf!=0&&
n>=NNb&&
n<=(NNb+NNnTB))
continue;
276 cout<<
"n: "<<
n<<
"Bg index"<<(
n-sig_entries)<<endl;
278 NNcc=(double)background.
netcc[1];
279 NNrho=(
double)background.
rho[0];
280 for (
int i=0;
i<nINP;
i++){
285 double output=mlp->Evaluate(0,params);
287 if(out>0.6) bg_05=bg_05+1;
288 cout<<
"rho: "<<background.
rho[0]<<
" cc "<<background.
netcc[1]<<
" out: "<<out<<endl;
291 cout<<
"riempito bg"<<endl;
293 Realentries=NNTree2->GetEntries();
296 cout<<
"chiuso file"<<endl;
300 cout<<
"dopo richiamo funzione"<<endl;
303 cout<<
" Sig with out<06: "<<sig_05<<endl;
304 cout<<
" Bg with out>06: "<<bg_05<<endl;
305 cout<<
" Realentries "<<Realentries<<endl;
311 name.ReplaceAll(
"outfile/",
"");
312 TFile* fTEST =TFile::Open(ifile.Data());
313 TTree* NNTree2=(TTree*)fTEST->Get(
"Parameters");
318 NNTree2->SetBranchAddress(
"ANNout",&out);
319 NNTree2->SetBranchAddress(
"cc",&cc);
320 NNTree2->SetBranchAddress(
"rho",&rho);
321 NNTree2->SetBranchAddress(
"#TestSig",&nSi);
322 NNTree2->GetEntry(0);
324 cout<<
"nSig: "<<nSig<<
" nSi: "<<nSi<<endl;
326 cout<<
"dentro funzione dopodef"<<endl;
328 double* rhoSig[ncurve];
329 for (
int i=0;
i<ncurve;
i++) rhoSig[
i]=
new double[nSig];
330 cout<<
"dopo def rhoSig"<<endl;
334 int const nBg=NNTree2->GetEntries()-nSig;
335 for (
int i=0;
i<ncurve;
i++) {
337 for (
int j=0;
j<nSig;
j++) rhoSig[
i][
j]=0.;
339 cout<<
"dopo def rhoSig"<<endl;
340 double* rhoBg[ncurve];
341 for (
int i=0;
i<ncurve;
i++) rhoBg[
i]=
new double[nBg];
344 for (
int i=0;
i<ncurve;
i++) {
346 for (
int j=0;
j<nBg;
j++) rhoBg[
i][
j]=0.;
348 cout<<
"dopo def rhoBg"<<endl;
350 for (
int i=0;
i<
nCC;
i++) ccTh[
i]=0.;
352 for (
int i=0;
i<
nANN;
i++) NNTh[
i]=0.;
358 cout<<NNTree2->GetEntries()<<endl;
359 for(
int n=0;
n<NNTree2->GetEntries();
n++){
361 NNTree2->GetEntry(
n);
362 cout<<
"rho "<<rho<<
" cc "<<cc<<
" out "<<out<<endl;
365 if(cc<ccTh[
i])
continue;
367 if(
j==0) NNTh[
j]=-1000.;
373 if(out<NNTh[
j])
continue;
376 NBg[i*nANN+
j]= NBg[i*nANN+
j]+1;
377 while(rhoBg[i*nANN+j][ni]!=0)ni=ni+1;
378 rhoBg[i*nANN+
j][ni]=
rho;
383 NSig[i*nANN+
j]= NSig[i*nANN+
j]+1;
384 while(rhoSig[i*nANN+j][ni]!=0)ni=ni+1;
385 rhoSig[i*nANN+
j][ni]=
rho;
393 cout<<
"dopo riempimento variabili"<<endl;
395 for (
int i=0;
i<ncurve;
i++) indexS[
i]=
new int[nSig];
398 for (
int y=0;
y<ncurve;
y++) {
411 TMath::Sort(nSig,rhoSig[
y],indexS[y],
false);
413 for (
int k=0;
k<nSig;
k++) {
418 int ij=indexS[
y][
k-1];
420 if(rhoSig[y][ii]!=0) {
423 if(rhoSig[y][ii]!=rhoSig[y][ij]) gS[
y]->SetPoint(igS_p++,rhoSig[y][ii],yy);
424 cout<<
"igS"<<igS<<
" x "<<rhoSig[
y][ii]<<
" y: "<<yy<<endl;
429 if(rhoSig[y][ii]!=0){
431 if(rhoSig[y][ii]>100)
continue;
432 gS[
y]->SetPoint(0,rhoSig[y][ii],yy);
443 for (
int i=0;
i<ncurve;
i++) indexB[
i]=
new int[nBg];
446 for (
int y=0;
y<ncurve;
y++) {
450 TMath::Sort(nBg,rhoBg[
y],indexB[y],
false);
452 for (
int k=0;
k<nBg;
k++) {
456 int ij=indexB[
y][
k-1];
458 if(rhoBg[y][ii]!=0) {
461 if(rhoBg[y][ii]!=rhoBg[y][ij]) gB[
y]->SetPoint(igB_p++,rhoBg[y][ii],yy);
467 if(rhoBg[y][ii]!=0) {
469 gB[
y]->SetPoint(0,rhoBg[y][ii],yy);
476 cout<<
"dopo inserimento puntiB"<<endl;
479 TCanvas* cS=
new TCanvas(
"Efficiency_vs_rho",
"Efficiency_vs_rho",0,0,1200,700);
481 cS->cd(1)->SetLogy();
482 TMultiGraph* mg1=
new TMultiGraph();
484 gS[0]->SetMarkerColor(2);
485 gS[0]->SetLineColor(2);
486 mg1->SetTitle(
"cc=0.5;rho;#Events");
487 if(gS[1]->GetN()!=0) mg1->Add(gS[1]);
488 if(gS[0]->GetN()!=0) mg1->Add(gS[0]);
491 gS[
h]->SetMarkerColor(3);
492 gS[
h]->SetLineColor(3);
494 if(gS[
h]->GetN()!=0) mg1->Add(gS[
h]);
498 cS->cd(2)->SetLogy();
499 TMultiGraph* mg2=
new TMultiGraph();
501 gS[
nANN]->SetMarkerColor(2);
502 gS[
nANN]->SetLineColor(2);
503 mg2->SetTitle(
"cc=0.55;rho;#Events");
504 if(gS[nANN]->GetN()!=0) mg2->Add(gS[nANN]);
506 gS[nANN+
h]->SetMarkerColor(3);
507 gS[nANN+
h]->SetLineColor(3);
509 if(gS[nANN+
h]->GetN()!=0) mg2->Add(gS[nANN+
h]);
513 cS->cd(3)->SetLogy();
514 TMultiGraph* mg3=
new TMultiGraph();
516 gS[nANN*2]->SetMarkerColor(2);
517 gS[nANN*2]->SetLineColor(2);
518 mg3->SetTitle(
"cc=0.6;rho;#Events");
519 if(gS[nANN*2+1]->GetN()!=0) mg3->Add(gS[nANN*2+1]);
520 if(gS[nANN*2]->GetN()!=0) mg3->Add(gS[nANN*2]);
522 gS[2*nANN+
h]->SetMarkerColor(3);
523 gS[2*nANN+
h]->SetLineColor(3);
525 if(gS[2*nANN+
h]->GetN()!=0) mg3->Add(gS[2*nANN+
h]);
529 cS->cd(4)->SetLogy();
530 TMultiGraph* mg4=
new TMultiGraph();
532 gS[nANN*3]->SetMarkerColor(2);
533 gS[nANN*3]->SetLineColor(2);
534 mg4->SetTitle(
"cc=0.65;rho;#Events");
535 if(gS[nANN*3]->GetN()!=0) mg4->Add(gS[nANN*3]);
538 gS[3*nANN+
h]->SetMarkerColor(3);
539 gS[3*nANN+
h]->SetLineColor(3);
541 if(gS[3*nANN+
h]->GetN()!=0) mg4->Add(gS[3*nANN+
h]);
545 cout<<
"nuovo canv"<<endl;
546 TCanvas* cB=
new TCanvas(
"Number_vs_rho",
"Number_vs_rho",0,0,1200,700);
548 cB->cd(1)->SetLogy();
549 TMultiGraph* mg1B=
new TMultiGraph();
551 gB[0]->SetMarkerColor(2);
552 gB[0]->SetLineColor(2);
553 mg1B->SetTitle(
"cc=0.5;rho;#Events");
554 if(gB[0]->GetN()!=0) mg1B->Add(gB[0]);
556 gB[
h]->SetMarkerColor(3);
557 gB[
h]->SetLineColor(3);
559 if(gB[
h]->GetN()!=0) mg1B->Add(gB[
h]);
563 cB->cd(2)->SetLogy();
564 TMultiGraph* mg2B=
new TMultiGraph();
566 gB[
nANN]->SetMarkerColor(2);
567 gB[
nANN]->SetLineColor(2);
568 mg2B->SetTitle(
"cc=0.55;rho;#Events");
569 if(gB[nANN]->GetN()!=0) mg2B->Add(gB[nANN]);
571 gB[nANN+
h]->SetMarkerColor(3);
572 gB[nANN+
h]->SetLineColor(3);
574 if(gB[nANN+
h]->GetN()!=0) mg2B->Add(gB[nANN+
h]);
578 cB->cd(3)->SetLogy();
579 TMultiGraph* mg3B=
new TMultiGraph();
581 gB[2*
nANN]->SetMarkerColor(2);
582 gB[2*
nANN]->SetLineColor(2);
583 mg3B->SetTitle(
"cc=0.6;rho;#Events");
584 if(gB[2*nANN]->GetN()!=0) mg3B->Add(gB[2*nANN]);
586 gB[2*nANN+
h]->SetMarkerColor(3);
587 gB[2*nANN+
h]->SetLineColor(3);
589 if(gB[2*nANN+
h]->GetN()!=0) mg3B->Add(gB[2*nANN+
h]);
593 cB->cd(4)->SetLogy();
594 TMultiGraph* mg4B=
new TMultiGraph();
596 gB[3*
nANN]->SetMarkerColor(2);
597 gB[3*
nANN]->SetLineColor(2);
598 mg4B->SetTitle(
"cc=0.65;rho;#Events");
599 if(gB[3*nANN]->GetN()!=0) mg4B->Add(gB[3*nANN]);
601 gB[3*nANN+
h]->SetMarkerColor(3);
602 gB[3*nANN+
h]->SetLineColor(3);
604 if(gB[3*nANN+
h]->GetN()!=0) mg4B->Add(gB[3*nANN+
h]);
616 char CnameS2root[1024];
617 char CnameB2root[1024];
618 CnameS.ReplaceAll(
".root",
".png");
619 CnameB.ReplaceAll(
".root",
".png");
620 sprintf(CnameS2,
"logN_rho/logN_rho_S_dANN%1.2f_%s",
deltaANN,CnameS.Data());
621 sprintf(CnameB2,
"logN_rho/logN_rho_B_dANN%1.2f_%s",
deltaANN,CnameB.Data());
622 sprintf(CnameS2root,
"logN_rho/logN_rho_S_dANN%1.2f_%s",
deltaANN,CnameSroot.Data());
623 sprintf(CnameB2root,
"logN_rho/logN_rho_B_dANN%1.2f_%s",
deltaANN,CnameBroot.Data());
626 cS->Print(CnameS2root);
627 cB->Print(CnameB2root);
635 name.ReplaceAll(
"outfile/",
"");
636 TFile* fTEST =TFile::Open(ifile.Data());
637 TTree* NNTree2=(TTree*)fTEST->Get(
"Parameters");
642 NNTree2->SetBranchAddress(
"ANNout",&out);
643 NNTree2->SetBranchAddress(
"cc",&cc);
644 NNTree2->SetBranchAddress(
"rho",&rho);
645 NNTree2->SetBranchAddress(
"#TestSig",&nSi);
646 NNTree2->GetEntry(0);
652 double* ANNSig[ncurve2];
653 for (
int i=0;
i<ncurve2;
i++) ANNSig[
i]=
new double[nSig];
657 int const nBg=NNTree2->GetEntries()-nSig;
658 for (
int i=0;
i<ncurve2;
i++) {
660 for (
int j=0;
j<nSig;
j++) ANNSig[
i][
j]=0.;
662 double* ANNBg[ncurve2];
663 for (
int i=0;
i<ncurve2;
i++) ANNBg[
i]=
new double[nBg];
667 for (
int i=0;
i<ncurve2;
i++) {
669 for (
int j=0;
j<nBg;
j++) ANNBg[
i][
j]=0.;
672 for (
int i=0;
i<
nCC;
i++) ccTh[
i]=0.;
674 for (
int i=0;
i<
nRHO;
i++) rhoTh[
i]=0.;
681 for(
int n=0;
n<NNTree2->GetEntries();
n++){
682 NNTree2->GetEntry(
n);
683 cout<<
"rho "<<rho<<
" cc "<<cc<<
" out "<<out<<endl;
686 if(cc<ccTh[
i])
continue;
689 if(rho<rhoTh[
j])
continue;
692 NBg[i*nRHO+
j]= NBg[i*nRHO+
j]+1;
693 while(ANNBg[i*nRHO+j][ni]!=0)ni=ni+1;
694 ANNBg[i*nRHO+
j][ni]=
out;
698 NSig[i*nRHO+
j]= NSig[i*nRHO+
j]+1;
699 while(ANNSig[i*nRHO+j][ni]!=0)ni=ni+1;
700 ANNSig[i*nRHO+
j][ni]=
out;
707 int* indexS[ncurve2];
708 for (
int i=0;
i<ncurve2;
i++) indexS[
i]=
new int[nSig];
710 TGraph * gS[ncurve2];
711 for (
int y=0;
y<ncurve2;
y++) {
715 TMath::Sort(nSig,ANNSig[
y],indexS[y],
false);
716 for (
int k=0;
k<nSig;
k++) {
721 int ij=indexS[
y][
k-1];
723 if(ANNSig[y][ii]!=0) {
726 if(ANNSig[y][ii]!=ANNSig[y][ij]) gS[
y]->SetPoint(igS_p++,ANNSig[y][ii],yy);
733 if(ANNSig[y][ii]!=0){
735 gS[
y]->SetPoint(0,ANNSig[y][ii],yy);
746 int* indexB[ncurve2];
747 for (
int i=0;
i<ncurve2;
i++) indexB[
i]=
new int[nBg];
749 TGraph * gB[ncurve2];
750 for (
int y=0;
y<ncurve2;
y++) {
754 TMath::Sort(nBg,ANNBg[
y],indexB[y],
false);
756 for (
int k=0;
k<nBg;
k++) {
760 int ij=indexB[
y][
k-1];
762 if(ANNBg[y][ii]!=0) {
765 if(ANNBg[y][ii]!=ANNBg[y][ij]) gB[
y]->SetPoint(igB_p++,ANNBg[y][ii],yy);
771 if(ANNBg[y][ii]!=0) {
773 gB[
y]->SetPoint(0,ANNBg[y][ii],yy);
782 TCanvas* cS=
new TCanvas(
"Efficiency_vs_ANN",
"Efficiency_vs_ANN",0,0,1200,700);
786 TMultiGraph* mg1=
new TMultiGraph();
787 mg1->SetTitle(
"cc=0.5;ANN;#Events");
789 gS[
h]->SetLineColor(4);
790 if(gS[
h]->GetN()!=0) mg1->Add(gS[
h]);
795 TMultiGraph* mg2=
new TMultiGraph();
796 mg2->SetTitle(
"cc=0.55;ANN;#Events");
798 gS[nRHO+
h]->SetLineColor(4);
799 if(gS[nRHO+
h]->GetN()!=0) mg2->Add(gS[nRHO+
h]);
805 TMultiGraph* mg3=
new TMultiGraph();
806 mg3->SetTitle(
"cc=0.6;ANN;#Events");
808 gS[2*nRHO+
h]->SetLineColor(4);
809 if(gS[2*nRHO+
h]->GetN()!=0) mg3->Add(gS[2*nRHO+
h]);
814 TMultiGraph* mg4=
new TMultiGraph();
815 mg4->SetTitle(
"cc=0.65;ANN;#Events");
817 gS[3*nRHO+
h]->SetLineColor(4);
818 if(gS[3*nRHO+
h]->GetN()!=0) mg4->Add(gS[3*nRHO+
h]);
822 TCanvas* cB=
new TCanvas(
"Number_vs_ANN",
"Number_vs_ANN",0,0,1200,700);
824 cB->cd(1)->SetLogy();
825 TMultiGraph* mg1B=
new TMultiGraph();
826 mg1B->SetTitle(
"cc=0.5;ANN;#Events");
828 gB[
h]->SetLineColor(4);
829 if(gB[
h]->GetN()!=0) mg1B->Add(gB[
h]);
832 cB->cd(2)->SetLogy();
833 TMultiGraph* mg2B=
new TMultiGraph();
834 mg2B->SetTitle(
"cc=0.55;ANN;#Events");
836 gB[nRHO+
h]->SetLineColor(4);
837 if(gB[nRHO+
h]->GetN()!=0) mg2B->Add(gB[nRHO+
h]);
840 cB->cd(3)->SetLogy();
841 TMultiGraph* mg3B=
new TMultiGraph();
842 mg3B->SetTitle(
"cc=0.6;ANN;#Events");
844 gB[2*nRHO+
h]->SetLineColor(4);
845 if(gB[2*nRHO+
h]->GetN()!=0) mg3B->Add(gB[2*nRHO+
h]);
848 cB->cd(4)->SetLogy();
849 TMultiGraph* mg4B=
new TMultiGraph();
850 mg4B->SetTitle(
"cc=0.6;ANN;#Events");
852 gB[3*nRHO+
h]->SetLineColor(4);
853 if(gB[3*nRHO+
h]->GetN()!=0) mg4B->Add(gB[3*nRHO+
h]);
865 char CnameS2root[1024];
866 char CnameB2root[1024];
867 CnameS.ReplaceAll(
".root",
".png");
868 CnameB.ReplaceAll(
".root",
".png");
869 sprintf(CnameS2,
"ANNthres/N_ANN_S_%s",CnameS.Data());
870 sprintf(CnameB2,
"ANNthres/N_ANN_B_%s",CnameB.Data());
871 sprintf(CnameS2root,
"ANNthres/N_ANN_S_%s",CnameSroot.Data());
872 sprintf(CnameB2root,
"ANNthres/N_ANN_B_%s",CnameBroot.Data());
875 cS->Print(CnameS2root);
876 cB->Print(CnameB2root);
888 name.ReplaceAll(
"outfile/",
"");
889 TFile* fTEST =TFile::Open(ifile.Data());
890 TTree* NNTree2=(TTree*)fTEST->Get(
"Parameters");
895 NNTree2->SetBranchAddress(
"ANNout",&out);
896 NNTree2->SetBranchAddress(
"cc",&cc);
897 NNTree2->SetBranchAddress(
"rho",&rho);
898 NNTree2->SetBranchAddress(
"#TestSig",&nSi);
899 NNTree2->GetEntry(0);
901 cout<<
"nSig: "<<nSig<<
" nSi: "<<nSi<<endl;
902 int const nBg=NNTree2->GetEntries()-nSig;
931 for (
int i=0;
i<nSig;
i++){
939 for (
int i=0;
i<nBg;
i++){
945 for (
int n = 0;
n <NNTree2->GetEntries();
n++){
946 NNTree2->GetEntry(
n);
974 for (
int k=0;
k<nBg;
k++)indexB[
k]=0;
975 TMath::Sort(nBg,aRHOB,indexB,
false);
979 double ANN1th[
nth+1];
981 for (
int k=0;
k<(
nth);
k++){
990 for (
int k=0;
k<(
nth);
k++){
991 ANN1th[
k+1]=0.+dANNth*(
k+1);
992 cc1th[
k]=0.5+dccth*
k;
1007 for (
int a=0;
a<nBg;
a++){
1008 cout<<
"tutti: a "<<
a<<
" rho "<<aRHOB[
a]<<
" aNNB "<<aANNB[
a]<<
" cc "<<aCCB[
a]<<endl;
1009 if(aRHOB[
a]<
RHOth)
continue;
1010 cout<<
"dopo rho: a "<<
a<<
" rho "<<aRHOB[
a]<<
" aNNB "<<aANNB[
a]<<
" cc "<<aCCB[
a]<<endl;
1011 for (
int b=0;b<
nth;b++){
1012 if(aCCB[
a]<cc1th[b])
continue;
1013 cout<<
"dopo cc: "<<cc1th[b]<<
"a "<<
a<<
" rho "<<aRHOB[
a]<<
" aNNB "<<aANNB[
a]<<
" cc "<<aCCB[
a]<<endl;
1014 for(
int c=0;
c<nth+1;
c++){
1015 if(aANNB[
a]<ANN1th[
c])
continue;
1016 cout<<
"dopo ANN:"<<ANN1th[
c]<<
" a "<<
a<<
" rho "<<aRHOB[
a]<<
" aNNB "<<aANNB[
a]<<
" cc "<<aCCB[
a]<<endl;
1017 Z[b*(nth+1)+c]=Z[b*(nth+1)+
c]+1;
1018 cout<<
"Zindex: "<<b*(nth+1)+c<<
" Z: "<<Z[b*(nth+1)+
c]<<endl;
1029 for (
int d=0; d<
nth+1;d++) bins[d]=0.;
1031 for (
int d=0; d<nth+1;d++) bins[d+1]=0.+(d+1)*dANNth;
1032 TH2D* hglitch=
new TH2D(
"Survived_glitches(RHO_cut:5.5)",
"Survived_glitches(RHO_cut:5.5)",nth,0.5,1,nth+1,bins);
1034 hglitch->SetDirectory(0);
1036 for (
int b=0;b<
nth;b++)
for(
int c=0;
c<nth+1;
c++) {
1038 hglitch->Fill(cc1th[b]+dccth/2.,bins[
c]+dANNth/2.,Z[(nth+1)*b+
c]);
1039 cout<<
" Zindex: "<<(nth+1)*b+
c<<
" x: "<<cc1th[b]<<
" y: "<<ANN1th[
c]<<
" z: "<<Z[(nth+1)*b+
c]<<endl;
1042 gB[0]=
new TGraph(nBg,aCCB,aRHOB);
1043 gS[0]=
new TGraph(nSig,aCCS,aRHOS);
1044 gB[1]=
new TGraph(nBg,aCCB,aANNB);
1045 gS[1]=
new TGraph(nSig,aCCS,aANNS);
1046 gB[2]=
new TGraph(nBg,aRHOB,aANNB);
1047 gS[2]=
new TGraph(nSig,aRHOS,aANNS);
1049 for(
int i=0;
i<3;
i++){
1050 gS[
i]->SetMarkerColor(2);
1051 gB[
i]->SetMarkerColor(4);
1052 gS[
i]->SetMarkerStyle(6);
1053 gB[
i]->SetMarkerStyle(7);
1056 TCanvas*
c=
new TCanvas(
"Plots",
"Plots",0,0,1200,700);
1059 TMultiGraph* mg1=
new TMultiGraph();
1060 mg1->SetTitle(
"cc_rho;cc;rho");
1061 if(gB[0]->GetN()!=0) mg1->Add(gB[0]);
1062 if(gS[0]->GetN()!=0) mg1->Add(gS[0]);
1065 TMultiGraph* mg2=
new TMultiGraph();
1066 mg2->SetTitle(
"cc_ANNoutput;cc;ANNoutput");
1067 if(gB[1]->GetN()!=0) mg2->Add(gB[1]);
1068 if(gS[1]->GetN()!=0) mg2->Add(gS[1]);
1071 TMultiGraph* mg3=
new TMultiGraph();
1072 mg3->SetTitle(
"rho_ANNoutput;rho;ANNoutput");
1073 if(gB[2]->GetN()!=0) mg3->Add(gB[2]);
1074 if(gS[2]->GetN()!=0) mg3->Add(gS[2]);
1077 TText*
text=
new TText(0.37,0.0,
"no cuts on ANN");
1078 hglitch->SetStats(0);
1079 hglitch->GetXaxis()->SetTitle(
"cc");
1080 hglitch->GetYaxis()->SetTitle(
"ANNoutput");
1081 hglitch->GetZaxis()->SetTitle(
"count");
1082 hglitch->Draw(
"colz");
1089 char Cname2root[1024];
1090 Cname.ReplaceAll(
".root",
".png");
1091 sprintf(Cname2,
"CC_RHO_ANN_Plots/CC_RHO_ANN_Plots_%s",Cname.Data());
1092 sprintf(Cname2root,
"CC_RHO_ANN_Plots/CC_RHO_ANN_Plots_%s",Cnameroot.Data());
1094 c->Print(Cname2root);
1096 TCanvas*
c2=
new TCanvas(
"Plots_Bkg_on_Sig",
"Plots_Bkg_on_Sig",0,0,1200,700);
1100 gB[0]->Draw(
"p,same");
1104 gB[1]->Draw(
"p,same");
1107 gB[2]->Draw(
"p,same");
1109 TText* text2=
new TText(0.37,0.0,
"no cuts on ANN");
1110 hglitch->SetStats(0);
1111 hglitch->GetXaxis()->SetTitle(
"cc");
1112 hglitch->GetYaxis()->SetTitle(
"ANNoutput");
1113 hglitch->GetZaxis()->SetTitle(
"count");
1114 hglitch->Draw(
"colz");
1120 char Cname2_2[1024];
1121 char Cname2_2root[1024];
1122 Cname_2.ReplaceAll(
".root",
".png");
1123 sprintf(Cname2_2,
"CC_RHO_ANN_Plots/CC_RHO_ANN_Plots_BoS_%s",Cname_2.Data());
1124 sprintf(Cname2_2root,
"CC_RHO_ANN_Plots/CC_RHO_ANN_Plots_BoS_%s",Cname_2root.Data());
1125 c2->SaveAs(Cname2_2);
1126 c2->Print(Cname2_2root);
void graph(TString ifile)
Float_t * rho
biased null statistics
wavearray< double > a(hp.size())
cout<< endl;cout<< "ts size = "<< ts.size()<< " ts rate = "<< ts.rate()<< endl;tf.Forward(ts, wdm);int levels=tf.getLevel();cout<< "tf size = "<< tf.size()<< endl;double dF=tf.resolution();double dT=1./(2 *dF);cout<< "rate(hz) : "<< RATE<< "\ layers : "<< nLAYERS<< "\ dF(hz) : "<< dF<< "\ dT(ms) : "<< dT *1000.<< endl;int itime=TIME_PIXEL_INDEX;int ifreq=FREQ_PIXEL_INDEX;int index=(levels+1) *itime+ifreq;double time=itime *dT;double freq=(ifreq >0) ? ifreq *dF :dF/4;cout<< endl;cout<< "PIXEL TIME = "<< time<< " sec "<< endl;cout<< "PIXEL FREQ = "<< freq<< " Hz "<< endl;cout<< endl;wavearray< double > x
cout<< fileout<< endl;out.open(fileout, ios::out);if(!out.good()) { cout<< "cwb_mkhtml_cbc.C : Error Opening File : "<< fileout<< endl;exit(1);} out<< "<br><br><< endl; MakePlotsHtmlCellTable(&out, "Sensitive distance for M< sub > total</sub > bins
void Plots(TString ifile)
cout<< "SNR "<< xsnr<< endl;wavearray< double > f
cout<< "Injected signals: "<< mdc.GetEntries()<< endl;cout<< "Injected signals in histogram factor_events_inj: "<< NEVTS<< endl;float myifar, ecor, m1, m2, netcc[3], neted, penalty;float rho[2];float chirp[6];float range[2];float frequency[2];float iSNR[3], sSNR[3];sim.SetBranchAddress("mass", mass);sim.SetBranchAddress("factor", &factor);sim.SetBranchAddress("range", range);sim.SetBranchAddress("chirp", chirp);sim.SetBranchAddress("rho", rho);sim.SetBranchAddress("netcc", netcc);sim.SetBranchAddress("neted", &neted);sim.SetBranchAddress("ifar", &myifar);sim.SetBranchAddress("ecor", &ecor);sim.SetBranchAddress("penalty", &penalty);sim.SetBranchAddress("time", mytime);sim.SetBranchAddress("iSNR", iSNR);sim.SetBranchAddress("sSNR", sSNR);sim.SetBranchAddress("spin", spin);sim.SetBranchAddress("frequency", frequency);float **volume=new float *[NBINS_mass1];float **volume_first_shell=new float *[NBINS_mass1];float **radius=new float *[NBINS_mass1];float **error_volume=new float *[NBINS_mass1];float **error_volume_first_shell=new float *[NBINS_mass1];float **error_radius=new float *[NBINS_mass1];for(int i=0;i< NBINS_mass1;i++) { volume[i]=new float[NBINS_mass2];volume_first_shell[i]=new float[NBINS_mass2];radius[i]=new float[NBINS_mass2];error_volume[i]=new float[NBINS_mass2];error_volume_first_shell[i]=new float[NBINS_mass2];error_radius[i]=new float[NBINS_mass2];for(int j=0;j< NBINS_mass2;j++) { volume[i][j]=0.;volume_first_shell[i][j]=0.;radius[i][j]=0.;error_volume[i][j]=0.;error_volume_first_shell[i][j]=0.;error_radius[i][j]=0.;} } float **spin_mtot_volume=new float *[NBINS_MTOT+1];float **spin_mtot_radius=new float *[NBINS_MTOT+1];float **error_spin_mtot_volume=new float *[NBINS_MTOT+1];float **error_spin_mtot_radius=new float *[NBINS_MTOT+1];for(int i=0;i< NBINS_MTOT+1;i++) { spin_mtot_volume[i]=new float[NBINS_SPIN+1];spin_mtot_radius[i]=new float[NBINS_SPIN+1];error_spin_mtot_volume[i]=new float[NBINS_SPIN+1];error_spin_mtot_radius[i]=new float[NBINS_SPIN+1];for(int j=0;j< NBINS_SPIN+1;j++) { spin_mtot_volume[i][j]=0.;error_spin_mtot_volume[i][j]=0.;spin_mtot_radius[i][j]=0.;error_spin_mtot_radius[i][j]=0.;} } char fname[1024];sprintf(fname, "%s/recovered_signals.txt", netdir);ofstream fev;fev.open(fname, std::ofstream::out);sprintf(line, "#GPS@L1 FAR[Hz] eFAR[Hz] Pval " "ePval factor rho frequency iSNR sSNR \");fev<< line<< endl;ofstream *fev_single=new ofstream[nfactor];for(int l=1;l< nfactor+1;l++) { sprintf(fname, "%s/recovered_signals_%d.txt", netdir, l);fev_single[l - 1].open(fname, std::ofstream::out);fev_single[l - 1]<< line<< endl;} double Vrho[RHO_NBINS], eVrho[RHO_NBINS], Rrho[RHO_NBINS], eRrho[RHO_NBINS], Trho[RHO_NBINS];for(int i=0;i< RHO_NBINS;i++) { Vrho[i]=0.;eVrho[i]=0.;Rrho[i]=0.;eRrho[i]=0.;Trho[i]=RHO_MIN+i *RHO_BIN;} double dV, dV1, dV_spin_mtot, nevts, internal_volume;int nT;int countv=0;int cnt=0;int cnt2=0;int cntfreq=0;bool bcut=false;double liveTot=sim.GetMaximum("ifar");double BKG_LIVETIME_yr=liveTot/CYS;double BKG_LIVETIME_Myr=BKG_LIVETIME_yr/(1.e6);cout.precision(14);cout<< "Total live time ---> background
Float_t * netcc
effective correlated SNR
void Annth(TString ifile)
sprintf(tfres,"(1/%g)x(%g) (sec)x(Hz)", 2 *df, df)
void Test0_dANN_var(TString NN_FILE, TString TEST_FILE, int TS=0, int TB=0, int s=0, int b=0, int uf=0)