Logo Coherent WaveBurst  
Reference Guide 		Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JoinCutANN_Mchirp_ROCcurves_rho2.C
Go to the documentation of this file.
1 #define XIFO 4
2 
3 #pragma GCC system_header
4 
5 #include "cwb.hh"
6 #include "config.hh"
7 #include "network.hh"
8 #include "wavearray.hh"
9 #include "TString.h"
10 #include "TObjArray.h"
11 #include "TObjString.h"
12 #include "TPaletteAxis.h"
13 #include "TMultiLayerPerceptron.h"
14 #include "TMLPAnalyzer.h"
15 #include "TRandom.h"
16 #include "TComplex.h"
17 #include "TMath.h"
18 #include "mdc.hh"
19 #include "watplot.hh"
20 #include <vector>
21 
22 
23 //#define nINP 64
24 #define Ncc3D 50
25 #define maxCC 1
26 #define minCC 0
27 #define nIFO 3
28 #define nRHO 3
29 #define nANN 10
30 #define deltaANN 0.1
31 #define deltacc 0.05
32 #define deltarho 1
33 #define nCC 4
34 #define NPIX 20
35 #define RHOth 5.5
36 #define nth 10
37 #define rhomin 6
38 #define ccmin 0.6
39 #define iMc 1
40 #define NMc_inj 50
41 #define minMc_inj 0
42 #define maxMc_inj 100
43 #define NMc_r 50
44 #define maxMc_r 40
45 #define minMc_r -5
46 #define NANN_av 50
47 #define minANN_av -0.2
48 #define maxANN_av 1.5
49 
50 //definition parameters for SNR_plots
51 #define N_ANNth 50
52 #define N_RHOth 50
53 #define ANNth_max 1.5
54 #define ANNth_min -0.5
55 #define RHOth_max 10.0
56 #define RHOth_min 4.5
57 
58 #define x_asym 0.00009999999999999999999999999999
59 #define nANNth 6
60 #define min_yasym -0.1
61 #define max_yasym 0.6
62 #define max_Mcth 2.5
63 #define min_Mcth 1.5
64 #define npoints 50
65 #define max_k 10
66 #define min_k 0.01
67 
68 #define cutMc1 0.5
69 #define cutMc2 1.5
70 #define cutANN1 0.15
71 #define cutANN2 0.25
72 #define maxrho 20.
73 #define minrho 0.
74 
75 
76 using namespace std;
77 void graph(TString ifile);
78 void Annth(TString ifile, int &FAth_n, int &FDth_n, double &thresFA);
79 void Mcth(TString ifile, int &McFAth_n, int &McFDth_n, double &McthresFA);
80 void PlotsAv_cc(TString ifile);
81 void PlotsAv_Mc(TString ifile);
82 void SNR_plots(TString ifile);
83 void CoG_plots(TString ifile);
84 void CutMcANN(TString ifile);
85 //void Test0(TString NN_FILE,TString TEST_FILE,TString ofile, int TS=0, int TB=0, int s=0, int b=0, int uf=0){
86 
87 //NN_FILE-NN_FILE8: name of the files where the interested ANNs are saved (!!they must become from a directory whose name ends in "NN").
88 //TEST_FILE: name of the file which contains the events you are going to analyse (!!they must become from a directory whose name ends in "nnTREE").
89 //TS, TB: numbers of events to test for each class.
90 //s, b: index values of the first events considered for both the classes. If b < #SIG_events, b is newly set to b=b+#SIG_events
91 //uf !=0 means the Test file of events is the same used to train the network, in this case ni==0. It means that the events used for the training are analysed only by the other ANNs, in this way we test the events with at least NNn-1 networks. If uf==0 the Test and Train files are different, in this case ni!=0 and all the events can beanalysed by all the nNN networks.
92 //consider_all=1(0): the events not used for the training procedure are (aren't) considered for the test.
93 //ni: =0 excludes from the test the training set of events
94 void JoinCutANN_Mchirp_ROCcurves_rho2(TString NN_FILE,TString NN_FILE2, TString NN_FILE3, TString NN_FILE4, TString NN_FILE5, TString NN_FILE6, TString NN_FILE7, TString NN_FILE8, TString TEST_FILE, TString NOMEtot_S, int TS=0, int TB=0, int s=0, int b=0, int uf=1, int consider_all=0, int av_on_nNN=0){
95 int ni=0;
96 double FAdes=0.01;
97 if (uf==0&&av_on_nNN!=0) cout<<"all events were not used for the training procedure"<<endl;
98 if(uf==0) {ni=1; consider_all=1; av_on_nNN=0;}
99 else ni=0; //events used for the training not considered for the ANN trained on them
100 
101 int TS0=0;
102 TS0=TS;
103 int TB0=0;
104 TB0=TB;
105 //if(uf!=0&&consider_all==0) ni=1;
106 
107 //COUNT THE ANNs
108 int n_NN=0;
109 if (NN_FILE.CompareTo("")){
110  n_NN=n_NN+1;
111  }
112 if (NN_FILE2.CompareTo("")){
113  n_NN=n_NN+1;
114  }
115 if (NN_FILE3.CompareTo("")){
116  n_NN=n_NN+1;
117  }
118 if (NN_FILE4.CompareTo("")){
119  n_NN=n_NN+1;
120  }
121 if (NN_FILE5.CompareTo("")){
122  n_NN=n_NN+1;
123  }
124 if (NN_FILE6.CompareTo("")){
125  n_NN=n_NN+1;
126  }
127 if (NN_FILE7.CompareTo("")){
128  n_NN=n_NN+1;
129  }
130 if (NN_FILE8.CompareTo("")){
131  n_NN=n_NN+1;
132  }
133 int const nNN=n_NN;
134 
135 // int p=0;
136  char NNi[nNN][1024];
137  char NNi2[nNN][1024];
138  sprintf(NNi2[0],"%s",NN_FILE.Data());
139  if(nNN>1) sprintf(NNi2[1],"%s",NN_FILE2.Data());
140  if(nNN>2) sprintf(NNi2[2],"%s",NN_FILE3.Data());
141  if(nNN>3) sprintf(NNi2[3],"%s",NN_FILE4.Data());
142  if(nNN>4) sprintf(NNi2[4],"%s",NN_FILE5.Data());
143  if(nNN>5) sprintf(NNi2[5],"%s",NN_FILE6.Data());
144  if(nNN>6) sprintf(NNi2[6],"%s",NN_FILE7.Data());
145  if(nNN>7) sprintf(NNi2[7],"%s",NN_FILE8.Data());
146 
147 if(nNN==1) {av_on_nNN=1; consider_all=0; ni=0;}
148 //saving the ANN names
149  for (int u=0;u<nNN;u++){
150  int p=0;
151  while (NNi2[u][p]){
152  if (NNi2[u][p]=='N'){
153  if((NNi2[u][p+1]=='N')&&(NNi2[u][p+2]=='\/')) {
154  int hh=p+3;
155  while (NNi2[u][hh]!='\0'){NNi[u][hh-p-3]=NNi2[u][hh];hh=hh+1;}
156  break;
157  }
158  }
159  p=p+1;
160  }
161 }
162 
163 //saving the Test-file name
164 int q=0;
165  char Filei[1024];
166  char Filei2[1024];
167  sprintf(Filei2,"%s",TEST_FILE.Data());
168  while (Filei2[q]){
169  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]=='\/')) {
170  int hh=q+7;
171  while (Filei2[hh]!='\0') {Filei[hh-q-7]=Filei2[hh];hh=hh+1;
172  }
173  for (int h0=hh-q-7;h0<1024;h0++) Filei[h0]='\0';
174  break;
175  }
176  q=q+1;
177  }
178 
179  cout<<Filei<<" original String: "<<Filei2<<endl;
180 
181 
182  TString NNi0[nNN];
183  NNi0[0]=NNi[0];
184  if(nNN>1) NNi0[1]=NNi[1];
185  if(nNN>2) NNi0[2]=NNi[2];
186  if(nNN>3) NNi0[3]=NNi[3];
187  if(nNN>4) NNi0[4]=NNi[4];
188  if(nNN>5) NNi0[5]=NNi[5];
189  if(nNN>6) NNi0[6]=NNi[6];
190  if(nNN>7) NNi0[7]=NNi[7];
191 
192  TString Filei0(Filei);
193 
194 //removing the extensions
195 for (int u=0;u<nNN;u++){
196  NNi0[u].ReplaceAll(".root","");
197  }
198  Filei0.ReplaceAll(".root","");
199 
200 //Definitions of some parameters
201  TFile *fnet[nNN];
202  TMultiLayerPerceptron* mlp[nNN];
203  TTree* infot[nNN];
204  int TotBg[nNN];
205  int FD0[nNN]; //FD=False Dismissal
206  int FA0[nNN]; //FA=False Alarm
207 
208  for (int yy=0; yy<nNN; yy++){TotBg[yy]=0;FA0[yy]=0;FD0[yy]=0;}
209 
210  int NNs[nNN]; //first SIG-event considered for ANN trainings
211  int NNb[nNN]; //first BKG-event considered for ANN trainings
212  int NNnTS[nNN]; //number of SIG-events considered for ANN trainings
213  int NNnTB[nNN]; //number of BKG-events considered for ANN trainings
214 
215 int min_NNb=5000000;
216 int min_NNs=5000000;
217 int max_NNb=0;
218 int max_NNs=0;
219 
220 //Extracting ANNs and information from the files
221 //if (uf!=0&&ni==0&&(TS!=0||TB!=0)){
222  for (int u=0;u<nNN;u++){
223  fnet[u]=TFile::Open(NNi2[u]);
224  mlp[u]=(TMultiLayerPerceptron*)fnet[u]->Get("TMultiLayerPerceptron");
225 
226  cout<<"dopo TMLP"<<endl;
227  if(mlp[u]==NULL) {cout << "Error getting mlp" << endl;exit(1);}
228  infot[u]=(TTree*)fnet[u]->Get("info");
229 
230  infot[u]->SetBranchAddress("Rand_start_Sig",&NNs[u]);
231  infot[u]->SetBranchAddress("Rand_start_Bg",&NNb[u]);
232  infot[u]->SetBranchAddress("#trainSig",&NNnTS[u]);
233  infot[u]->SetBranchAddress("#trainBg",&NNnTB[u]);
234  infot[u]->GetEntry(0);
235  cout<<"n "<<u<<" b "<<NNb[u]<<" min_NNb "<<min_NNb<<endl;
236 
237  if(NNs[u]<min_NNs) min_NNs=NNs[u];
238  if(NNb[u]<min_NNb) min_NNb=NNb[u];
239  if(NNs[u]>max_NNs) max_NNs=NNs[u];
240  if(NNb[u]>max_NNb) max_NNb=NNb[u];
241  cout<<"n "<<u<<" s "<<NNs[u]<<" b "<<NNb[u]<<" s fin "<<NNs[u]+NNnTS[u]<<" b fin "<<NNb[u]+NNnTB[u]<<" nTS "<<NNnTS[u]<<" nTB "<<NNnTB[u]<<endl;
242  }
243 
244 //Changing b value to its sum to the first event considered for the training
245 //if (uf!=0&&ni==0&&(b<min_NNb)&&(TS!=0||TB!=0)){ b+=min_NNb;cout<<" change in b value to have BKG events: "<<b<<endl;}
246 /*if (uf!=0&&ni==0&&(b<min_NNb)){ b+=min_NNb;cout<<" change in b value to have BKG events: "<<b<<endl;}
247 if (uf!=0&&ni==0&&(s<min_NNs)){ s+=min_NNs;cout<<" change in s value to have SIG events: "<<s<<endl;}*/
248 //if (uf!=0&&consider_all==0&&av_on_nNN==0&&(b<min_NNb)){ b+=min_NNb;cout<<" change in b value to have BKG events: "<<b<<endl;}
249 //if (uf!=0&&consider_all==0&&av_on_nNN==0&&(s<min_NNs)){ s+=min_NNs;cout<<" change in s value to have SIG events: "<<s<<endl;}
250 
251 //if (uf!=0&&av_on_nNN!=0&&(b>min_NNb)&&(b<))
252 cout<<" min_NNb "<<min_NNb<<" b "<<b<<endl;
253 cout<<"Def. tree e mlp"<<endl;
254 
255 //opening the Test-file
256  TFile* fTEST =TFile::Open(TEST_FILE.Data());
257 cout<<"Def. tree e mlp"<<endl;
258  TTree* NNTree=(TTree*)fTEST->Get("nnTree");
259 cout<<"Def. tree e mlp"<<endl;
260  int entries=NNTree->GetEntries();
261 // int entries=0;
262 
263  cout<<"entries: "<<entries<<endl;
264  cout<<" NNTree->GetEntries() "<<NNTree->GetEntries()<<endl;
265 //exit(0);
266 //Extract information necessary for the tree which defines the mlp
267  int ndim;
268  int ninp;
269  int y;
270  int entriesTot;
271  int bg_entries;
272  int sig_entries;
273 
274  NNTree->SetBranchAddress("#Entries_type",&entriesTot);
275  NNTree->SetBranchAddress("Matrix_dim",&ndim);
276  NNTree->SetBranchAddress("#inputs",&ninp);
277  NNTree->SetBranchAddress("amplitude_mode",&y);
278 
279  NNTree->GetEntry(0);
280  sig_entries=entriesTot;
281  NNTree->GetEntry(entries-1);
282  bg_entries=entriesTot;
283  cout<<"sig "<<sig_entries<<" bg "<<bg_entries<<endl;
284  int const NDIM=ndim;
285  int const nINP=ninp;
286  cout<<"NDIM: "<<NDIM<<endl;
287  cout<<"nINP: "<<nINP<<endl;
288  cout<<"sig e: "<<sig_entries<<endl;
289  cout<<"bg e: "<<bg_entries<<endl;
290 
291 //defining the sample which has less events
292  int minevents=0;
293  if (sig_entries>bg_entries) minevents=bg_entries;
294  else minevents=sig_entries;
295 
296  //if(b==0) b=sig_entries;
297  if(TB==0 && TS==0){
298  TS=minevents;
299  TB=minevents;
300  }
301 cout<<"TS "<<TS<<" TB "<<TB<<" minevents "<<minevents<<endl;
302 //defining compatible parameters
303  if (b<sig_entries){
304  int a =b;
305  b+=sig_entries;
306  cout<<"Error: Bg index<sig_entries: new set of b parameter-> b="<<b<<" instead of b="<<a<<endl;
307  //exit(0);
308 }
309  if(s>sig_entries){
310  int a=s;
311  s-=sig_entries;
312  cout<<"Error: Sig index>sig_entries: new set of s parameter-> s="<<b<<" instead of s="<<a<<endl;
313 }
314 if((TS>sig_entries-s||TB>(bg_entries-(b-sig_entries)))&&(TS==TB)) {TS=minevents-s;TB=minevents-(b-sig_entries);
315  if (TS<TB) TB=TS;
316  else TS=TB;
317  cout<<"Error:S>sig_entries or TB>bg_entries, to maintain ugual number of analysed events TS=TB="<<TB<<endl;
318  }
319  if((TS>sig_entries-s||TB>bg_entries-(b-sig_entries))&&(TS!=TB)) {
320  if(TS>sig_entries) TS=sig_entries-s;
321  if (TB>bg_entries) TB=bg_entries-(b-sig_entries);
322  cout<<"Error:S>sig_entries or TB>bg_entries, new TS and TB values are thus define-> TS="<<TS<<" TB="<<TB<<endl;
323  }
324 
325 
326  char NOMEtot[1024];
327  sprintf(NOMEtot,"%s",NOMEtot_S.Data());
328  cout<<"output name: "<<NOMEtot<<endl;
329 
330 //extracting original files
331  TString SIG_FILE;
332  TString BG_FILE;
333  char FILE_NAME[516];
334  NNTree->SetBranchAddress("Files_name",&FILE_NAME);
335  NNTree->GetEntry(0);
336  SIG_FILE=FILE_NAME;
337  NNTree->GetEntry(entries-1);
338  BG_FILE=FILE_NAME;
339  cout<<"fine ifdef RHO_CC"<<endl;
340 
341 //extracting other information
342  TChain sigTree("waveburst");//search the Tree "waveburst" in files
343  sigTree.Add(SIG_FILE.Data());
344  netevent signal(&sigTree,nIFO);
345  int sig_entries2 = signal.GetEntries();
346  cout << "sig entries2 : " << sig_entries2 << endl;
347 
348  TChain bgTree("waveburst");
349  bgTree.Add(BG_FILE.Data());
350  netevent background(&bgTree,nIFO);
351  int bg_entries2 = background.GetEntries();
352  cout << "bg entries2 : " << bg_entries2 << endl;
353 
354  cout<<"b: "<<b<<endl;
355  cout<<"s: "<<s<<endl;
356 
357  // add leaf
358  Float_t x[nINP];
359  for (int jj=0; jj<nINP;jj++) x[jj]=0.;
360  char ilabel[nINP][16];
361 
362  // define a branch for each input
363  for(int i=0;i<nINP;i++) {
364  sprintf(ilabel[i],"x%i",i+1);
365  NNTree->SetBranchAddress(ilabel[i], &x[i]);
366  }
367 
368 char ofile[1024];
369 sprintf(ofile,"average_file/%s.root",NOMEtot);
370 TFile*f =new TFile(ofile,"RECREATE");
371  TTree* NNTree2=new TTree("Parameters","Parameters");
372  NNTree2->SetDirectory(f);
373 
374 // defining useful information
375  double average=0.;
376  double out[nNN];
377  for(int y=0; y<nNN; y++) out[y]=0.;
378  double Mc=0.;
379  double NNcc=0.;
380  double NNrho=0.;
381  double Mc_inj=0.;
382  double std=0.;
383  int ev_ind=0;
384  int CoG=0;
385  // int ev_ind;
386  int index_ev=0;
387  char TestFile[1024];
388 //adding new information to the original tree
389  NNTree2->Branch("Center_of_Gravity",&CoG,"Center_of_Gravity/I");
390  NNTree2->Branch("Average",&average,"Average/D");
391  NNTree2->Branch("StandardDevaition",&std,"StandardDeviation/D");
392  NNTree2->Branch("cc",&NNcc,"cc/D");
393  NNTree2->Branch("Mchirp",&Mc,"Mchirp/D");
394  NNTree2->Branch("Mchirp_injected",&Mc_inj,"Mchirp_injected/D");
395  NNTree2->Branch("rho",&NNrho,"rho/D");
396  NNTree2->Branch("index",&index_ev,"index/I");
397 
398 // NNTree2->Branch("event_index",&ev_ind,"event_index/I");
399  NNTree2->Branch("Test_file",&TestFile,"Test_file/C");
400 
401  for(int u=0;u<nNN;u++){
402  char NNoutl[512];
403  char NNoutl2[512];
404  char NNnamel[512];
405  char NNnamel2[512];
406  sprintf(NNoutl,"NNout%i",u);
407  sprintf(NNoutl2,"NNout%i/D",u);
408  sprintf(NNnamel,"NNname%i",u);
409  sprintf(NNnamel2,"NNname%i/C",u);
410  NNTree2->Branch(NNoutl,&out[u],NNoutl2);
411  NNTree2->Branch(NNnamel,&NNi[u],NNnamel2);
412  }
413 
414  char Testf[1024];
415  NNTree2->Branch("TestFile",&Testf,"TestFile/C");
416  sprintf(Testf,"%s",TEST_FILE.Data());
417  int nTestS=0;
418 
419 
420 //defining how many examples to use
421  int scount=0;
422  if(uf==0) nTestS=TS;
423  else {
424  //for(int n=s;n<s+TS;n++) {
425  for(int n=s;n<sig_entries;n++) {
426  int countNN=0;
427  for(int u=0;u<nNN;u++){
428  if (uf!=0&&n>NNs[u]&&n<(NNs[u]+NNnTS[u])) countNN=countNN+1;
429  }
430 
431  // if(consider_all==0&&av_on_nNN==0&&countNN==0)continue; //skip if no network is trained on the considered event
432  // if(consider_all==0&&av_on_nNN!=0&&countNN!=0) {sigskip=sigskip+1; continue; }
433 
434  if(countNN==0&&av_on_nNN!=0) scount=scount+1;
435  else{
436  //if(nNN==1&&countNN==1&&av_on_nNN!=0) scount=scount+1;
437  if(countNN==1&&consider_all==0&&av_on_nNN==0) scount=scount+1;
438  if(countNN<2&&consider_all!=0) scount=scount+1; //if you are indifferently intereseted in average on nNN or nNN-1
439  if(countNN>1){cout<<"Error: training non independent"<<n<<endl; exit(0);} //exit if the ANNs are training on dependent set of events
440  cout<<"n "<<n<<"; countNN "<<countNN<<"; scount "<< scount<<"; consider_all "<<consider_all<<"; av_on_nNN "<<av_on_nNN<<endl;
441  if(scount==TS) break;}
442  }
443  nTestS=scount;
444 }
445 cout<<" nTestS "<< nTestS<<" TS "<<TS<<endl;
446 //exit(0);
447 int B_eff=0;
448 int FA=0;
449 int bg_05=0;
450 int cont_su=0;
451 int cont_su5=0;
452 int cont_su4=0;
453 int cont_su3=0;
454 
455 int nTestB=0;
456  int bcount=0;
457  if(uf==0) nTestB=TB;
458  else {
459  //for(int n=b;n<b+TB;n++) {
460  for(int n=b;n<sig_entries+bg_entries;n++) {
461  int countNN=0;
462  for(int u=0;u<nNN;u++){
463  if (uf!=0&&n>NNb[u]&&n<(NNb[u]+NNnTB[u])) countNN=countNN+1;
464  }
465  if(countNN==0&&av_on_nNN!=0) bcount=bcount+1;
466  //if(nNN==1&&countNN==1&&av_on_nNN!=0) bcount=bcount+1;
467  if(countNN==1&&consider_all==0&&av_on_nNN==0) bcount=bcount+1;
468  if(countNN<2&&consider_all!=0) bcount=bcount+1; //if you are indifferently intereseted in average on nNN or nNN-1
469  if(countNN>1){cout<<"Error: training non independent"<<n<<endl; exit(0);} //exit if the ANNs are training on dependent set of events
470  if(bcount==TB) break;
471  }
472  nTestB=bcount;
473 }
474 
475 if(uf!=0) nTestS=nTestS-1;
476 if(uf!=0) nTestB=nTestB-1;
477 //TB=TB-1;
478 //TS=TS-1;
479 cout<<"TS "<<TS<<" TB "<<TB<<" nTestS "<<nTestS<<" nTestB "<<nTestB<<endl;
480 if(TS0==0&&TB0==0) {
481  if(nTestB>nTestS) {nTestB=nTestS; TB=nTestS; TS=nTestS;}
482  else {nTestS=nTestB; TB=nTestB; TS=nTestB;}
483  }
484 
485 if(nTestS<TS) TS= nTestS;
486 else nTestS=TS;
487 
488 // NNTree2->Branch("#TestSig",&TS,"#TestSig/I");
489  NNTree2->Branch("#TestSig",&TS,"#TestSig/I");
490  //NNTree2->Branch("#TestSig",&nTestS,"#TestSig/I");
491  cout<<"nTestS: "<<nTestS<<" TS: "<<TS<<endl;
492  cout<<"dopo def tree"<<endl;
493 cout<<"s test "<<s<< " s+TS "<<s+TS<<" b "<<b<<" b+ TB "<<b+TB<<" nTestS "<<nTestS<<" TB "<<TB<<" nTestB "<<nTestB<<endl;
494 double params[nINP];
495  int sig_05=0;
496  for (int i=0; i<nINP; i++) params[i]=0.;
497 
498  int S_eff=0;
499  int FD=0;
500 int sigskip=0;
501 
502  TString txt_originaldata(NOMEtot_S);
503  txt_originaldata.ReplaceAll(".root","_originalData.txt");
504  //txt_originaldata.ReplaceAll(".root","_originalData.txt");
505  char txt_originaldata_c[1024];
506  sprintf(txt_originaldata_c,"txt_files/%s.txt",txt_originaldata.Data());
507  ofstream od_txt(txt_originaldata_c);
508 
509 
510  //for(int n=s;n<s+TS;n++) {
511  for(int n=s;n<sig_entries;n++) {
512 double ti_0=0;
513 double fi_0=0;
514 
515  index_ev=n;
516  average=0.;
517  int indNN=-1;
518  int countNN=0;
519  for(int u=0;u<nNN;u++){
520  if (uf!=0&&n>=NNs[u]&&n<(NNs[u]+NNnTS[u])) { indNN=u;countNN=countNN+1;}
521  }
522  if (countNN>1){cout<<"Error: training non independent"; exit(0);}
523  if(nNN==1){cout<<"out==Averege:only 1NN is considered"<<endl;}
524  if(consider_all==0&&av_on_nNN==0&&countNN==0)continue; //skip if no network is trained on the considered event
525  if(consider_all==0&&av_on_nNN!=0&&countNN!=0) {sigskip=sigskip+1; continue; }
526  NNTree->GetEntry(n);
527  signal.GetEntry(n);
528  char line_od[1024];
529 
530  // exracting parameters
531  NNcc=(double)signal.netcc[0];
532  NNrho=(double)signal.rho[1];
533  Mc=(double)signal.chirp[iMc];
534  Mc_inj=(double)signal.chirp[0];
535 
536  sprintf(TestFile,"%s",TEST_FILE.Data());
537 
538  sprintf(line_od,"%i \n",n);
539  od_txt<<line_od;
540  char line_od2[1024];
541  // evaluating the event
542  for (int i=0; i<nINP; i++){
543  params[i]=x[i];
544 
545  int ti_i=0;
546  int fi_i=0;
547  ti_i=i/NDIM;
548  fi_i=i-ti_i*NDIM;
549  ti_0+=(double)x[i]*ti_i;
550  fi_0+=(double)x[i]*fi_i;
551  sprintf(line_od2,"%s %1.5f",line_od2, x[i]);
552  }
553  sprintf(line_od2,"%s \n",line_od2);
554 
555  if (ti_0<fi_0) CoG=1;
556  else CoG=0;
557 
558 //cout<<" ti_0 "<<ti_0<<" fi_0 "<<fi_0<<" CoG "<<CoG<<endl;
559  od_txt<<line_od2;
560  for (int u=0;u<nNN;u++) {
561 // cout<<" u "<<u<<" nNN "<<nNN<<endl;
562  double output=mlp[u]->Evaluate(0,params);
563 // cout<<output<<endl;
564  out[u]=output;
565  }
566 
567 // cout<<"Dopo param"<<endl;
568 
569 //calculating the average and the standard deviation
570 
571  int c_nNN0=0;
572  for (int u=0;u<nNN;u++){
573  if((u!=indNN&&ni==0)||uf==0) {average=out[u]+average;std+=out[u]*out[u];c_nNN0=c_nNN0+1; if(out[u]<0.6) FD0[u]=FD0[u]+1;}
574  }
575 
576  average=average/(c_nNN0);
577  if((c_nNN0)!=1)std=pow((std/(nNN-1-countNN)-average*average),0.5);
578  else std=-1;
579 
580  char line_od3[1024];
581  sprintf(line_od3,"ev_ind %i average %1.5f Mc %1.5f Mc_inj %1.5f n %i \n",index_ev, average,Mc,Mc_inj,n);
582  od_txt<<line_od3;
583 
584  if (average<0.6) FD=FD+1;
585  NNTree2->Fill();
586  S_eff=S_eff+1;
587  if(S_eff==1) s=n;
588  if(S_eff==TS) break;
589  cout<<" cc "<<NNcc<<" NNrho "<<NNrho<<" Mc_inj "<< Mc_inj<<endl;
590 }
591 
592 od_txt.close();
593 //exit(0);
594 cout<<txt_originaldata_c<<endl;
595 //exit(0);
596 //closing the cycle on SIG events
597 cout<<" S_eff "<<S_eff<<" s "<<s<<" b "<<b<<" sigskip "<<sigskip<<endl;
598 //exit(0);
599 //defining BKG parameters
600 
601 int bgskip=0;
602 //opening the cycle on BKG events
603  for(int n=b;n<sig_entries+bg_entries;n++) {
604 double ti_0=0;
605 double fi_0=0;
606  index_ev=0;
607  index_ev=n;
608  average=0.;
609  int indNN=-1;
610  int countNN=0;
611  for(int u=0;u<nNN;u++){
612  cout<<" uf "<<uf<<" n "<<n<<" u "<<u<<" NN b[u] "<<NNb[u]<<" NNb[u]+NNnTB[u]"<<NNb[u]+NNnTB[u]<<endl;
613  if (uf!=0&&n>=NNb[u]&&n<(NNb[u]+NNnTB[u])) { indNN=u;countNN=countNN+1; }
614  }
615  if (countNN>1){cout<<"Error: training non independent"; exit(0);}
616  if(nNN==1){cout<<"out==Averege:only 1NN is considered"<<endl;}
617  cout<<"consider_all "<<consider_all<<" av_on_nNN "<<av_on_nNN<<" countNN "<<countNN<<" indNN "<<indNN<<endl;
618  if(consider_all==0&&av_on_nNN==0&&countNN==0)continue; //skip if no network is trained on the considered event
619  if(consider_all==0&&av_on_nNN!=0&&countNN!=0) { cout<<" n skipped "<<n<<endl; bgskip=bgskip+1; continue;}
620  NNTree->GetEntry(n);
621  cout<<"BKG->n: "<<n<<"Bg index"<<(n-sig_entries)<<endl;
622  background.GetEntry(n-sig_entries);
623  //NNcc=(double)background.netcc[1];
624  NNcc=(double)background.netcc[0];
625  //NNrho=(double)background.rho[0];
626  NNrho=(double)background.rho[1];
627  Mc=(double)background.chirp[iMc];
628  Mc_inj=(double)background.chirp[0];
629  /*int hn=0;
630  ev_ind=0;
631  hn=n;
632  ev_ind=hn;*/
633  sprintf(TestFile,"%s",TEST_FILE.Data());
634 
635 for (int i=0; i<nINP; i++){
636  params[i]=x[i];
637  int ti_i=0;
638  int fi_i=0;
639  ti_i=i/NDIM;
640  fi_i=i-ti_i*NDIM;
641  ti_0+=(double)x[i]*ti_i;
642  fi_0+=(double)x[i]*fi_i;
643 
644  }
645  if (ti_0<fi_0) CoG=1;
646  else CoG=0;
647 
648  for(int u=0;u<nNN;u++) {
649  double output=mlp[u]->Evaluate(0,params);
650  cout<<output<<endl;
651  out[u]=output;
652  }
653 
654  int c_nNN0=0;
655  int nnsu=0;
656  for (int u=0;u<nNN;u++){
657  if((u!=indNN&&ni==0)||uf==0) {average=out[u]+average;std+=out[u]*out[u];c_nNN0=c_nNN0+1; TotBg[u]=TotBg[u]+1;if(out[u]>0.6) {FA0[u]=FA0[u]+1;nnsu=nnsu+1;}}
658  }
659 
660 
661  average=average/(c_nNN0);
662  if((c_nNN0)!=1)std=pow((std/(nNN-1-countNN)-average*average),0.5);
663  else std=-1;
664 
665  if(nnsu==6) cont_su=cont_su+1;
666  if(nnsu>4) cont_su5=cont_su5+1;
667  if(nnsu>3) cont_su4=cont_su4+1;
668  if(nnsu>2) cont_su3=cont_su3+1;
669 
670 //calculating the final average value
671 
672 
673 //saving false alarms
674  if(average>0.6) FA=FA+1;
675 
676  NNTree2->Fill();
677  B_eff=B_eff+1;
678  if(B_eff==1) b=n;
679  if(B_eff==TB) break;
680  }
681 //cout<<"bkg filled"<<endl;
682 
683 NNTree2->Write();
684 f->Close();
685 cout<<"B_eff "<<B_eff<<" TB "<<TB<<endl;
686 //cout<<"closed file"<<endl;
687 cout<<ofile<<endl;
688 
689 int cont=0;
690 
691 double freq_c[nNN];
692 double meanf=0.;
693 double dev=0.;
694 
695 for (int yy=0; yy<nNN; yy++){
696  freq_c[yy]=0.;
697 // cout<<" FA0[yy] "<<FA0[yy]<<" FD0[yy] "<<FD0[yy]<<" FD "<<FD<<" TotBg[yy] "<<TotBg[yy]<<endl;
698  if(TotBg[yy]!=0){freq_c[yy]=(double)FA0[yy]/TotBg[yy];
699  //cout<<" freq_c[yy] "<<freq_c[yy]<<endl;
700  }
701  if(freq_c[yy]!=0) {meanf=freq_c[yy]+meanf;cont=cont+1;
702  //cout<<" cont "<<cont<<endl;
703  }
704  dev=freq_c[yy]*freq_c[yy]+dev;
705 }
706 
707 if(cont==0) cout<<"Error cont==0"<<endl;
708 if(cont!=0) meanf=meanf/cont;
709 dev=0.;
710 for (int yy=0; yy<nNN; yy++){
711 //cout<<" dev "<<dev<<endl;
712 //cout<<"freq_c[yy] "<<freq_c[yy]<<" meanf "<<meanf<<endl;
713 if(freq_c[yy]!=0)dev+=(freq_c[yy]-meanf)*(freq_c[yy]-meanf);
714 }
715 double McFAdes=FAdes;
716 int FAth_n=FAdes*B_eff;
717 int McFAth_n=FAdes*B_eff;
718 cout<<"FAth_n "<<FAth_n<<" FAdes*B_eff "<<FAdes*B_eff<<endl;
719 int FDth_n=-1;
720 int McFDth_n=-1;
721 double thresFA=-1;
722 double McthresFA=-1;
723 //Annth(ofile, FAth_n, FDth_n, thresFA);
724 //Mcth(ofile, McFAth_n, McFDth_n, McthresFA);
725 //graph(ofile);
726 //PlotsAv_cc(ofile);
727 //PlotsAv_Mc(ofile);
728 //SNR_plots(ofile);
729 //CoG_plots(ofile);
730 CutMcANN(ofile);
731 if(cont!=0) dev=pow(dev/(cont-1),0.5);
732 cout<<"meanf "<<meanf<<" dev "<<dev<<endl;
733 cout<<" FA average "<<FA<<endl;
734 cout<<" cont_su "<<cont_su<<" cont_su5 "<<cont_su5<<" cont_su4 "<<cont_su4<<" cont_su3 "<<cont_su3<<endl;
735 
736 cout<<" S_eff "<<S_eff<<" B_eff "<<B_eff<<" s "<<s<<" b "<<b<<" bgskip "<<bgskip<<" sigskip "<<sigskip<<endl;
737 cout<<" FA0[1] "<<FA0[1]<<" TotBg[1] "<<TotBg[1]<<endl;
738 cout<<" FD0[1] "<<FD0[1]<<" TotBg[1] "<<TotBg[1]<<endl;
739 //double S_eff2=(double
740 double FDdes=(double)FDth_n/S_eff;
741 double McFDdes=(double)McFDth_n/S_eff;
742 FAdes=(double)FAth_n/B_eff;
743 McFAdes=(double)McFAth_n/B_eff;
744 cout<<"FAth_n "<<FAth_n<<" FAdes: "<<FAdes<<" ANN_av_FDdes "<<FDdes<<" ANN_av_FDth_n "<<FDth_n<<" ANN_av_thresFA "<<thresFA<<endl;
745 cout<<"McFAth_n "<<McFAth_n<<" McFAdes: "<<McFAdes<<" Mc_FDdes "<<McFDdes<<" McFDth_n "<<McFDth_n<<" McthresFA "<<McthresFA<<endl;
746 }
747 
748 void graph(TString ifile){
749  TString name(ifile);
750  name.ReplaceAll("average_file/","");
751  TFile* fTEST =TFile::Open(ifile.Data());
752  TTree* NNTree2=(TTree*)fTEST->Get("Parameters");
753  double av;
754  double cc;
755  double rho;
756  int nSi;
757  NNTree2->SetBranchAddress("Average",&av);
758  NNTree2->SetBranchAddress("cc",&cc);
759  NNTree2->SetBranchAddress("rho",&rho);
760  NNTree2->SetBranchAddress("#TestSig",&nSi);
761  NNTree2->GetEntry(0);
762  int const nSig=nSi;
763  cout<<"nSig: "<<nSig<<" nSi: "<<nSi<<endl;
764  cout<<" BG: "<<NNTree2->GetEntries()-nSi<<endl;
765  int const ncurve=nANN*nCC;
766  cout<<"dentro funzione dopodef"<<endl;
767 //LOG(NBg)vs RHO------------------------------------------
768  double* rhoSig[ncurve];
769  for (int i=0;i<ncurve;i++) rhoSig[i]=new double[nSig];
770  cout<<"dopo def rhoSig"<<endl;
771  //double rhoSig[20][10];
772  int NSig[ncurve];
773 // int nBg=0;
774  int const nBg=NNTree2->GetEntries()-nSig;
775  for (int i=0;i<ncurve;i++) {
776  NSig[i]=0;
777  for (int j=0;j<nSig;j++) rhoSig[i][j]=0.;
778  }
779  cout<<"dopo def rhoSig"<<endl;
780  double* rhoBg[ncurve];
781  for (int i=0;i<ncurve;i++) rhoBg[i]=new double[nBg];
782  //double rhoBg[20][10];
783  int NBg[ncurve];
784  for (int i=0;i<ncurve;i++) {
785  NBg[i]=0;
786  for (int j=0;j<nBg;j++) rhoBg[i][j]=0.;
787  }
788  cout<<"dopo def rhoBg"<<endl;
789  double ccTh[nCC];
790  for (int i=0;i<nCC;i++) ccTh[i]=0.;
791  double NNTh[nANN];
792  for (int i=0;i<nANN;i++) NNTh[i]=0.;
793  //double deltacc=0.;
794 // double deltaANN=0.;
795 // deltacc=0.2/nCC;
796  //deltaANN=0.6/(nANN-1);
797 
798  cout<<NNTree2->GetEntries()<<endl;
799  for(int n=0;n<NNTree2->GetEntries();n++){
800  cout<<n<<endl;
801  NNTree2->GetEntry(n);
802  cout<<"rho "<<rho<<" cc "<<cc<<" av "<<av<<endl;
803  for(int i=0; i<nCC;i++){
804  ccTh[i]=0.5+i*deltacc;
805  if(cc<ccTh[i]) continue;
806  for(int j=0; j<nANN;j++){
807  if(j==0) NNTh[j]=-1000.;
808  //else NNTh[j]=0.5+(j-1)*deltaANN;
809  //else NNTh[j]=0.1+(j-1)*deltaANN;
810  else NNTh[j]=0.+(j)*deltaANN;
811  //else NNTh[j]=0.+(j-1)*deltaANN/1000.;
812  //else NNTh[j]=1.;
813  if(av<NNTh[j]) continue;
814  int ni=0;
815  if(n>=nSig) {
816 
817  NBg[i*nANN+j]= NBg[i*nANN+j]+1;
818  while(rhoBg[i*nANN+j][ni]!=0)ni=ni+1;
819  rhoBg[i*nANN+j][ni]=rho;
820  // cout<<"rho: "<<rho<<" colonna "<<i*nANN+j<<" riga "<<ni<<endl;
821 // cout<<" soglia_cc "<<ccTh[i]<<" soglia_ANN "<<NNTh[j]<<endl;
822  }
823  else {
824  NSig[i*nANN+j]= NSig[i*nANN+j]+1;
825  while(rhoSig[i*nANN+j][ni]!=0)ni=ni+1;
826  rhoSig[i*nANN+j][ni]=rho;
827  // cout<<"rho: "<<rho<<" colonna "<<i*nANN+j<<" riga "<<ni<<endl;
828 // cout<<" soglia_cc "<<ccTh[i]<<" soglia_ANN "<<NNTh[j]<<endl;
829  }
830  // }
831  }
832  }
833  }
834  cout<<"dopo riempimento variabili"<<endl;
835  int* indexS[ncurve];
836  for (int i=0;i<ncurve;i++) indexS[i]=new int[nSig];
837 
838  TGraph * gS[ncurve];
839  for (int y=0;y<ncurve;y++) {
840  int igS=0;
841  int igS_p=0;
842  /* int indexS[nSig];
843  double rhoS[nSig];
844  for(int i=0;i<nSig;i++) {
845  rhoS[i]=0.;
846  indexS[i]=0;
847  }*/
848  // ig=1;
849  // for(int i=0;i<nSig;i++) rhoS[i]=rhoSig[y][i];
850  gS[y]=new TGraph();
851 // gS[y]->SetMarkerStyle(7);
852  TMath::Sort(nSig,rhoSig[y],indexS[y],false);
853 // cout<<"dopo Sort "<<y<<endl;
854  for (int k=0;k<nSig;k++) {
855  int ii=indexS[y][k];
856  int yy=0;
857  if (k>0){
858 // cout<<"k "<<k<<endl;
859  int ij=indexS[y][k-1];
860  //if(rhoSig[y][ii]!=0&&rhoSig[y][ii]!=rhoSig[y][ij]) {
861  if(rhoSig[y][ii]!=0) {
862  yy=NSig[y]-igS;
863  //gS[y]->SetPoint(igS,rhoSig[y][ii],yy);
864  if(rhoSig[y][ii]!=rhoSig[y][ij]) gS[y]->SetPoint(igS_p++,rhoSig[y][ii],yy);
865  cout<<"igS"<<igS<<" x "<<rhoSig[y][ii]<<" y: "<<yy<<endl;
866  igS=igS+1;
867  }
868  }
869  else {
870  if(rhoSig[y][ii]!=0){
871  yy=NSig[y]-igS;
872  gS[y]->SetPoint(0,rhoSig[y][ii],yy);
873  igS=igS+1;
874 // cout<<" x "<<rhoSig[y][ii]<<" y: "<<yy<<endl;
875  }
876 
877  }
878  }
879  }
880  // cout<<"dopo inserimento puntiiS"<<endl;
881 // cout<<ncurve<<endl;
882  int* indexB[ncurve];
883  for (int i=0;i<ncurve;i++) indexB[i]=new int[nBg];
884 
885  TGraph * gB[ncurve];
886  for (int y=0;y<ncurve;y++) {
887  int igB=0;
888  int igB_p=0;
889  gB[y]=new TGraph();
890  TMath::Sort(nBg,rhoBg[y],indexB[y],false);
891 // cout<<"dopo Sort "<<y<<endl;
892  for (int k=0;k<nBg;k++) {
893  int ii=indexB[y][k];
894  int yy=0;
895  if (k>0){
896  int ij=indexB[y][k-1];
897  //if(rhoBg[y][ii]!=0&&rhoBg[y][ii]!=rhoBg[y][ij]) {
898  if(rhoBg[y][ii]!=0) {
899  yy=NBg[y]-igB;
900  //gB[y]->SetPoint(igB,rhoBg[y][ii],yy);
901  if(rhoBg[y][ii]!=rhoBg[y][ij]) gB[y]->SetPoint(igB_p++,rhoBg[y][ii],yy);
902  igB=igB+1;
903 // cout<<"igB"<<igB<<" x "<<rhoBg[y][ii]<<" y: "<<yy<<endl;
904  }
905  }
906  else {
907  if(rhoBg[y][ii]!=0) {
908  yy=NBg[y]-igB;
909  gB[y]->SetPoint(0,rhoBg[y][ii],yy);
910  igB=igB+1;
911 // cout<<"igB"<<igB<<" x "<<rhoBg[y][ii]<<" y: "<<yy<<endl;
912  }
913  }
914  }
915  }
916  cout<<"dopo inserimento puntiB"<<endl;
917  //gB[1]->SetMarkerStyle(7);
918  //gB[1]->SetPoint(1,1,2);
919  TCanvas* cS=new TCanvas("Efficiency_vs_rho","Efficiency_vs_rho",0,0,1200,700);
920  cS->Divide(2,2);
921  //cS->cd(1)->SetLogy();
922  cS->cd(1);
923  TMultiGraph* mg1=new TMultiGraph();
924 // gS[0]->SetMarkerStyle(7);
925  gS[0]->SetMarkerColor(2);
926  gS[0]->SetLineColor(2);
927  mg1->SetTitle("cc=0.5;rho;#Events");
928  if(gS[0]->GetN()!=0) mg1->Add(gS[0]);
929 // gS[0]->Draw("apl");
930  for (int h=1;h<nANN;h++){
931  gS[h]->SetMarkerColor(3);
932  gS[h]->SetLineColor(3);
933 // gS[h]->SetMarkerStyle(h+1);
934  if(gS[h]->GetN()!=0) mg1->Add(gS[h]);
935  // gS[h]->Draw("apl,same");
936  }
937  mg1->Draw("apl");
938  //cS->cd(2)->SetLogy();
939  cS->cd(2);
940  TMultiGraph* mg2=new TMultiGraph();
941 // gS[nANN]->SetMarkerStyle(7);
942  gS[nANN]->SetMarkerColor(2);
943  gS[nANN]->SetLineColor(2);
944  mg2->SetTitle("cc=0.55;rho;#Events");
945  if(gS[nANN]->GetN()!=0) mg2->Add(gS[nANN]);
946  for (int h=1;h<nANN;h++){
947  gS[nANN+h]->SetMarkerColor(3);
948  gS[nANN+h]->SetLineColor(3);
949 // gS[nANN+h]->SetMarkerStyle(h+1);
950  if(gS[nANN+h]->GetN()!=0) mg2->Add(gS[nANN+h]);
951  // gS[nANN+h]->Draw("apl,same");
952  }
953  mg2->Draw("apl");
954  //cS->cd(3)->SetLogy();
955  cS->cd(3);
956  TMultiGraph* mg3=new TMultiGraph();
957 // gS[nANN*2]->SetMarkerStyle(7);
958  gS[nANN*2]->SetMarkerColor(2);
959  gS[nANN*2]->SetLineColor(2);
960  mg3->SetTitle("cc=0.6;rho;#Events");
961  if(gS[nANN*2]->GetN()!=0) mg3->Add(gS[nANN*2]);
962  for (int h=1;h<nANN;h++){
963  gS[2*nANN+h]->SetMarkerColor(3);
964  gS[2*nANN+h]->SetLineColor(3);
965 // gS[2*nANN+h]->SetMarkerStyle(h+1);
966  if(gS[2*nANN+h]->GetN()!=0) mg3->Add(gS[2*nANN+h]);
967  // gS[2*nANN+h]->Draw("apl,same");
968  }
969  mg3->Draw("apl");
970  //cS->cd(4)->SetLogy();
971  cS->cd(4);
972  TMultiGraph* mg4=new TMultiGraph();
973 // gS[nANN*3]->SetMarkerStyle(7);
974  gS[nANN*3]->SetMarkerColor(2);
975  gS[nANN*3]->SetLineColor(2);
976  mg4->SetTitle("cc=0.65;rho;#Events");
977  if(gS[nANN*3]->GetN()!=0) mg4->Add(gS[nANN*3]);
978  // gS[nANN*3]->Draw("apl");
979  for (int h=1;h<nANN;h++){
980  gS[3*nANN+h]->SetMarkerColor(3);
981  gS[3*nANN+h]->SetLineColor(3);
982 // gS[3*nANN+h]->SetMarkerStyle(h+1);
983  if(gS[3*nANN+h]->GetN()!=0) mg4->Add(gS[3*nANN+h]);
984  // gS[3*nANN+h]->Draw("apl,same");
985  }
986  mg4->Draw("apl");
987  cout<<"nuovo canv"<<endl;
988  TCanvas* cB=new TCanvas("Number_vs_rho","Number_vs_rho",0,0,1200,700);
989  cB->Divide(2,2);
990  cB->cd(1)->SetLogy();
991  TMultiGraph* mg1B=new TMultiGraph();
992 // gB[0]->SetMarkerStyle(7);
993  gB[0]->SetMarkerColor(2);
994  gB[0]->SetLineColor(2);
995  mg1B->SetTitle("cc=0.5;rho;#Events");
996  if(gB[0]->GetN()!=0) mg1B->Add(gB[0]);
997  for (int h=1;h<nANN;h++){
998  gB[h]->SetMarkerColor(3);
999  gB[h]->SetLineColor(3);
1000  // gB[h]>SetMarkerStyle(h+1);
1001  if(gB[h]->GetN()!=0) mg1B->Add(gB[h]);
1002  //gB[h]->Draw("apl,same");
1003  }
1004  mg1B->Draw("apl");
1005  cB->cd(2)->SetLogy();
1006  TMultiGraph* mg2B=new TMultiGraph();
1007  //gB[nANN]->SetMarkerStyle(7);
1008  gB[nANN]->SetMarkerColor(2);
1009  gB[nANN]->SetLineColor(2);
1010  mg2B->SetTitle("cc=0.55;rho;#Events");
1011  if(gB[nANN]->GetN()!=0) mg2B->Add(gB[nANN]);
1012  for (int h=1;h<nANN;h++){
1013  gB[nANN+h]->SetMarkerColor(3);
1014  gB[nANN+h]->SetLineColor(3);
1015  //gB[nANN+h]>SetMarkerStyle(h+1);
1016  if(gB[nANN+h]->GetN()!=0) mg2B->Add(gB[nANN+h]);
1017  //gB[h]->Draw("apl,same");
1018  }
1019  mg2B->Draw("apl");
1020  cB->cd(3)->SetLogy();
1021  TMultiGraph* mg3B=new TMultiGraph();
1022 // gB[2*nANN]->SetMarkerStyle(7);
1023  gB[2*nANN]->SetMarkerColor(2);
1024  gB[2*nANN]->SetLineColor(2);
1025  mg3B->SetTitle("cc=0.6;rho;#Events");
1026  if(gB[2*nANN]->GetN()!=0) mg3B->Add(gB[2*nANN]);
1027  for (int h=1;h<nANN;h++){
1028  gB[2*nANN+h]->SetMarkerColor(3);
1029  gB[2*nANN+h]->SetLineColor(3);
1030 // gB[2*nANN+h]>SetMarkerStyle(h+1);
1031  if(gB[2*nANN+h]->GetN()!=0) mg3B->Add(gB[2*nANN+h]);
1032  //gB[h]->Draw("apl,same");
1033  }
1034  mg3B->Draw("apl");
1035  cB->cd(4)->SetLogy();
1036  TMultiGraph* mg4B=new TMultiGraph();
1037 // gB[3*nANN]->SetMarkerStyle(7);
1038  gB[3*nANN]->SetMarkerColor(2);
1039  gB[3*nANN]->SetLineColor(2);
1040  mg4B->SetTitle("cc=0.65;rho;#Events");
1041  if(gB[3*nANN]->GetN()!=0) mg4B->Add(gB[3*nANN]);
1042  for (int h=1;h<nANN;h++){
1043  gB[3*nANN+h]->SetMarkerColor(3);
1044  gB[3*nANN+h]->SetLineColor(3);
1045 // gB[3*nANN+h]>SetMarkerStyle(h+1);
1046  if(gB[3*nANN+h]->GetN()!=0) mg4B->Add(gB[3*nANN+h]);
1047  //gB[h]->Draw("apl,same");
1048  }
1049  mg4B->Draw("apl");
1050 
1051 // cout<<"dopo Draw()"<<endl;
1052  TString CnameS(name);
1053  TString CnameB(name);
1054  TString CnameSroot(name);
1055  TString CnameBroot(name);
1056  char CnameS2[1024];
1057  char CnameB2[1024];
1058  char CnameS2root[1024];
1059  char CnameB2root[1024];
1060  CnameS.ReplaceAll(".root",".png");
1061  CnameB.ReplaceAll(".root",".png");
1062  sprintf(CnameS2,"logN_rho_av/logN_rho_S_dANN%1.2f_%s",deltaANN,CnameS.Data());
1063  sprintf(CnameB2,"logN_rho_av/logN_rho_B_dANN%1.2f_%s",deltaANN,CnameB.Data());
1064  sprintf(CnameS2root,"logN_rho_av/logN_rho_S_dANN%1.2f_%s",deltaANN,CnameSroot.Data());
1065  sprintf(CnameB2root,"logN_rho_av/logN_rho_B_dANN%1.2f_%s",deltaANN,CnameBroot.Data());
1066  cS->SaveAs(CnameS2);
1067  cB->SaveAs(CnameB2);
1068  cS->Print(CnameS2root);
1069  cB->Print(CnameB2root);
1070 // cout<<"fine"<<endl;
1071 
1072 }
1073 
1074 
1075 void Annth(TString ifile, int &FAth_n, int &FDth_n, double &thresFA){
1076  TString name(ifile);
1077  name.ReplaceAll("average_file/","");
1078  TFile* fTEST =TFile::Open(ifile.Data());
1079  TTree* NNTree2=(TTree*)fTEST->Get("Parameters");
1080  double av;
1081  double cc;
1082  double rho;
1083  int nSi;
1084  NNTree2->SetBranchAddress("Average",&av);
1085  NNTree2->SetBranchAddress("cc",&cc);
1086  NNTree2->SetBranchAddress("rho",&rho);
1087  NNTree2->SetBranchAddress("#TestSig",&nSi);
1088  NNTree2->GetEntry(0);
1089  int const nSig=nSi;
1090  cout<<"nSig: "<<nSig<<" nSi: "<<nSi<<endl;
1091  int const ncurve2=nRHO*nCC;
1092 
1093 
1094  double* ANNSig[ncurve2];
1095  for (int i=0;i<ncurve2;i++) ANNSig[i]=new double[nSig];
1096  //double rhoSig[20][10];
1097  int NSig[ncurve2];
1098 // int nBg=0;
1099  int const nBg=NNTree2->GetEntries()-nSig;
1100  for (int i=0;i<ncurve2;i++) {
1101  NSig[i]=0;
1102  for (int j=0;j<nSig;j++) ANNSig[i][j]=0.;
1103  }
1104  double* ANNBg[ncurve2];
1105  for (int i=0;i<ncurve2;i++) ANNBg[i]=new double[nBg];
1106 
1107  //double rhoBg[20][10];
1108  int NBg[ncurve2];
1109  for (int i=0;i<ncurve2;i++) {
1110  NBg[i]=0;
1111  for (int j=0;j<nBg;j++) ANNBg[i][j]=0.;
1112  }
1113  double ccTh[nCC];
1114  for (int i=0;i<nCC;i++) ccTh[i]=0.;
1115  double rhoTh[nRHO];
1116  for (int i=0;i<nRHO;i++) rhoTh[i]=0.;
1117  // double deltacc=0.;
1118  // double deltarho=0.;
1119  // deltacc=0.2/nCC;
1120 // deltarho=1./(nRHO);
1121 
1122 //int thresFA=-1;
1123  for(int n=0;n<NNTree2->GetEntries();n++){
1124  NNTree2->GetEntry(n);
1125  cout<<"rho "<<rho<<" cc "<<cc<<" av "<<av<<endl;
1126  for(int i=0; i<nCC;i++){
1127  ccTh[i]=0.5+i*deltacc;
1128  ccTh[0]=0.;
1129  if(cc<ccTh[i]) continue;
1130  for(int j=0; j<nRHO;j++){
1131  rhoTh[j]=5+j*deltarho;
1132  rhoTh[0]=0.;
1133  if(rho<rhoTh[j]) continue;
1134  int ni=0;
1135  if(n>=nSig) {
1136  NBg[i*nRHO+j]= NBg[i*nRHO+j]+1;
1137  if(i*nRHO+j==0) cout<<" NBg[i*nRHO+j] "<<NBg[i*nRHO+j]<<endl;
1138  while(ANNBg[i*nRHO+j][ni]!=0)ni=ni+1;
1139  ANNBg[i*nRHO+j][ni]=av;
1140  // if (NBg[0]==FAth_n) thresFA=av;
1141  // cout<<" soglia_cc "<<ccTh[i]<<" soglia_RHO "<<rhoTh[j]<<endl;
1142  }
1143  else {
1144  NSig[i*nRHO+j]= NSig[i*nRHO+j]+1;
1145  if(i*nRHO+j==0) cout<<" NSig[i*nRHO+j] "<<NSig[i*nRHO+j]<<endl;
1146  while(ANNSig[i*nRHO+j][ni]!=0)ni=ni+1;
1147  ANNSig[i*nRHO+j][ni]=av;
1148  // cout<<" soglia_cc "<<ccTh[i]<<" soglia_RHO "<<rhoTh[j]<<endl;
1149  }
1150  }
1151  }
1152  }
1153  int* indexS[ncurve2];
1154  for (int i=0;i<ncurve2;i++) indexS[i]=new int[nSig];
1155 
1156  TGraph * gS[ncurve2];
1157  for (int y=0;y<ncurve2;y++) {
1158  int igS=0;
1159  int igS_p=0;
1160  gS[y]=new TGraph();
1161  TMath::Sort(nSig,ANNSig[y],indexS[y],false);
1162  for (int k=0;k<nSig;k++) {
1163  int ii=indexS[y][k];
1164  int yy=0;
1165  if (k>0){
1166  //cout<<"k "<<k<<endl;
1167  int ij=indexS[y][k-1];
1168  //if(ANNSig[y][ii]!=0&&ANNSig[y][ii]!=ANNSig[y][ij]) {
1169  if(ANNSig[y][ii]!=0) {
1170  yy=NSig[y]-igS;
1171  //gS[y]->SetPoint(igS,ANNSig[y][ii],yy);
1172  if(ANNSig[y][ii]!=ANNSig[y][ij]) gS[y]->SetPoint(igS_p++,ANNSig[y][ii],yy);
1173  // cout<<"igS"<<igS<<" x "<<ANNSig[y][ii]<<" y: "<<yy<<endl;
1174  igS=igS+1;
1175 
1176  }
1177  }
1178  else {
1179  if(ANNSig[y][ii]!=0){
1180  yy=NSig[y]-igS;
1181  gS[y]->SetPoint(0,ANNSig[y][ii],yy);
1182  igS=igS+1;
1183  // cout<<" x "<<ANNSig[y][ii]<<" y: "<<yy<<endl;
1184  }
1185 
1186  }
1187  }
1188  }
1189 
1190 
1191 int thres_ind=0;
1192  int* indexB[ncurve2];
1193  for (int i=0;i<ncurve2;i++) indexB[i]=new int[nBg];
1194 
1195  TGraph * gB[ncurve2];
1196  for (int y=0;y<ncurve2;y++) {
1197  int igB=0;
1198  int igB_p=0;
1199  gB[y]=new TGraph();
1200  TMath::Sort(nBg,ANNBg[y],indexB[y],false);
1201  // cout<<"dopo Sort "<<y<<endl;
1202  for (int k=0;k<nBg;k++) {
1203  int ii=indexB[y][k];
1204  int ik=indexB[y][k-1];
1205  int yy=0;
1206  if (k>0){
1207  int ij=indexB[y][k-1];
1208  //if(ANNBg[y][ii]!=0&&ANNBg[y][ii]!=ANNBg[y][ij]) {
1209  if(ANNBg[y][ii]!=0) {
1210  yy=NBg[y]-igB;
1211  if(y==0 && yy==FAth_n+1) {thresFA=ANNBg[y][ii];}
1212  if (y==0 && yy<FAth_n+1) if(ANNBg[y][ii]==ANNBg[y][ik]){FAth_n=FAth_n-1;}
1213  if(y==0) cout<<"yy "<<yy<<" y "<<y<<" FAth_n "<<FAth_n<<" thresFA "<<thresFA<<" ANNBg[y][ii] "<<ANNBg[y][ii]<<endl;
1214  //gB[y]->SetPoint(igB,ANNBg[y][ii],yy);
1215  if(ANNBg[y][ii]!=ANNBg[y][ij]) gB[y]->SetPoint(igB_p++,ANNBg[y][ii],yy);
1216  igB=igB+1;
1217  // cout<<"igB"<<igB<<" x "<<ANNBg[y][ii]<<" y: "<<yy<<endl;
1218  }
1219  }
1220  else {
1221  if(ANNBg[y][ii]!=0) {
1222  yy=NBg[y]-igB;
1223  gB[y]->SetPoint(0,ANNBg[y][ii],yy);
1224  igB=igB+1;
1225  // cout<<"igB"<<igB<<" x "<<ANNBg[y][ii]<<" y: "<<yy<<endl;
1226  }
1227  }
1228  }
1229  }
1230 //int FDth_n=0;
1231 //for(int i=0; i<nCC;i++){
1232 // for(int j=0; j<nRHO;j++){
1233 FDth_n=0;
1234  for(int ni=1;ni<nSig;ni++){
1235  cout<<"FAth_n "<<FAth_n<<" thresFA "<<thresFA<<" ANNSig[0][indexS[0][ni]] "<<ANNSig[0][indexS[0][ni]]<<" FDth_n "<<FDth_n<<" NSig[0]-(ni-1 )"<<NSig[0]-(ni-1)<<" NSig[0] "<<NSig[0]<<endl;
1236 // if (ANNSig[0][indexS[0][ni]]==thresFA || (ANNSig[0][indexS[0][ni-1]]<thresFA && ANNSig[0][indexS[0][ni]]>thresFA)) FDth_n=NSig[0]-(ni-1);
1237  if (ANNSig[0][indexS[0][ni]]<=thresFA) FDth_n=FDth_n+1;
1238  }
1239 // }
1240 //}
1241 //exit(0);
1242  TCanvas* cS=new TCanvas("Efficiency_vs_ANN","Efficiency_vs_ANN",0,0,1200,700);
1243  cS->Divide(2,2);
1244  //cS->cd(1)->SetLogy();
1245  cS->cd(1);
1246  TMultiGraph* mg1=new TMultiGraph();
1247  mg1->SetTitle("cc: no_cut;ANN;#Events");
1248  for (int h=0;h<nRHO;h++){
1249  gS[h]->SetLineColor(4);
1250  if(gS[h]->GetN()!=0) mg1->Add(gS[h]);
1251  }
1252  mg1->Draw("al");
1253  cS->cd(2);
1254  // cS->cd(2)->SetLogy();
1255  TMultiGraph* mg2=new TMultiGraph();
1256  mg2->SetTitle("cc=0.55;ANN;#Events");
1257  for (int h=0;h<nRHO;h++){
1258  gS[nRHO+h]->SetLineColor(4);
1259  if(gS[nRHO+h]->GetN()!=0) mg2->Add(gS[nRHO+h]);
1260  }
1261  mg2->Draw("al");
1262 
1263  //cS->cd(3)->SetLogy();
1264  cS->cd(3);
1265  TMultiGraph* mg3=new TMultiGraph();
1266  mg3->SetTitle("cc=0.6;ANN;#Events");
1267  for (int h=0;h<nRHO;h++){
1268  gS[2*nRHO+h]->SetLineColor(4);
1269  if(gS[2*nRHO+h]->GetN()!=0) mg3->Add(gS[2*nRHO+h]);
1270  }
1271  mg3->Draw("al");
1272  //cS->cd(4)->SetLogy();
1273  cS->cd(4);
1274  TMultiGraph* mg4=new TMultiGraph();
1275  mg4->SetTitle("cc=0.65;ANN;#Events");
1276  for (int h=0;h<nRHO;h++){
1277  gS[3*nRHO+h]->SetLineColor(4);
1278  if(gS[3*nRHO+h]->GetN()!=0) mg4->Add(gS[3*nRHO+h]);
1279  }
1280  mg4->Draw("al");
1281 
1282  TCanvas* cB=new TCanvas("Number_vs_ANN","Number_vs_ANN",0,0,1200,700);
1283  cB->Divide(2,2);
1284  cB->cd(1)->SetLogy();
1285  TMultiGraph* mg1B=new TMultiGraph();
1286  mg1B->SetTitle("cc: no_cut;ANN;#Events");
1287  for (int h=0;h<nRHO;h++){
1288  gB[h]->SetLineColor(4);
1289  if(gB[h]->GetN()!=0) mg1B->Add(gB[h]);
1290  }
1291  mg1B->Draw("al");
1292  cB->cd(2)->SetLogy();
1293  TMultiGraph* mg2B=new TMultiGraph();
1294  mg2B->SetTitle("cc=0.55;ANN;#Events");
1295  for (int h=0;h<nRHO;h++){
1296  gB[nRHO+h]->SetLineColor(4);
1297  if(gB[nRHO+h]->GetN()!=0) mg2B->Add(gB[nRHO+h]);
1298  }
1299  mg2B->Draw("al");
1300  cB->cd(3)->SetLogy();
1301  TMultiGraph* mg3B=new TMultiGraph();
1302  mg3B->SetTitle("cc=0.6;ANN;#Events");
1303  for (int h=0;h<nRHO;h++){
1304  gB[2*nRHO+h]->SetLineColor(4);
1305  if(gB[2*nRHO+h]->GetN()!=0) mg3B->Add(gB[2*nRHO+h]);
1306  }
1307  mg3B->Draw("al");
1308  cB->cd(4)->SetLogy();
1309  TMultiGraph* mg4B=new TMultiGraph();
1310  mg4B->SetTitle("cc=0.65;ANN;#Events");
1311  for (int h=0;h<nRHO;h++){
1312  gB[3*nRHO+h]->SetLineColor(4);
1313  if(gB[3*nRHO+h]->GetN()!=0) mg4B->Add(gB[3*nRHO+h]);
1314  }
1315  mg4B->Draw("al");
1316 
1317 
1318  //cout<<"dopo Draw()"<<endl;
1319  TString CnameS(name);
1320  TString CnameB(name);
1321  TString CnameSroot(name);
1322  TString CnameBroot(name);
1323  char CnameS2[1024];
1324  char CnameB2[1024];
1325  char CnameS2root[1024];
1326  char CnameB2root[1024];
1327  CnameS.ReplaceAll(".root",".png");
1328  CnameB.ReplaceAll(".root",".png");
1329  sprintf(CnameS2,"ANNthres_av/N_ANN_S_%s",CnameS.Data());
1330  sprintf(CnameB2,"ANNthres_av/N_ANN_B_%s",CnameB.Data());
1331  sprintf(CnameS2root,"ANNthres_av/N_ANN_S_%s",CnameSroot.Data());
1332  sprintf(CnameB2root,"ANNthres_av/N_ANN_B_%s",CnameBroot.Data());
1333  cS->SaveAs(CnameS2);
1334  cB->SaveAs(CnameB2);
1335  cS->Print(CnameS2root);
1336  cB->Print(CnameB2root);
1337  //cout<<"fine"<<endl;
1338 
1339 //CARTELLA Annth
1340 
1341 
1342 
1343 }
1344 
1345 void Mcth(TString ifile, int &McFAth_n, int &McFDth_n, double &McthresFA){
1346 TString name(ifile);
1347  name.ReplaceAll("outfile/","");
1348  name.ReplaceAll("average_file/","");
1349  TFile* fTEST =TFile::Open(ifile.Data());
1350  TTree* NNTree2=(TTree*)fTEST->Get("Parameters");
1351  //double outbis;
1352  double av;
1353  double out;
1354  double cc;
1355  double Mc;
1356  double Mc_inj;
1357  double rho;
1358  int nSi;
1359 // NNTree2->SetBranchAddress("ANNout",&out);
1360  NNTree2->SetBranchAddress("Mchirp",&Mc);
1361  NNTree2->SetBranchAddress("Average",&av);
1362  NNTree2->SetBranchAddress("cc",&cc);
1363  NNTree2->SetBranchAddress("rho",&rho);
1364  NNTree2->SetBranchAddress("Mchirp_injected",&Mc_inj);
1365  NNTree2->SetBranchAddress("#TestSig",&nSi);
1366  NNTree2->GetEntry(0);
1367  int const nSig=nSi;
1368  cout<<"nSig: "<<nSig<<" nSi: "<<nSi<<endl;
1369  int const ncurve2=nRHO*nCC;
1370 
1371  double* McSig[ncurve2];
1372  for (int i=0;i<ncurve2;i++) McSig[i]=new double[nSig];
1373  //double rhoSig[20][10];
1374  int NSig[ncurve2];
1375 // int nBg=0;
1376  int const nBg=NNTree2->GetEntries()-nSig;
1377  for (int i=0;i<ncurve2;i++) {
1378  NSig[i]=0;
1379  for (int j=0;j<nSig;j++) McSig[i][j]=0.;
1380  }
1381  double* McBg[ncurve2];
1382  for (int i=0;i<ncurve2;i++) McBg[i]=new double[nBg];
1383 
1384  //double rhoBg[20][10];
1385  int NBg[ncurve2];
1386  for (int i=0;i<ncurve2;i++) {
1387  NBg[i]=0;
1388  for (int j=0;j<nBg;j++) McBg[i][j]=0.;
1389  }
1390  double ccTh[nCC];
1391  for (int i=0;i<nCC;i++) ccTh[i]=0.;
1392  double rhoTh[nRHO];
1393  for (int i=0;i<nRHO;i++) rhoTh[i]=0.;
1394 // double deltacc=0.;
1395 // double deltarho=0.;
1396 // deltacc=0.2/nCC;
1397 // deltarho=1./(nRHO);
1398 
1399  for(int n=0;n<NNTree2->GetEntries();n++){
1400  NNTree2->GetEntry(n);
1401  cout<<"rho "<<rho<<" cc "<<cc<<" av "<<av<<endl;
1402  for(int i=0; i<nCC;i++){
1403  ccTh[i]=0.5+i*deltacc;
1404  ccTh[0]=0.;
1405  if(cc<ccTh[i]) continue;
1406  for(int j=0; j<nRHO;j++){
1407  rhoTh[j]=5+j*deltarho;
1408  rhoTh[0]=0.;
1409  if(rho<rhoTh[j]) continue;
1410  int ni=0;
1411  if(n>=nSig) {
1412  NBg[i*nRHO+j]= NBg[i*nRHO+j]+1;
1413  if(i*nRHO+j==0) cout<<" NBg[i*nRHO+j] "<<NBg[i*nRHO+j]<<endl;
1414  while(McBg[i*nRHO+j][ni]!=0)ni=ni+1;
1415  McBg[i*nRHO+j][ni]=Mc;
1416  // if (NBg[0]==FAth_n) thresFA=av;
1417  // cout<<" soglia_cc "<<ccTh[i]<<" soglia_RHO "<<rhoTh[j]<<endl;
1418  }
1419  else {
1420  NSig[i*nRHO+j]= NSig[i*nRHO+j]+1;
1421  if(i*nRHO+j==0) cout<<" NSig[i*nRHO+j] "<<NSig[i*nRHO+j]<<endl;
1422  while(McSig[i*nRHO+j][ni]!=0)ni=ni+1;
1423  McSig[i*nRHO+j][ni]=Mc;
1424  // cout<<" soglia_cc "<<ccTh[i]<<" soglia_RHO "<<rhoTh[j]<<endl;
1425  }
1426  }
1427  }
1428  }
1429 
1430  int* indexS[ncurve2];
1431  for (int i=0;i<ncurve2;i++) indexS[i]=new int[nSig];
1432 
1433  TGraph * gS[ncurve2];
1434  for (int y=0;y<ncurve2;y++) {
1435  int igS=0;
1436  int igS_p=0;
1437  gS[y]=new TGraph();
1438  TMath::Sort(nSig,McSig[y],indexS[y],false);
1439  for (int k=0;k<nSig;k++) {
1440  int ii=indexS[y][k];
1441  int yy=0;
1442  if (k>0){
1443  //cout<<"k "<<k<<endl;
1444  int ij=indexS[y][k-1];
1445  //if(ANNSig[y][ii]!=0&&ANNSig[y][ii]!=ANNSig[y][ij]) {
1446  if(McSig[y][ii]!=0) {
1447  yy=NSig[y]-igS;
1448  //gS[y]->SetPoint(igS,ANNSig[y][ii],yy);
1449  if(McSig[y][ii]!=McSig[y][ij]) gS[y]->SetPoint(igS_p++,McSig[y][ii],yy);
1450  // cout<<"igS"<<igS<<" x "<<ANNSig[y][ii]<<" y: "<<yy<<endl;
1451  igS=igS+1;
1452 
1453  }
1454  }
1455  else {
1456  if(McSig[y][ii]!=0){
1457  yy=NSig[y]-igS;
1458  gS[y]->SetPoint(0,McSig[y][ii],yy);
1459  igS=igS+1;
1460  // cout<<" x "<<ANNSig[y][ii]<<" y: "<<yy<<endl;
1461  }
1462 
1463  }
1464  }
1465  }
1466 
1467 int thres_ind=0;
1468  int* indexB[ncurve2];
1469  for (int i=0;i<ncurve2;i++) indexB[i]=new int[nBg];
1470 
1471  TGraph * gB[ncurve2];
1472  for (int y=0;y<ncurve2;y++) {
1473  int igB=0;
1474  int igB_p=0;
1475  gB[y]=new TGraph();
1476  TMath::Sort(nBg,McBg[y],indexB[y],false);
1477  // cout<<"dopo Sort "<<y<<endl;
1478  for (int k=0;k<nBg;k++) {
1479  int ii=indexB[y][k];
1480  int ik=indexB[y][k-1];
1481  int yy=0;
1482  if (k>0){
1483  int ij=indexB[y][k-1];
1484  //if(ANNBg[y][ii]!=0&&ANNBg[y][ii]!=ANNBg[y][ij]) {
1485  if(McBg[y][ii]!=0) {
1486  yy=NBg[y]-igB;
1487  if(y==0 && yy==McFAth_n+1) {McthresFA=McBg[y][ii];
1488  }
1489  if (y==0 && yy<McFAth_n+1) { cout<<"dentro primo if"<<endl;
1490  //if(McBg[y][ii]==McBg[y][ik]){McFAth_n=McFAth_n-1; cout<<"dentro secondo if"<<endl;}
1491  if(McBg[y][ii]>McBg[y][ik]){McFAth_n=McFAth_n;cout<<"maggiore ii "<<" diff"<<(McBg[y][ii]-McBg[y][ik])<<endl;}
1492  if(McBg[y][ii]<McBg[y][ik]){McFAth_n=McFAth_n;cout<<"minore ii"<<endl;}
1493  if(McBg[y][ii]==McBg[y][ik]){McFAth_n=McFAth_n-1;cout<<"uguale"<<endl;}
1494  //else { McFAth_n=McFAth_n-1;cout<<"dentro secondo if"<<endl;}
1495 }
1496  if(y==0) cout<<"ii "<<ii<<" yy "<<yy<<" y "<<y<<" McFAth_n "<<McFAth_n<<"McthresFA "<<McthresFA<<" McBg[y][ii] "<<McBg[y][ii]<<" McBg[y][indexB[y][k-1] "<<McBg[y][indexB[y][k-1]]<<endl;
1497  //gB[y]->SetPoint(igB,ANNBg[y][ii],yy);
1498  if(McBg[y][ii]!=McBg[y][ij]) gB[y]->SetPoint(igB_p++,McBg[y][ii],yy);
1499  igB=igB+1;
1500  // cout<<"igB"<<igB<<" x "<<ANNBg[y][ii]<<" y: "<<yy<<endl;
1501  }
1502  }
1503  else {
1504  if(McBg[y][ii]!=0) {
1505  yy=NBg[y]-igB;
1506  gB[y]->SetPoint(0,McBg[y][ii],yy);
1507  igB=igB+1;
1508  // cout<<"igB"<<igB<<" x "<<ANNBg[y][ii]<<" y: "<<yy<<endl;
1509  }
1510  }
1511  }
1512  }
1513 
1514 McFDth_n=0;
1515  for(int ni=1;ni<nSig;ni++){
1516  cout<<"McFAth_n "<<McFAth_n<<" McthresFA "<<McthresFA<<" McSig[0][indexS[0][ni]] "<<McSig[0][indexS[0][ni]]<<" McFDth_n "<<McFDth_n<<" NSig[0]-(ni-1 )"<<NSig[0]-(ni-1)<<" NSig[0] "<<NSig[0]<<endl;
1517 // if (ANNSig[0][indexS[0][ni]]==thresFA || (ANNSig[0][indexS[0][ni-1]]<thresFA && ANNSig[0][indexS[0][ni]]>thresFA)) FDth_n=NSig[0]-(ni-1);
1518  if (McSig[0][indexS[0][ni]]<=McthresFA) McFDth_n=McFDth_n+1;
1519  }
1520 
1521 
1522  TCanvas* cS=new TCanvas("Efficiency_vs_Mchirp","Efficiency_vs_Mchirp",0,0,1200,700);
1523  cS->Divide(2,2);
1524  //cS->cd(1)->SetLogy();
1525  cS->cd(1);
1526  TMultiGraph* mg1=new TMultiGraph();
1527  mg1->SetTitle("cc: no_cut;Mchirp;#Events");
1528  for (int h=0;h<nRHO;h++){
1529  gS[h]->SetLineColor(4);
1530  if(gS[h]->GetN()!=0) mg1->Add(gS[h]);
1531  }
1532  mg1->Draw("al");
1533  cS->cd(2);
1534  // cS->cd(2)->SetLogy();
1535  TMultiGraph* mg2=new TMultiGraph();
1536  mg2->SetTitle("cc=0.55;Mchirp;#Events");
1537  for (int h=0;h<nRHO;h++){
1538  gS[nRHO+h]->SetLineColor(4);
1539  if(gS[nRHO+h]->GetN()!=0) mg2->Add(gS[nRHO+h]);
1540  }
1541  mg2->Draw("al");
1542 
1543  //cS->cd(3)->SetLogy();
1544  cS->cd(3);
1545  TMultiGraph* mg3=new TMultiGraph();
1546  mg3->SetTitle("cc=0.6;Mchirp;#Events");
1547  for (int h=0;h<nRHO;h++){
1548  gS[2*nRHO+h]->SetLineColor(4);
1549  if(gS[2*nRHO+h]->GetN()!=0) mg3->Add(gS[2*nRHO+h]);
1550  }
1551  mg3->Draw("al");
1552  //cS->cd(4)->SetLogy();
1553  cS->cd(4);
1554  TMultiGraph* mg4=new TMultiGraph();
1555  mg4->SetTitle("cc=0.65;Mchirp;#Events");
1556  for (int h=0;h<nRHO;h++){
1557  gS[3*nRHO+h]->SetLineColor(4);
1558  if(gS[3*nRHO+h]->GetN()!=0) mg4->Add(gS[3*nRHO+h]);
1559  }
1560  mg4->Draw("al");
1561 
1562  TCanvas* cB=new TCanvas("Number_vs_Mchirp","Number_vs_Mchirp",0,0,1200,700);
1563  cB->Divide(2,2);
1564  cB->cd(1)->SetLogy();
1565  TMultiGraph* mg1B=new TMultiGraph();
1566  mg1B->SetTitle("cc: no_cut;Mchirp;#Events");
1567  for (int h=0;h<nRHO;h++){
1568  gB[h]->SetLineColor(4);
1569  if(gB[h]->GetN()!=0) mg1B->Add(gB[h]);
1570  }
1571  mg1B->Draw("al");
1572  cB->cd(2)->SetLogy();
1573  TMultiGraph* mg2B=new TMultiGraph();
1574  mg2B->SetTitle("cc=0.55;Mchirp;#Events");
1575  for (int h=0;h<nRHO;h++){
1576  gB[nRHO+h]->SetLineColor(4);
1577  if(gB[nRHO+h]->GetN()!=0) mg2B->Add(gB[nRHO+h]);
1578  }
1579  mg2B->Draw("al");
1580  cB->cd(3)->SetLogy();
1581  TMultiGraph* mg3B=new TMultiGraph();
1582  mg3B->SetTitle("cc=0.6;Mchirp;#Events");
1583  for (int h=0;h<nRHO;h++){
1584  gB[2*nRHO+h]->SetLineColor(4);
1585  if(gB[2*nRHO+h]->GetN()!=0) mg3B->Add(gB[2*nRHO+h]);
1586  }
1587  mg3B->Draw("al");
1588  cB->cd(4)->SetLogy();
1589  TMultiGraph* mg4B=new TMultiGraph();
1590  mg4B->SetTitle("cc=0.6;Mchirp;#Events");
1591  for (int h=0;h<nRHO;h++){
1592  gB[3*nRHO+h]->SetLineColor(4);
1593  if(gB[3*nRHO+h]->GetN()!=0) mg4B->Add(gB[3*nRHO+h]);
1594  }
1595  mg4B->Draw("al");
1596 
1597 
1598  TString CnameS(name);
1599  TString CnameB(name);
1600  TString CnameSroot(name);
1601  TString CnameBroot(name);
1602  char CnameS2[1024];
1603  char CnameB2[1024];
1604  char CnameS2root[1024];
1605  char CnameB2root[1024];
1606  CnameS.ReplaceAll(".root",".png");
1607  CnameB.ReplaceAll(".root",".png");
1608  sprintf(CnameS2,"Mcthres/N_Mc_S_%s",CnameS.Data());
1609  sprintf(CnameB2,"Mcthres/N_Mc_B_%s",CnameB.Data());
1610  sprintf(CnameS2root,"Mcthres/N_Mc_S_%s",CnameSroot.Data());
1611  sprintf(CnameB2root,"Mcthres/N_Mc_B_%s",CnameBroot.Data());
1612  cS->SaveAs(CnameS2);
1613  cB->SaveAs(CnameB2);
1614  cS->Print(CnameS2root);
1615  cB->Print(CnameB2root);
1616 
1617 
1618 }
1619 
1620 void PlotsAv_cc(TString ifile){
1621  TString name(ifile);
1622  name.ReplaceAll("outfile/","");
1623  name.ReplaceAll("average_file/","");
1624  TFile* fTEST =TFile::Open(ifile.Data());
1625  TTree* NNTree2=(TTree*)fTEST->Get("Parameters");
1626  //double outbis;
1627  double av;
1628  double out;
1629  double cc;
1630  double rho;
1631  int nSi;
1632 // NNTree2->SetBranchAddress("ANNout",&out);
1633  NNTree2->SetBranchAddress("Average",&av);
1634  NNTree2->SetBranchAddress("cc",&cc);
1635  NNTree2->SetBranchAddress("rho",&rho);
1636  NNTree2->SetBranchAddress("#TestSig",&nSi);
1637  NNTree2->GetEntry(0);
1638  int const nSig=nSi;
1639  cout<<"nSig: "<<nSig<<" nSi: "<<nSi<<endl;
1640  int const nBg=NNTree2->GetEntries()-nSig;
1641  cout<<"nBg: "<<nBg<<" Entries: "<<NNTree2->GetEntries()<<endl;
1642 
1643  TGraph * gS[3];
1644  TGraph * gB[3];
1645  gB[0]=new TGraph();
1646  gS[0]=new TGraph();
1647 // cout<<"nSig: "<<nSig<<endl;
1648  for (int n = 0; n <NNTree2->GetEntries(); n++){
1649  NNTree2->GetEntry(n);
1650  if(n<nSig) {
1651  gS[0]->SetPoint(n,cc,av);
1652  cout<<"Sig_graph1: x="<<cc<<" y: "<<av<<endl;
1653  }
1654  else {
1655  gB[0]->SetPoint(n-nSig,cc,av);
1656  cout<<"Bg_graph1: x="<<cc<<" y: "<<av<<endl;
1657  }
1658  }
1659 
1660 
1661  /*for (int a=0;a<nBg;a++){
1662  if(aRHOB[a]>=RHOth){
1663  for (int b=0;b<nth;b++){
1664  if(aCCB[a]>=cc1th[b]){
1665  for(int c=0;c<nth+1;c++){
1666  if(aANNB[a]>=ANN1th[c]) Z[b*nth+c]=Z[b*nth+c]+1;
1667  if(c=10) cout<<" ANN "<<ANN1th[c]<<" Zcount: "<<Z[b*nth+c]<<" aNNB "<<aANNB[a]<<" a "<<a<<" b "<<b<<endl;
1668  }
1669  }
1670  }
1671  }
1672  }*/
1673  //gS[0]->GetHistogram()->GetXaxis()->SetTitle("Mchirp");
1674  //gB[0]->GetHistogram()->GetYaxis()->SetTitle("Average on ANN ouputs");
1675  //gS[0]->GetHistogram()->GetXaxis()->SetTitle("Mchirp");
1676  //gB[0]->GetHistogram()->GetYaxis()->SetTitle("Average on ANN ouputs");
1677  gS[0]->SetMarkerColor(2);
1678  gB[0]->SetMarkerColor(4);
1679  gS[0]->SetMarkerStyle(6);
1680  gB[0]->SetMarkerStyle(7);
1681 
1682  TCanvas* c=new TCanvas("Plots","Plots",0,0,1200,700);
1683  c->Divide(1,2);
1684  c->cd(1);
1685  TMultiGraph* mg1=new TMultiGraph();
1686  mg1->SetTitle("Av_cc");
1687  if(gB[0]->GetN()!=0) mg1->Add(gB[0]);
1688  if(gS[0]->GetN()!=0) mg1->Add(gS[0]);
1689  mg1->Draw("ap");
1690  mg1->GetHistogram()->GetXaxis()->SetTitle("cc");
1691  mg1->GetHistogram()->GetYaxis()->SetTitle("Average on ANN ouputs");
1692 // mg1->Draw("ap");
1693  c->cd(2);
1694  TMultiGraph* mg2=new TMultiGraph();
1695  mg2->SetTitle("Av_cc");
1696  if(gS[0]->GetN()!=0) mg2->Add(gS[0]);
1697  if(gB[0]->GetN()!=0) mg2->Add(gB[0]);
1698  mg2->Draw("ap");
1699  mg2->GetHistogram()->GetXaxis()->SetTitle("cc");
1700  mg2->GetHistogram()->GetYaxis()->SetTitle("Average on ANN ouputs");
1701 
1702  cout<<" name "<<name<<endl;
1703  TString Cname(name);
1704  TString Cnameroot(name);
1705  char Cname2[1024];
1706  char Cname2root[1024];
1707  Cname.ReplaceAll(".root",".png");
1708  sprintf(Cname2,"average_png/out_Plots_%s",Cname.Data());
1709  sprintf(Cname2root,"average_png/out_Plots_%s",Cnameroot.Data());
1710  c->SaveAs(Cname2);
1711  c->Print(Cname2root);
1712 /*
1713  TCanvas* c2=new TCanvas("Plots_Bkg_on_Sig","Plots_Bkg_on_Sig",0,0,1200,700);
1714  c2->Divide(2,2);
1715  c2->cd(1);
1716  gS[0]->Draw("ap");
1717  gB[0]->Draw("p,same");
1718 
1719  c2->cd(2);
1720  gS[1]->Draw("ap");
1721  gB[1]->Draw("p,same");
1722  c2->cd(3);
1723  gS[2]->Draw("ap");
1724  gB[2]->Draw("p,same");
1725  c2->cd(4);
1726  TText* text2=new TText(0.37,0.0,"no cuts on ANN");
1727  hglitch->SetStats(0);
1728  hglitch->GetXaxis()->SetTitle("cc");
1729  hglitch->GetYaxis()->SetTitle("ANNoutput");
1730  hglitch->GetZaxis()->SetTitle("count");
1731  hglitch->Draw("colz");
1732  text2->Draw();
1733  gPad->SetLogz();
1734 
1735  TString Cname_2(name);
1736  TString Cname_2root(name);
1737  char Cname2_2[1024];
1738  char Cname2_2root[1024];
1739  Cname_2.ReplaceAll(".root",".png");
1740  sprintf(Cname2_2,"CC_RHO_ANN_Plots/CC_RHO_ANN_Plots_BoS_%s",Cname_2.Data());
1741  sprintf(Cname2_2root,"CC_RHO_ANN_Plots/CC_RHO_ANN_Plots_BoS_%s",Cname_2root.Data());
1742  c2->SaveAs(Cname2_2);
1743  c2->Print(Cname2_2root);
1744 
1745 //CARTELLA CCi_RHO_ANN_Plots */
1746 }
1747 void PlotsAv_Mc(TString ifile){
1748  TString name(ifile);
1749  name.ReplaceAll("outfile/","");
1750  name.ReplaceAll("average_file/","");
1751  TFile* fTEST =TFile::Open(ifile.Data());
1752  TTree* NNTree2=(TTree*)fTEST->Get("Parameters");
1753  //double outbis;
1754  double av;
1755  double out;
1756  double cc;
1757  double Mc;
1758  double Mc_inj;
1759  double rho;
1760  int nSi;
1761 // NNTree2->SetBranchAddress("ANNout",&out);
1762  NNTree2->SetBranchAddress("Mchirp",&Mc);
1763  NNTree2->SetBranchAddress("Average",&av);
1764  NNTree2->SetBranchAddress("cc",&cc);
1765  NNTree2->SetBranchAddress("rho",&rho);
1766  NNTree2->SetBranchAddress("Mchirp_injected",&Mc_inj);
1767  NNTree2->SetBranchAddress("#TestSig",&nSi);
1768  NNTree2->GetEntry(0);
1769  int const nSig=nSi;
1770  cout<<"nSig: "<<nSig<<" nSi: "<<nSi<<endl;
1771  //exit(0);
1772  int const nBg=NNTree2->GetEntries()-nSig;
1773 
1774  TGraph * gS[3];
1775  TGraph * gB[3];
1776  gB[0]=new TGraph();
1777  gS[0]=new TGraph();
1778 
1779 
1780  TCanvas* ANN_Mc_c3D=new TCanvas("ANN_average_vs_Mc_reconstructed","ANN_average_vs_Mc_reconstructed",0,0,1200,700);
1781  ANN_Mc_c3D->Divide(1,2);
1782  TH2D* ANN_Mc_SIG3D=new TH2D("ANN_average_vs_Mc_recostructed_sig","ANN_average_vs_Mc_recostructed_sig",NMc_r,minMc_r,maxMc_r,NANN_av,minANN_av,maxANN_av);
1783  TH2D* ANN_Mc_BKG3D=new TH2D("ANN_average_vs_Mc_recostructed_bkg","ANN_average_vs_Mc_recostructed_bkg",NMc_r,minMc_r,maxMc_r,NANN_av,minANN_av,maxANN_av);
1784  // Mc_inj_g->SetPoint(0,rhoSig[y][ii],yy);
1785  ANN_Mc_BKG3D->SetDirectory(0);
1786  ANN_Mc_BKG3D->SetStats(0);
1787  ANN_Mc_SIG3D->SetDirectory(0);
1788  ANN_Mc_SIG3D->SetStats(0);
1789 
1790  TCanvas* ANN_cc_c3D=new TCanvas("ANN_average_vs_cc","ANN_average_vs_cc",0,0,1200,700);
1791  ANN_cc_c3D->Divide(1,2);
1792  TH2D* ANN_cc_SIG3D=new TH2D("ANN_average_vs_cc_SIG","ANN_average_vs_cc_SIG",Ncc3D, minCC, maxCC,NANN_av,minANN_av,maxANN_av);
1793  TH2D* ANN_cc_BKG3D=new TH2D("ANN_average_vs_cc_BKG","ANN_average_vs_cc_BKG",Ncc3D, minCC, maxCC,NANN_av,minANN_av,maxANN_av);
1794  ANN_cc_SIG3D->SetDirectory(0);
1795  ANN_cc_SIG3D->SetStats(0);
1796  ANN_cc_BKG3D->SetDirectory(0);
1797  ANN_cc_BKG3D->SetStats(0);
1798  double deltacc3D=0;
1799  deltacc3D=(maxCC-minCC)/(Ncc3D);
1800 
1801 
1802  TCanvas* Mc_inj_c=new TCanvas("ANN_vs_Mc_injected","ANN_vs_Mc_injected",0,0,1200,700);
1803  Mc_inj_c->Divide(1,2);
1804  TH2D* Mc_inj_gMc=new TH2D("Mc_injected_Mc_recostructed","Mc_injected_recostructed",NMc_inj,minMc_inj,maxMc_inj,NMc_r,minMc_r,maxMc_r);
1805  TH2D* Mc_inj_g=new TH2D("Mc_injected_comparison","Mc_injected_comparison",NMc_inj,minMc_inj,maxMc_inj,NANN_av,minANN_av,maxANN_av);
1806  // Mc_inj_g->SetPoint(0,rhoSig[y][ii],yy);
1807  Mc_inj_gMc->SetDirectory(0);
1808  Mc_inj_gMc->SetStats(0);
1809  Mc_inj_g->SetDirectory(0);
1810  Mc_inj_g->SetStats(0);
1811  double deltaMc_inj=0;
1812  deltaMc_inj=(maxMc_inj-minMc_inj)/(NMc_inj);
1813  double deltaANN_av_M=0;
1814  deltaANN_av_M=(maxANN_av-minANN_av)/(NANN_av);
1815  double deltaMc_r=0;
1816  deltaMc_r=(maxMc_r-minMc_r)/(NMc_r);
1817 
1818 
1819  TString txt_name0(name);
1820  txt_name0.ReplaceAll(".root",".txt");
1821  char txt_name[1024];
1822  sprintf(txt_name,"txt_files/%s",txt_name0.Data());
1823  ofstream file_txt_out(txt_name);
1824 
1825 
1826 // cout<<"nSig: "<<nSig<<endl;
1827  int aa1=0;
1828  int bb1=0;
1829  for (int n = 0; n <NNTree2->GetEntries(); n++){
1830  NNTree2->GetEntry(n);
1831  char line[2048];
1832  if(n<nSig&&n!=nSig) {
1833  Mc_inj_gMc->Fill(Mc_inj,Mc);
1834  Mc_inj_g->Fill(Mc_inj,av);
1835  gS[0]->SetPoint(n,Mc,av);
1836  cout<<" n "<<n<<" Mc "<<Mc<<" av "<<av<<" M_inj "<<Mc_inj<<endl;
1837  cout<<"Sig_graph1: x="<<cc<<" y: "<<av<<endl;
1838  aa1=aa1+1;
1839  ANN_Mc_SIG3D->Fill(Mc,av);
1840  ANN_cc_SIG3D->Fill(cc,av);
1841  sprintf(line,"SIG:Average_%1.5f; Mc_rec_%1.5f; Mc_inj_%1.5f; event %i; count_%i; nSig_%i\n",av,Mc,Mc_inj,n,aa1,nSig);
1842  file_txt_out<<line;
1843  cout<<" n "<<n<<" nSig "<<nSig<<endl;
1844  }
1845  else {
1846  bb1=bb1+1;
1847  gB[0]->SetPoint(n-nSig,Mc,av);
1848  ANN_Mc_BKG3D->Fill(Mc,av);
1849  ANN_cc_BKG3D->Fill(cc,av);
1850  // cout<<"Bg_graph1: x="<<cc<<" y: "<<av<<endl;
1851  sprintf(line,"BKG:Average %1.5f; Mc_rec %1.5f; Mc_inj_%1.5f; event %i; count_%i\n",av,Mc,Mc_inj,n,bb1);
1852  file_txt_out<<line;
1853 
1854  }
1855  }
1856  //exit(0);
1857 
1858 
1859  gS[0]->SetMarkerColor(2);
1860  gB[0]->SetMarkerColor(4);
1861  gS[0]->SetMarkerStyle(6);
1862  gB[0]->SetMarkerStyle(7);
1863 
1864  ANN_Mc_c3D->cd(1);
1865  ANN_Mc_SIG3D->GetXaxis()->SetTitle("SIG:Mchirp_reconstructed");
1866  ANN_Mc_SIG3D->GetYaxis()->SetTitle("SIG:ANN_average");
1867  ANN_Mc_SIG3D->GetZaxis()->SetTitle("count");
1868  ANN_Mc_SIG3D->Draw("colz");
1869  ANN_Mc_c3D->cd(2);
1870  ANN_Mc_BKG3D->GetXaxis()->SetTitle("BKG:Mchirp_reconstructed");
1871  ANN_Mc_BKG3D->GetYaxis()->SetTitle("BKG:ANN_average");
1872  ANN_Mc_BKG3D->GetZaxis()->SetTitle("count");
1873  ANN_Mc_BKG3D->Draw("colz");
1874 
1875 
1876  TString ANN_Mc_name3D(name);
1877  TString ANN_Mc_rootname3D(name);
1878  char ANN_Mc_cname23D[1024];
1879  char ANN_Mc_crootname23D[1024];
1880  ANN_Mc_name3D.ReplaceAll(".root",".png");
1881  sprintf(ANN_Mc_cname23D,"average_png/nSig%i_nBg%i_3DMc_ANNav_Plots_%s",nSig,nBg,ANN_Mc_name3D.Data());
1882  sprintf(ANN_Mc_crootname23D,"average_png/nSig%i_nBg%i_3DMc_ANNav_Plots_%s",nSig,nBg,ANN_Mc_rootname3D.Data());
1883  ANN_Mc_c3D->SaveAs(ANN_Mc_cname23D);
1884  ANN_Mc_c3D->Print(ANN_Mc_crootname23D);
1885 
1886 
1887  ANN_cc_c3D->cd(1);
1888  ANN_cc_SIG3D->GetXaxis()->SetTitle("SIG:cc");
1889  ANN_cc_SIG3D->GetYaxis()->SetTitle("SIG:ANN_average");
1890  ANN_cc_SIG3D->GetZaxis()->SetTitle("count");
1891  ANN_cc_SIG3D->Draw("colz");
1892  ANN_cc_c3D->cd(2);
1893  ANN_cc_BKG3D->GetXaxis()->SetTitle("BKG:cc");
1894  ANN_cc_BKG3D->GetYaxis()->SetTitle("BKG:ANN_average");
1895  ANN_cc_BKG3D->GetZaxis()->SetTitle("count");
1896  ANN_cc_BKG3D->Draw("colz");
1897 
1898 
1899  TString ANN_cc_cname3D(name);
1900  TString ANN_cc_crootname3D(name);
1901  char ANN_cc_cname23D[1024];
1902  char ANN_cc_crootname23D[1024];
1903  ANN_cc_cname3D.ReplaceAll(".root",".png");
1904  sprintf(ANN_cc_cname23D,"average_png/nSig%i_nBg%i_3DANNav_cc_Plots_%s",nSig,nBg,ANN_cc_cname3D.Data());
1905  sprintf(ANN_cc_crootname23D,"average_png/nSig%i_nBg%i_3DANNav_cc_Plots_%s",nSig,nBg,ANN_cc_crootname3D.Data());
1906  ANN_cc_c3D->SaveAs(ANN_cc_cname23D);
1907  ANN_cc_c3D->Print(ANN_cc_crootname23D);
1908 
1909  Mc_inj_c->cd(1);
1910  Mc_inj_g->GetXaxis()->SetTitle("Mchirp_injected");
1911  Mc_inj_g->GetYaxis()->SetTitle("ANN_average");
1912  Mc_inj_g->GetZaxis()->SetTitle("count");
1913  Mc_inj_g->Draw("colz");
1914  Mc_inj_c->cd(2);
1915  Mc_inj_gMc->GetXaxis()->SetTitle("Mchirp_injected");
1916  Mc_inj_gMc->GetYaxis()->SetTitle("Mchirp_reconstructed");
1917  Mc_inj_gMc->GetZaxis()->SetTitle("count");
1918  Mc_inj_gMc->Draw("colz");
1919 
1920 
1921  TString Mc_inj_cname(name);
1922  TString Mc_inj_crootname(name);
1923  char Mc_inj_cname2[1024];
1924  char Mc_inj_crootname2[1024];
1925  Mc_inj_cname.ReplaceAll(".root",".png");
1926  sprintf(Mc_inj_cname2,"average_png/nSig%i_nBg%i_Mc_inj_Plots_%s",nSig,nBg,Mc_inj_cname.Data());
1927  sprintf(Mc_inj_crootname2,"average_png/nSig%i_nBg%i_Mc_inj_Plots_%s",nSig,nBg,Mc_inj_crootname.Data());
1928  Mc_inj_c->SaveAs(Mc_inj_cname2);
1929  Mc_inj_c->Print(Mc_inj_crootname2);
1930 
1931 
1932  TCanvas* c=new TCanvas("Plots","Plots",0,0,1200,700);
1933  c->Divide(1,2);
1934  c->cd(1);
1935  TMultiGraph* mg1=new TMultiGraph();
1936  mg1->SetTitle("Av_Mc");
1937  if(gB[0]->GetN()!=0) mg1->Add(gB[0]);
1938  if(gS[0]->GetN()!=0) mg1->Add(gS[0]);
1939  mg1->Draw("ap");
1940  mg1->GetHistogram()->GetXaxis()->SetTitle("Mchirp");
1941  mg1->GetHistogram()->GetYaxis()->SetTitle("Average on ANN ouputs");
1942  c->cd(2);
1943  TMultiGraph* mg2=new TMultiGraph();
1944  mg2->SetTitle("Av_Mc");
1945  if(gS[0]->GetN()!=0) mg2->Add(gS[0]);
1946  if(gB[0]->GetN()!=0) mg2->Add(gB[0]);
1947  mg2->Draw("ap");
1948  mg2->GetHistogram()->GetXaxis()->SetTitle("Mchirp");
1949  mg2->GetHistogram()->GetYaxis()->SetTitle("Average on ANN ouputs");
1950 
1951  cout<<" name "<<name<<endl;
1952  TString Cname(name);
1953  TString Cnameroot(name);
1954  char Cname2[1024];
1955  char Cname2root[1024];
1956  Cname.ReplaceAll(".root",".png");
1957  sprintf(Cname2,"average_png/nSig%i_nBg%i_Mc_Plots_%s",nSig,nBg,Cname.Data());
1958  sprintf(Cname2root,"average_png/nSig%i_nBg%i_Mc_Plots_%s",nSig,nBg,Cnameroot.Data());
1959  cout<<"Cname2 "<<Cname2<<endl;
1960  c->SaveAs(Cname2);
1961  c->Print(Cname2root);
1962  cout<<" gS[0]->GetN() "<<gS[0]->GetN()<<endl;
1963  //cout<<"entries "<<mg1->GetHistogram()->GetEntries()<<endl;
1964  cout<<" nSig "<<nSig<<" nBg "<<nBg<<" aa1 "<<aa1<<endl;
1965 }
1966 
1967 
1968 void SNR_plots(TString ifile){
1969 
1970  TString name(ifile);
1971  name.ReplaceAll("outfile/","");
1972  name.ReplaceAll("average_file/","");
1973  TFile* fTEST =TFile::Open(ifile.Data());
1974  TTree* NNTree2=(TTree*)fTEST->Get("Parameters");
1975  double av;
1976  double out;
1977  double cc;
1978  double Mc;
1979  double Mc_inj;
1980  double rho;
1981  int nSi;
1982  NNTree2->SetBranchAddress("Mchirp",&Mc);
1983  NNTree2->SetBranchAddress("Average",&av);
1984  NNTree2->SetBranchAddress("cc",&cc);
1985  NNTree2->SetBranchAddress("rho",&rho);
1986  NNTree2->SetBranchAddress("Mchirp_injected",&Mc_inj);
1987  NNTree2->SetBranchAddress("#TestSig",&nSi);
1988  NNTree2->GetEntry(0);
1989  int const nSig=nSi;
1990  cout<<"nSig: "<<nSig<<" nSi: "<<nSi<<endl;
1991  int const nBg=NNTree2->GetEntries()-nSig;
1992  int const nTot=NNTree2->GetEntries();
1993  int const n_points=N_ANNth*N_RHOth;
1994  int nBgSur[3][n_points];
1995  int nSigSur[3][n_points];
1996  double zax[3][n_points];
1997 
1998  for (int i=0; i<3; i++) for (int u=0; u<n_points; u++) {nBgSur[i][u]=0; nSigSur[i][u]=0;zax[i][u]=0.;}
1999 
2000  double pANNth[N_ANNth];
2001  double pRHOth[N_RHOth];
2002 
2003  double pCCth[3];
2004 
2005  pCCth[0]=0.;
2006  pCCth[1]=0.;
2007  pCCth[2]=0.;
2008 
2009  pCCth[0]=0.5;
2010  pCCth[1]=0.6;
2011  pCCth[2]=0.7;
2012 
2013  double deltaANNth=0.;
2014  double deltaRHOth=0.;
2015 
2016  deltaANNth=(double)(ANNth_max-ANNth_min)/N_ANNth;
2017  deltaRHOth=(double)(RHOth_max-RHOth_min)/N_RHOth;
2018 
2019 
2020  TCanvas* SNR_c=new TCanvas("SIG-BKG_reduction","SIG-BKG_reduction",0,0,1200,700);
2021  SNR_c->Divide(1,3);
2022 
2023  TH2D* SNR3Dcc0=new TH2D("cc=0.5","cc=0.5",N_ANNth,ANNth_min,ANNth_max,N_RHOth,RHOth_min,RHOth_max);
2024  TH2D* SNR3Dcc1=new TH2D("cc=0.6","cc=0.6",N_ANNth,ANNth_min,ANNth_max,N_RHOth,RHOth_min,RHOth_max);
2025  TH2D* SNR3Dcc2=new TH2D("cc=0.7","cc=0.7",N_ANNth,ANNth_min,ANNth_max,N_RHOth,RHOth_min,RHOth_max);
2026  SNR3Dcc0->SetDirectory(0);
2027  SNR3Dcc0->SetStats(0);
2028  SNR3Dcc1->SetDirectory(0);
2029  SNR3Dcc1->SetStats(0);
2030  SNR3Dcc2->SetDirectory(0);
2031  SNR3Dcc2->SetStats(0);
2032 
2033  for (int u=0; u<N_ANNth; u++) { pANNth[u]=0.; pANNth[u]=ANNth_min+u*deltaANNth; cout<<" pANNth[u] "<<pANNth[u]<<" u "<<u<<endl;}
2034  for (int u=0; u<N_RHOth; u++) { pRHOth[u]=0.; pRHOth[u]=RHOth_min+u*deltaRHOth; cout<<" pRHOth[u] "<<pRHOth[u]<<" u "<<u<<endl;}
2035 
2036 cout<<"fino def cose"<<endl;
2037 /*
2038  for (int u=0; u<nTot; u++){
2039  NNTree2->GetEntry(u);
2040  cout<<" cc "<<cc<<" av "<<av<<" rho "<< rho<<endl;
2041  for (int y=2; y>-1; y--){
2042  // cout<<" cc "<<cc<<" pCCth "<<pCCth[y]<<" y "<<y<<endl;
2043  if(cc>pCCth[y]) {
2044  for (int i=N_ANNth-1; i>-1; i--){
2045  // cout<<" av "<<av<<" pANNth[i] "<<pANNth[i]<<" i "<<i<<endl;
2046  if(av>pANNth[i]){
2047  for (int j=N_RHOth-1; j>-1; j++){
2048  //cout<<" rho "<<rho<<" pRHOth[j] "<<pRHOth[j]<<" j "<<j<<" nSig "<<nSig<<" u "<<u<<endl;
2049  if(rho>pRHOth[j]){
2050  for(int yy=y; yy>-1;yy--) for (int ii=i; ii>-1; ii--) for (int jj=j; jj>-1;jj--){
2051  if(u<nSig) nSigSur[yy][jj+N_ANNth*ii]=nSigSur[yy][jj+N_ANNth*ii]+1;
2052  else nBgSur[yy][jj+N_ANNth*ii]=nBgSur[yy][jj+N_ANNth*ii]+1;
2053  // cout<<" yy "<<yy<<" ii "<<ii<<" jj "<<jj<<" nBgSur[yy][jj+N_ANNth*ii] "<<nBgSur[yy][jj+N_ANNth*ii]<<" nSigSur[yy][jj+N_ANNth*ii] "<<nSigSur[yy][jj+N_ANNth*ii]<<endl;
2054  }
2055  break;
2056  }
2057  else continue;
2058  }
2059  break;
2060  }
2061  else continue;
2062  }
2063  break;
2064  }
2065  else continue;
2066  }
2067  }
2068 //exit(0);
2069 cout<<"dopo ciclo"<<endl;
2070 */
2071 for (int u=0; u<nTot; u++){
2072  NNTree2->GetEntry(u);
2073  for (int y=0; y<3; y++){
2074  if(cc>pCCth[y]){
2075  for (int i=0; i<N_ANNth; i++){
2076  if(av>pANNth[i]) {
2077  for (int j=0; j<N_RHOth; j++){
2078  if(rho>pRHOth[j]){
2079  if(u<nSig) nSigSur[y][j+N_ANNth*i]=nSigSur[y][j+N_ANNth*i]+1;
2080  else nBgSur[y][j+N_ANNth*i]=nBgSur[y][j+N_ANNth*i]+1;
2081  }
2082  else continue;
2083  }
2084  }
2085  else continue;
2086  }
2087  }
2088  else continue;
2089  }
2090  }
2091 
2092 
2093 //exit(0);
2094 int nBgSur0=0;
2095 nBgSur0=nBgSur[0][0];
2096 int nSigSur0=0;
2097 nSigSur0=nSigSur[0][0];
2098 if(nBgSur0==0){nBgSur0=nBg;}
2099 if(nSigSur0==0){nSigSur0=nSig;}
2100 int count_zax=0;
2101 cout<<"dopo ciclo"<<" nSigSur0 "<<nSigSur0<<" nBgSur0 "<<nBgSur0<<endl;
2102 for (int u=0; u<3; u++) for (int j=0; j<N_RHOth; j++) for (int i=0; i<N_ANNth; i++) {
2103  //cout<<" u "<<u<<" j "<<j<<" i "<<i<<" nSigSur[u][i+j*N_ANNth] "<<nSigSur[u][i+j*N_ANNth]<<endl;
2104  //if(nSigSur[u][i+j*N_ANNth]!=0){zax[u][i+j*N_ANNth]=(nBgSur[u][i+j*N_ANNth]/nBgSur0)/(nSigSur[u][i+j*N_ANNth]/nSigSur0);};
2105  if(nSigSur[u][i+j*N_ANNth]==0) {
2106  if (nSigSur[u][i+j*N_ANNth]==0 && nBgSur[u][i+j*N_ANNth]==0) zax[u][i+j*N_ANNth]=-1;
2107  else zax[u][i+j*N_ANNth]=-2;
2108  }
2109  else {
2110  double num=0.;
2111  double den=0.;
2112  num=(double)nBgSur[u][i+j*N_ANNth]/nBgSur0;
2113  den=(double)nSigSur0/nSigSur[u][i+j*N_ANNth];
2114  zax[u][i+j*N_ANNth]=(double)num/den; count_zax=count_zax+1;
2115  cout<<(double)(nBgSur[u][i+j*N_ANNth]/nBgSur0)*(nSigSur0/nSigSur[u][i+j*N_ANNth])<<endl;
2116  }
2117 // cout<<" nBgSur[yy][jj+N_ANNth*ii] "<<nBgSur[u][j+N_ANNth*i]<<" nSigSur[y][j+N_ANNth*i] "<<nSigSur[u][j+N_ANNth*i]<<" zax[u][i+j*N_ANNth] "<<zax[u][i+j*N_ANNth]<<endl;
2118  //else zax[u][i+j*N_ANNth]=(nBgSur[u][i+j*N_ANNth]/nBgSur0)/(nSigSur[u][i+j*N_ANNth]/nSigSur0);
2119  if(u==2){
2120  /*cout<<" zax[0][i+j*N_ANNth] "<<zax[0][i+j*N_ANNth]<<endl;
2121  cout<<" zax[1][i+j*N_ANNth] "<<zax[1][i+j*N_ANNth]<<endl;
2122  cout<<" zax[2][i+j*N_ANNth] "<<zax[2][i+j*N_ANNth]<<endl;*/
2123  SNR3Dcc0->Fill(pANNth[i]+ deltaANNth/2.,pRHOth[j]+ deltaRHOth/2.,zax[0][i+j*N_ANNth]);
2124  SNR3Dcc1->Fill(pANNth[i]+ deltaANNth/2.,pRHOth[j]+ deltaRHOth/2.,zax[1][i+j*N_ANNth]);
2125  SNR3Dcc2->Fill(pANNth[i]+ deltaANNth/2.,pRHOth[j]+ deltaRHOth/2.,zax[2][i+j*N_ANNth]);
2126  }
2127 }
2128 for(int yy=0; yy<3;yy++) {
2129  cout<<" pCCth "<<pCCth[yy]<<endl;
2130  for (int ii=N_ANNth-1; ii>-1; ii--) {
2131  cout<<" pANNth[i] "<<pANNth[ii]<<endl;
2132  for (int jj=N_RHOth-1; jj>-1;jj--){
2133  cout<<" pRHOth[jj] "<<pRHOth[jj]<<endl;
2134  cout<<" nBgSur[yy][ii+N_ANNth*jj] "<<nBgSur[yy][ii+N_ANNth*jj]<<endl;
2135  cout<<" nSigSur[yy][ii+N_ANNth*jj] "<<nSigSur[yy][ii+N_ANNth*jj]<<endl;
2136  cout<<" zax[0][i+j*N_ANNth]"<<zax[yy][ii+jj*N_ANNth]<<endl;
2137  cout<<" (double)(nBgSur[u][i+j*N_ANNth]/nBgSur0)"<<(double)nBgSur[yy][ii+jj*N_ANNth]/nBgSur0<<endl;
2138  //cout<<" yy "<<yy<<" ii "<<ii<<" jj "<<jj<<" nBgSur[yy][jj+N_ANNth*ii] "<<nBgSur[yy][jj+N_ANNth*ii]<<" nSigSur[yy][jj+N_ANNth*ii] "<<nSigSur[yy][jj+N_ANNth*ii]<<endl;
2139  }
2140  }
2141  }
2142 cout<<"coutn_zax "<<count_zax<<endl;
2143  SNR_c->cd(1);
2144  SNR3Dcc0->GetXaxis()->SetTitle("ANN_average_thresholds");
2145  SNR3Dcc0->GetYaxis()->SetTitle("RHO_thresholds");
2146  SNR3Dcc0->GetZaxis()->SetTitle("normBKG/normSIG");
2147  SNR3Dcc0->Draw("colz");
2148 
2149  SNR_c->cd(2);
2150  SNR3Dcc1->GetXaxis()->SetTitle("ANN_average_thresholds");
2151  SNR3Dcc1->GetYaxis()->SetTitle("RHO_thresholds");
2152  SNR3Dcc1->GetZaxis()->SetTitle("normBKG/normSIG");
2153  SNR3Dcc1->Draw("colz");
2154 
2155  SNR_c->cd(3);
2156  SNR3Dcc2->GetXaxis()->SetTitle("ANN_average_thresholds");
2157  SNR3Dcc2->GetYaxis()->SetTitle("RHO_thresholds");
2158  SNR3Dcc2->GetZaxis()->SetTitle("normBKG/normSIG");
2159  SNR3Dcc2->Draw("colz");
2160 
2161  TString SNR_n(name);
2162  TString SNR_nroot(name);
2163  char SNR_c_n[1024];
2164  char SNR_c_nroot[1024];
2165  SNR_n.ReplaceAll(".root",".png");
2166  sprintf(SNR_c_n,"SNR_plots/nSig%i_nBg%i_invSNR_Plots_%s",nSig,nBg,SNR_n.Data());
2167  sprintf(SNR_c_nroot,"SNR_plots/nSig%i_nBg%i_invSNR_Plots_%s",nSig,nBg,SNR_nroot.Data());
2168  SNR_c->SaveAs(SNR_c_n);
2169  SNR_c->Print(SNR_c_nroot);
2170 
2171 
2172 }
2173 
2174 void CoG_plots(TString ifile){
2175 
2176  TString name(ifile);
2177  name.ReplaceAll("outfile/","");
2178  name.ReplaceAll("average_file/","");
2179  TFile* fTEST =TFile::Open(ifile.Data());
2180  TTree* NNTree2=(TTree*)fTEST->Get("Parameters");
2181  double av;
2182  double out;
2183  double cc;
2184  double Mc;
2185  double Mc_inj;
2186  double rho;
2187  int nSi;
2188  int CoG;
2189  NNTree2->SetBranchAddress("Center_of_Gravity",&CoG);
2190  NNTree2->SetBranchAddress("Mchirp",&Mc);
2191  NNTree2->SetBranchAddress("Average",&av);
2192  NNTree2->SetBranchAddress("cc",&cc);
2193  NNTree2->SetBranchAddress("rho",&rho);
2194  NNTree2->SetBranchAddress("Mchirp_injected",&Mc_inj);
2195  NNTree2->SetBranchAddress("#TestSig",&nSi);
2196  NNTree2->GetEntry(0);
2197  int const nSig=nSi;
2198  cout<<"nSig: "<<nSig<<" nSi: "<<nSi<<endl;
2199  int const nTot=NNTree2->GetEntries();
2200  int const nBg=nTot-nSig;
2201 
2202  TCanvas* CoG_c=new TCanvas("SIG-BKG_reduction","SIG-BKG_reduction",0,0,1200,700);
2203 
2204  TH1D* imp=new TH1D("Center_of_Gravity_VS_ANN","Center_of_Gravity_VS_ANN",N_ANNth,ANNth_min,ANNth_max);
2205  TH1D* CoG1Sig=new TH1D("Center_of_Gravity_VS_ANN","Center_of_Gravity_VS_ANN",N_ANNth,ANNth_min,ANNth_max);
2206  TH1D* CoG1Bg=new TH1D("Center_of_Gravity_VS_ANN","Center_of_Gravity_VS_ANN",N_ANNth,ANNth_min,ANNth_max);
2207  TH1D* CoG0Sig=new TH1D("Center_of_Gravity_VS_ANN","Center_of_Gravity_VS_ANN",N_ANNth,ANNth_min,ANNth_max);
2208  TH1D* CoG0Bg=new TH1D("Center_of_Gravity_VS_ANN","Center_of_Gravity_VS_ANN",N_ANNth,ANNth_min,ANNth_max);
2209 
2210  CoG1Sig->SetDirectory(0);
2211  CoG1Sig->SetStats(0);
2212  CoG0Sig->SetDirectory(0);
2213  CoG0Sig->SetStats(0);
2214  CoG1Bg->SetDirectory(0);
2215  CoG1Bg->SetStats(0);
2216  CoG0Bg->SetDirectory(0);
2217  CoG0Bg->SetStats(0);
2218 
2219  for (int n=0; n<nTot; n++){
2220  NNTree2->GetEntry(n);
2221  if (n<nSig) {
2222  if(CoG==0) {CoG0Sig->Fill(av);}
2223  else {CoG1Sig->Fill(av);}
2224  }
2225  else {
2226  if(CoG==0) {CoG0Bg->Fill(av);}
2227  else {CoG1Bg->Fill(av);}
2228  }
2229  cout<<" CoG "<<CoG<<endl;
2230  }
2231 
2232  imp->SetDirectory(0);
2233  imp->SetStats(0);
2234  int max_n=0;
2235  // if(nBg>nSig) max_n=nBg;
2236  // else max_n=nSig;
2237 
2238  int max[4];
2239  for(int i=0; i<4; i++) max[i]=0;
2240 
2241  max[0]=CoG1Sig->GetMaximum();
2242  max[1]=CoG0Sig->GetMaximum();
2243  max[2]=CoG1Bg->GetMaximum();
2244  max[3]=CoG0Bg->GetMaximum();
2245 
2246  max_n=max[0];
2247  for(int i=0; i<4; i++) if(max_n<max[i]) max_n=max[i] ;
2248 
2249  imp->Fill(ANNth_min,max_n);
2250  imp->SetLineColor(10);
2251  imp->SetBins(N_ANNth,ANNth_min,ANNth_max);
2252  imp->GetXaxis()->SetTitle("ANN");
2253  imp->GetYaxis()->SetTitle("Count");
2254  imp->Draw();
2255 
2256  CoG1Sig->SetFillStyle(3018);
2257  CoG0Sig->SetFillStyle(3018);
2258  CoG0Bg->SetFillStyle(3004);
2259  CoG1Bg->SetFillStyle(3004);
2260  CoG1Sig->SetFillColor(2);
2261  CoG0Sig->SetFillColor(6);
2262  CoG0Bg->SetFillColor(4);
2263  CoG1Bg->SetFillColor(9);
2264  CoG1Sig->SetLineColor(2);
2265  CoG0Sig->SetLineColor(6);
2266  CoG0Bg->SetLineColor(4);
2267  CoG1Bg->SetLineColor(9);
2268  // CoG1Sig->SetBins(N_ANNth,ANNth_min,ANNth_max);
2269  // CoG1Sig->GetXaxis()->SetTitle("ANN");
2270  // CoG1Sig->GetYaxis()->SetTitle("Count");
2271  CoG1Sig->Draw("same");
2272  CoG0Sig->Draw("same");
2273  CoG0Bg->Draw("same");
2274  CoG1Bg->Draw("same");
2275 
2276  TString CoG_n(name);
2277  TString CoG_nroot(name);
2278  char CoG_c_n[1024];
2279  char CoG_c_nroot[1024];
2280  CoG_n.ReplaceAll(".root",".png");
2281  sprintf(CoG_c_n,"CoG_plots/nSig%i_nBg%i_invSNR_Plots_%s",nSig,nBg,CoG_n.Data());
2282  sprintf(CoG_c_nroot,"CoG_plots/nSig%i_nBg%i_invSNR_Plots_%s",nSig,nBg,CoG_nroot.Data());
2283  CoG_c->SaveAs(CoG_c_n);
2284  CoG_c->Print(CoG_c_nroot);
2285 
2286 }
2287 
2288 void CutMcANN(TString ifile){
2289 
2290  TString name(ifile);
2291  name.ReplaceAll("outfile/","");
2292  name.ReplaceAll("average_file/","");
2293  TFile* fTEST =TFile::Open(ifile.Data());
2294  TTree* NNTree2=(TTree*)fTEST->Get("Parameters");
2295  double av;
2296  double out;
2297  double cc;
2298  double Mc;
2299  double Mc_inj;
2300  double rho;
2301  int nSi;
2302  int CoG;
2303  NNTree2->SetBranchAddress("Center_of_Gravity",&CoG);
2304  NNTree2->SetBranchAddress("Mchirp",&Mc);
2305  NNTree2->SetBranchAddress("Average",&av);
2306  NNTree2->SetBranchAddress("cc",&cc);
2307  NNTree2->SetBranchAddress("rho",&rho);
2308  NNTree2->SetBranchAddress("Mchirp_injected",&Mc_inj);
2309  NNTree2->SetBranchAddress("#TestSig",&nSi);
2310  NNTree2->GetEntry(0);
2311  int const nSig=nSi;
2312  cout<<"nSig: "<<nSig<<" nSi: "<<nSi<<endl;
2313  int const nTot=NNTree2->GetEntries();
2314  int const nBg=nTot-nSig;
2315 
2316 
2317 TGraph *Mc_g[nANNth];
2318 TGraph *ANNMc_g[nANNth];
2319 TGraph * rho_g[7];
2320 
2321 double Mccuts[7];
2322 double ANNcuts[7];
2323 double TA_rho[7][npoints];
2324 double FA_rho[7][npoints];
2325 double FA_Mc[nANNth+1][npoints];
2326 double FA_ANNMc[nANNth][npoints];
2327 double TA_Mc[nANNth+1][npoints];
2328 double TA_ANNMc[nANNth][npoints];
2329 double th_Mc[npoints];
2330 double th_k[npoints];
2331 double th_rho[npoints];
2332 
2333 for(int i=0; i<7; i++){
2334 Mccuts[i]=0.0;
2335 ANNcuts[i]=0.0;
2336 }
2337 ANNcuts[0]=-1.0;
2338 ANNcuts[1]=cutANN2;
2339 ANNcuts[2]=-1.0;
2340 Mccuts[2]=cutMc2;
2341 ANNcuts[3]=cutANN1;
2342 Mccuts[3]=cutMc1;
2343 Mccuts[4]=cutMc2;
2344 ANNcuts[4]=cutANN2;
2345 ANNcuts[5]=cutANN2;
2346 Mccuts[5]=cutMc1;
2347 Mccuts[6]=cutMc2;
2348 ANNcuts[6]=cutANN1;
2349 
2350 for(int i=0; i<7; i++){
2351  for (int j=0; j<npoints; j++){
2352  TA_rho[i][j]=0.;
2353  FA_rho[i][j]=0.;
2354  }
2355 }
2356 double delta_rho=0.;
2357 delta_rho=(double)(maxrho-minrho)/npoints;
2358 
2359 for(int i=0; i<npoints; i++){
2360 th_rho[i]=0.;
2361 th_rho[i]=minrho+i*delta_rho;
2362 }
2363 double yasym[nANNth];
2364 //Mc_g[0]=new TGraph();
2365 for (int i=0; i<nANNth; i++){
2366 // Mc_g[i+1]=new TGraph();
2367 // ANNMc_g[i]=new TGraph();
2368  for (int j=0; j<npoints; j++){
2369  FA_Mc[i+1][j]=0.;
2370  FA_ANNMc[i][j]=0.;
2371  TA_Mc[i+1][j]=0.;
2372  TA_ANNMc[i][j]=0.;
2373  }
2374 }
2375 
2376 //ANN thresholds
2377 double delta_yasym=0.;
2378 delta_yasym=(double)(max_yasym-min_yasym)/nANNth;
2379 
2380 cout<<" delta_yasym "<<delta_yasym<<endl;
2381 for (int i=0; i<nANNth; i++){
2382  yasym[i]=0.;
2383  yasym[i]=min_yasym+i*delta_yasym;
2384  cout<<" yasym[i] "<<yasym[i]<<" i "<<i<<endl;
2385 }
2386 
2387 
2388 //thresholds on the parameters necessary to built the TA-FA curves
2389 double delta_Mc=0.;
2390 //delta_Mc=(max_Mcth-min_Mcth)/nANNth;
2391 delta_Mc=(double)(max_Mcth-min_Mcth)/npoints;
2392 
2393 double delta_logk=0.;
2394 double min_logk=0.;
2395 double max_logk=0.;
2396 
2397 min_logk=log10(min_k);
2398 max_logk=log10(max_k);
2399 delta_logk=(log10(max_k)-log10(min_k))/npoints;
2400 
2401 for (int i=0; i<npoints; i++){
2402  FA_Mc[0][i]=0.;
2403  TA_Mc[0][i]=0.;
2404  th_Mc[i]=min_Mcth+i*delta_Mc;
2405  th_k[i]=pow(10,min_logk+i*delta_logk);
2406  cout<<" th_k[i] "<<th_k[i]<<" th_Mc[i] "<<th_Mc[i]<<endl;
2407  }
2408 
2409 //FA and TA definitions in function of rho
2410 int countS=0;
2411 int countB=0;
2412 for(int n=0; n<nTot; n++){
2413 NNTree2->GetEntry(n);
2414  if(cc<0.5) continue;
2415  if(n<nSig) countS=countS+1;
2416  else countB=countB+1;
2417  for(int i=0; i<7; i++){
2418  if(Mc>Mccuts[i]&&av>ANNcuts[i]) {
2419  for (int j=0; j<npoints;j++) {
2420  if(rho>th_rho[j]&&n<nSig) TA_rho[i][j]=TA_rho[i][j]+1;
2421  if(rho>th_rho[j]&&(n>nSig||n==nSig)) FA_rho[i][j]=FA_rho[i][j]+1;
2422  }
2423  }
2424  }
2425 }
2426 
2427 
2428 for(int i=0; i<7; i++){
2429  for (int j=0; j<npoints;j++) {
2430  if(i==6) cout<<"TA_rho[i][j] "<<TA_rho[i][j]<<" countS "<<countS<<" FA_rho[i][j] "<<FA_rho[i][j]<<" countB "<<countB<<endl;
2431  FA_rho[i][j]=(double)FA_rho[i][j]/countB;
2432  TA_rho[i][j]=(double)TA_rho[i][j]/countS;
2433  if(FA_rho[i][j]==0) FA_rho[i][j]=0.00001;
2434  if(i==6) cout<<"TA_rho[i][j] "<<TA_rho[i][j]<<" FA_rho[i][j] "<<FA_rho[i][j]<<endl;
2435  }
2436  }
2437 
2438 
2439 // FA and TA definitions
2440 for (int n=0; n<nTot; n++){
2441  NNTree2->GetEntry(n);
2442  if(n<nSig){
2443  for(int i=0; i<nANNth; i++){
2444  for(int j=0; j<npoints; j++){
2445  if(Mc>th_Mc[j]&&i==0)TA_Mc[0][j]=TA_Mc[0][j]+1;
2446  if(Mc>th_Mc[j]&&av>yasym[i]) TA_Mc[i+1][j]=TA_Mc[i+1][j]+1;
2447  if(av>((th_k[j]/(Mc-x_asym))+yasym[i]))TA_ANNMc[i][j]=TA_Mc[i+1][j]+1;
2448 // if(i==0&&j==0)cout<<" av "<<av<<" Mc "<<Mc<<" th_k[j] "<<th_k[j]<<" x_asym "<<x_asym<<" asym[i] "<<yasym[i]<<" ((th_k[j]/(Mc-x_asym))+yasym[i])) "<<(th_k[j]/(Mc-x_asym))+yasym[i]<<" TA_ANNMc[i][j] "<<TA_ANNMc[i][j]<<endl;
2449  }
2450  }
2451  }
2452 
2453  else {
2454  for(int i=0; i<nANNth; i++){
2455  for(int j=0; j<npoints; j++){
2456  if(Mc>th_Mc[j]&&i==0)FA_Mc[0][j]=FA_Mc[0][j]+1;
2457  if(Mc>th_Mc[j]&&av>yasym[i]) FA_Mc[i+1][j]=FA_Mc[i+1][j]+1;
2458  if(av>(th_k[j]/(Mc-x_asym)+yasym[i]))FA_ANNMc[i][j]=FA_Mc[i+1][j]+1;
2459  }
2460  }
2461  }
2462 
2463 }
2464 //exit(0);
2465 
2466 for(int j=0; j<npoints; j++){
2467  TA_Mc[0][j]=TA_Mc[0][j]/nSig;
2468  FA_Mc[0][j]=FA_Mc[0][j]/nBg;
2469  for(int i=0; i<nANNth; i++){
2470  TA_Mc[i+1][j]=(double)TA_Mc[i+1][j]/nSig;
2471  FA_Mc[i+1][j]=(double)FA_Mc[i+1][j]/nBg;
2472  TA_ANNMc[i][j]=(double)TA_ANNMc[i][j]/nSig;
2473  FA_ANNMc[i][j]=(double)FA_ANNMc[i][j]/nBg;
2474  }
2475 }
2476 
2477 Mc_g[0]=new TGraph(npoints,FA_Mc[0],TA_Mc[0]);
2478 TGraph* g0=new TGraph();
2479 g0->SetPoint(0.00001,1,0);
2480 for (int i=0; i<7; i++){
2481  rho_g[i]=new TGraph(npoints,FA_rho[i],TA_rho[i]);
2482  }
2483  TCanvas* canrho=new TCanvas("Efficiency_vs_FalseAlarm","Efficiency_vs_Alarm",0,0,1200,700);
2484  canrho->SetLogx();
2485  TLegend *legrho = new TLegend(.75, .2, .99, .4);
2486  TMultiGraph* mg0=new TMultiGraph();
2487  if(g0->GetN()!=0) mg0->Add(g0);
2488  g0->SetMarkerColor(0);
2489 
2490  mg0->SetTitle("Comparison;FalseAlarms;Efficiency");
2491  for (int i=0;i<7;i++){
2492  rho_g[i]->SetLineColor(2+i);
2493  rho_g[i]->SetMarkerColor(2+i);
2494  rho_g[i]->SetMarkerStyle(i+20);
2495  //if((i+2)!=6&&(i+2)!=7&&(i+2)!=1) ANNMc_g[i]->SetMarkerSize(0.5);
2496  rho_g[i]->SetMarkerSize(0.5);
2497  // char def[1024];
2498  char defANNMc_rho[1024];
2499  sprintf(defANNMc_rho,"Mc_%1.2f & ANNcut_%1.2f",Mccuts[i],ANNcuts[i]);
2500  if(rho_g[i]->GetN()!=0) legrho->AddEntry(rho_g[i],defANNMc_rho,"pl");
2501  if(rho_g[i]->GetN()!=0) mg0->Add(rho_g[i]);
2502  // if(ANNMc_g[0]->GetN()!=0) mg1->Add(ANNMc_g[0]);
2503  }
2504  mg0->Draw("apl");
2505  legrho->Draw();
2506  TString rhoCname(name);
2507  TString rhoCnamer(name);
2508  rhoCname.ReplaceAll(".root",".png");
2509  char rhoCname2[1024];
2510  char rhoCname2r[1024];
2511  sprintf(rhoCname2,"Cuts_caROC/TAvsFA_changingRHO%s",rhoCname.Data());
2512  sprintf(rhoCname2r,"Cuts_caROC/TAvsFA_changingRHO%s",rhoCnamer.Data());
2513  canrho->SaveAs(rhoCname2);
2514  canrho->Print(rhoCname2r);
2515 
2516 
2517 for(int j=0; j<npoints; j++) cout<<" th_Mc[j] "<<th_Mc[j]<<" j "<<j<<" TA_Mc[0][j] "<<TA_Mc[0][j]<<" FA_Mc[0][j] "<<FA_Mc[0][j]<<endl;
2518 
2519 
2520 for (int i=0; i<nANNth; i++){
2521  Mc_g[i+1]=new TGraph(npoints,FA_Mc[i+1],TA_Mc[i+1]);
2522  ANNMc_g[i]=new TGraph(npoints,FA_ANNMc[i],TA_ANNMc[i]);
2523 }
2524 
2525  TCanvas* can=new TCanvas("Efficiency_vs_FalseAlarm","Efficiency_vs_Alarm",0,0,1200,700);
2526  can->SetLogx();
2527  TLegend *leg = new TLegend(.75, .4, .99, .65);
2528  //TLegend *leg =new TLengend(.75, .80, .95, .95);
2529  TMultiGraph* mg1=new TMultiGraph();
2530  mg1->SetTitle("Comparison;FalseAlarms;Efficiency");
2531  Mc_g[0]->SetLineColor(3);
2532  Mc_g[0]->SetMarkerStyle(1);
2533  for (int i=0;i<nANNth;i++){
2534  Mc_g[i+1]->SetLineColor(4);
2535  //if((i+2)!=6&&(i+1)!=7&&(i+2)!=1)Mc_g[i+1]->SetMarkerSize(0.5);
2536  Mc_g[i+1]->SetMarkerSize(0.5);
2537  Mc_g[i+1]->SetMarkerStyle(i+20);
2538  Mc_g[i+1]->SetMarkerColor(4);
2539  ANNMc_g[i]->SetLineColor(2);
2540  ANNMc_g[i]->SetMarkerColor(2);
2541  ANNMc_g[i]->SetMarkerStyle(i+20);
2542  //if((i+2)!=6&&(i+2)!=7&&(i+2)!=1) ANNMc_g[i]->SetMarkerSize(0.5);
2543  ANNMc_g[i]->SetMarkerSize(0.5);
2544  char defMc[1024];
2545  char defANNMc[1024];
2546  sprintf(defMc,"Mc & ANNcut_%1.2f",yasym[i]);
2547  sprintf(defANNMc,"Hyperbola_cut,y_asymtoto_%1.2f","yasym[i]");
2548  if(Mc_g[i+1]->GetN()!=0) leg->AddEntry(Mc_g[i+1],defMc,"pl");
2549  if(Mc_g[i+1]->GetN()!=0) mg1->Add(Mc_g[i+1]);
2550  if(ANNMc_g[i]->GetN()!=0) leg->AddEntry(ANNMc_g[i],defANNMc,"pl");
2551  if(ANNMc_g[i]->GetN()!=0) mg1->Add(ANNMc_g[i]);
2552  // if(ANNMc_g[0]->GetN()!=0) mg1->Add(ANNMc_g[0]);
2553  }
2554  if(Mc_g[0]->GetN()!=0) mg1->Add(Mc_g[0]);
2555  if(Mc_g[0]->GetN()!=0) leg->AddEntry(Mc_g[0],"Mc, no ANNcut","pl");
2556  mg1->Draw("apl");
2557  leg->Draw();
2558 // ANNMc_g[0]= new TGraph(npoints, )
2559 // cout<<" TA_ANNMc[0][0] "<<TA_ANNMc[0][0]<<" FA_ANNMc[0][0] "<<FA_ANNMc[0][0]<<endl;
2560  TString Cname(name);
2561  TString Cnamer(name);
2562  Cname.ReplaceAll(".root",".png");
2563  char Cname2[1024];
2564  char Cname2r[1024];
2565  sprintf(Cname2,"Cuts_caROC/TAvsFA_%s",Cname.Data());
2566  sprintf(Cname2r,"Cuts_caROC/TAvsFA_%s",Cnamer.Data());
2567  can->SaveAs(Cname2);
2568  can->Print(Cname2r);
2569 
2570 }
2571 
2572 
cutANN1
#define cutANN1
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:70
rho
double rho
Definition: cwb_report_prod_2.C:149
nANNth
#define nANNth
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:59
netevent::rho
Float_t * rho
biased null statistics
Definition: netevent.hh:93
CoG_plots
void CoG_plots(TString ifile)
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:2174
cwb_online.f
tuple f
Definition: cwb_online.py:91
Ncc3D
#define Ncc3D
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:24
a
wavearray< double > a(hp.size())
name
par[0] name
Definition: CWB_Plugin_HEN_SIM_Config.C:87
n
int n
Definition: cwb_net.C:10
deltaANN
#define deltaANN
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:30
TString
TString("c")
Definition: cwb_report_skymap.C:111
RHOth_max
#define RHOth_max
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:55
leg
TLegend * leg
Definition: compare_bkg.C:246
x
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<< "\t layers : "<< nLAYERS<< "\t dF(hz) : "<< dF<< "\t 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
Definition: BaseFunctionWDM.C:51
netevent::GetEntry
Int_t GetEntry(Int_t)
Definition: netevent.cc:387
nRHO
#define nRHO
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:28
min_yasym
#define min_yasym
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:60
Get
cout<<"Comparison bkg events above threshold: "<< entriesx<< endl;double *rhox=cnet.GetV1();comph=(TH1D *) gDirectory-> Get("hcomp")
Definition: compare_bkg.C:201
minrho
#define minrho
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:73
CutMcANN
void CutMcANN(TString ifile)
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:2288
std
STL namespace.
cutANN2
#define cutANN2
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:71
maxCC
#define maxCC
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:25
NMc_r
#define NMc_r
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:43
max_Mcth
#define max_Mcth
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:62
ofile
char ofile[512]
Definition: cwb_condor_rescue.C:40
j
int j
Definition: cwb_net.C:10
i
i drho i
Definition: cwb_epparameters.C:85
deltarho
#define deltarho
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:32
TB
CWB::Toolbox TB
Definition: ComputeSNR.C:5
deltacc
#define deltacc
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:31
minMc_inj
#define minMc_inj
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:41
out
ofstream out
Definition: cwb_merge.C:196
nIFO
#define nIFO
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:27
NMc_inj
#define NMc_inj
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:40
N_ANNth
#define N_ANNth
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:51
NDIM
#define NDIM
Definition: ClusterToFrameNew01_dtmin_approx_cent.C:31
h
wavearray< double > h
Definition: Regression_H1.C:25
minMc_r
#define minMc_r
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:45
Annth
void Annth(TString ifile, int &FAth_n, int &FDth_n, double &thresFA)
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:1075
background
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 \n");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
Definition: cbc_plots_sim4_ifar.C:1070
cutMc1
#define cutMc1
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:68
PlotsAv_cc
void PlotsAv_cc(TString ifile)
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:1620
output
char output[256]
Definition: cwb_report_prod_2.C:876
k
int k
Definition: cwb_report_prod_2.C:150
netevent::netcc
Float_t * netcc
effective correlated SNR
Definition: netevent.hh:94
nANN
#define nANN
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:29
x_asym
#define x_asym
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:58
RHOth_min
#define RHOth_min
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:56
ifile
TFile * ifile
Definition: cwb_dump_config.C:23
num
int num
Definition: check_slag_bkg.C:71
yy
wavearray< double > yy
Definition: TestFrame5.C:12
minANN_av
#define minANN_av
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:47
minCC
#define minCC
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:26
cutMc2
#define cutMc2
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:69
NANN_av
#define NANN_av
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:46
s
s s
Definition: cwb_net.C:137
maxMc_r
#define maxMc_r
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:44
npoints
#define npoints
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:64
line
TLine * line
Definition: compare_bkg.C:482
graph
void graph(TString ifile)
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:748
min_Mcth
#define min_Mcth
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:63
ANNth_max
#define ANNth_max
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:53
maxrho
#define maxrho
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:72
iMc
#define iMc
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:39
max_yasym
#define max_yasym
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:61
maxMc_inj
#define maxMc_inj
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:42
Mcth
void Mcth(TString ifile, int &McFAth_n, int &McFDth_n, double &McthresFA)
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:1345
sprintf
sprintf(tfres,"(1/%g)x(%g) (sec)x(Hz)", 2 *df, df)
min_k
#define min_k
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:66
netevent::GetEntries
Int_t GetEntries()
Definition: netevent.cc:381
nCC
#define nCC
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:33
max_k
#define max_k
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:65
ANNth_min
#define ANNth_min
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:54
maxANN_av
#define maxANN_av
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:48
PlotsAv_Mc
void PlotsAv_Mc(TString ifile)
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:1747
y
wavearray< double > y
Definition: Test10.C:31
SNR_plots
void SNR_plots(TString ifile)
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:1968
netevent::chirp
Float_t * chirp
range to source: [0/1]-rec/inj
Definition: netevent.hh:109
exit
exit(0)
JoinCutANN_Mchirp_ROCcurves_rho2
void JoinCutANN_Mchirp_ROCcurves_rho2(TString NN_FILE, TString NN_FILE2, TString NN_FILE3, TString NN_FILE4, TString NN_FILE5, TString NN_FILE6, TString NN_FILE7, TString NN_FILE8, TString TEST_FILE, TString NOMEtot_S, int TS=0, int TB=0, int s=0, int b=0, int uf=1, int consider_all=0, int av_on_nNN=0)
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:94
N_RHOth
#define N_RHOth
Definition: JoinCutANN_Mchirp_ROCcurves_rho2.C:52