Logo Coherent WaveBurst  
Reference Guide 		Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Test0_dANN_var.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 nIFO 3
25 #define nRHO 3
26 #define nANN 10
27 #define deltaANN 0.1
28 #define nCC 4
29 #define NPIX 20
30 #define RHOth 5.5
31 #define nth 10
32 using namespace std;
33 void graph(TString ifile);
34 void Annth(TString ifile);
35 void Plots(TString ifile);
36 //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){
37 void Test0_dANN_var(TString NN_FILE,TString TEST_FILE, int TS=0, int TB=0, int s=0, int b=0, int uf=0){
38 
39  int p=0;
40  char NNi[1024];
41  char NNi2[1024];
42  sprintf(NNi2,"%s",NN_FILE.Data());
43  while (NNi2[p]){
44  //cout<<NNi2[p]<<endl;
45  if (NNi2[p]=='N'){
46 // cout<<p<<endl;
47  if((NNi2[p+1]=='N')&&(NNi2[p+2]=='\/')) {
48 // cout<<" dentro if"<<endl;
49  int hh=p+3;
50  while (NNi2[hh]!='\0'){NNi[hh-p-3]=NNi2[hh];hh=hh+1;}
51  break;
52  }
53  }
54  p=p+1;
55  }
56 
57  int q=0;
58  char Filei[1024];
59 // for (int i=0;i<1024;i++)Filei[i]='';
60  char Filei2[1024];
61  sprintf(Filei2,"%s",TEST_FILE.Data());
62  while (Filei2[q]){
63 // cout<<Filei2[q]<<endl;
64  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]=='\/')) {
65  int hh=q+7;
66  //while (Filei2[hh]!='\0') {Filei[hh-p-8]=Filei2[hh];hh=hh+1;cout<<Filei2[hh]<<endl;}
67  while (Filei2[hh]!='\0') {Filei[hh-p-7]=Filei2[hh];hh=hh+1;
68  //cout<<Filei2[hh]<<endl;
69  }
70  for (int h0=hh-p-7;h0<1024;h0++) Filei[h0]='\0';
71  break;
72  }
73  q=q+1;
74  }
75 
76  cout<<NNi<<" original String: "<<NNi2<<endl;
77  cout<<Filei<<" original String: "<<Filei2<<endl;
78 
79  TString NNi0(NNi);
80  TString Filei0(Filei);
81 
82  NNi0.ReplaceAll(".root","");
83  Filei0.ReplaceAll(".root","");
84 
85 
86  TFile* fnet =TFile::Open(NN_FILE.Data());
87  TMultiLayerPerceptron *mlp = (TMultiLayerPerceptron*)fnet->Get("TMultiLayerPerceptron");
88  if(mlp==NULL) {cout << "Error getting mlp" << endl;exit(1);}
89  TTree* infot=(TTree*)fnet->Get("info");
90  int NNs;
91  int NNb;
92  int NNnTS;
93  int NNnTB;
94  infot->SetBranchAddress("Rand_start_Sig",&NNs);
95  infot->SetBranchAddress("Rand_start_Bg",&NNb);
96  infot->SetBranchAddress("#trainSig",&NNnTS);
97  infot->SetBranchAddress("#trainBg",&NNnTB);
98  infot->GetEntry(0);
99 
100  TFile* fTEST =TFile::Open(TEST_FILE.Data());
101  TTree* NNTree=(TTree*)fTEST->Get("nnTree");
102  int entries=NNTree->GetEntries();
103  cout<<"entries: "<<entries<<endl;
104 
105 //estraggo le info che servono anche al tree di def dell'mlp
106  int ndim;
107  int ninp;
108  int y;
109  int entriesTot;
110  int bg_entries;
111  int sig_entries;
112 
113  NNTree->SetBranchAddress("#Entries_type",&entriesTot);
114  NNTree->SetBranchAddress("Matrix_dim",&ndim);
115  NNTree->SetBranchAddress("#inputs",&ninp);
116  NNTree->SetBranchAddress("amplitude_mode",&y);
117 
118  NNTree->GetEntry(0);
119  sig_entries=entriesTot;
120  NNTree->GetEntry(entries-1);
121  bg_entries=entriesTot;
122 
123  int const NDIM=ndim;
124  int const nINP=ninp;
125 
126  cout<<"NDIM "<<NDIM<<endl;
127  cout<<"nINP"<<nINP<<endl;
128  cout<<"sig e "<<sig_entries<<endl;
129  cout<<"bg e "<<bg_entries<<endl;
130  int minevents=0;
131  if (sig_entries>bg_entries) minevents=bg_entries;
132  else minevents=sig_entries;
133 
134  if(b==0) b=sig_entries;
135  if(TB==0 && TS==0){
136  TS=minevents;
137  TB=minevents;
138  }
139  if (b<sig_entries){
140  cout<<"Error: Bg index<sig_entries"<<endl;
141  exit(0);
142  }
143  if((TS>sig_entries||TB>bg_entries)&&(TS==TB)) {TS=minevents-s;TB=minevents-b+sig_entries;}
144  if((TS>sig_entries||TB>bg_entries)&&(TS!=TB)) {TS=sig_entries-s;TB=bg_entries-b+sig_entries;}
145 
146 
147  char NOMEtot[1024];
148  sprintf(NOMEtot,"NN_%s_FILE_%s_TS%i_TB%i_st%i_bt%i_dANN%1.2f_uf%i",NNi0.Data(),Filei0.Data(),TS,TB,s,b,deltaANN,uf);
149  cout<<"nome: "<<NOMEtot<<endl;
150 
151  TString SIG_FILE;
152  TString BG_FILE;
153  char FILE_NAME[516];
154  NNTree->SetBranchAddress("Files_name",&FILE_NAME);
155  NNTree->GetEntry(0);
156  SIG_FILE=FILE_NAME;
157  NNTree->GetEntry(entries-1);
158  BG_FILE=FILE_NAME;
159  cout<<"fine ifdef RHO_CC"<<endl;
160 
161  TChain sigTree("waveburst");//cerca il Tree "waveburst nei file
162  sigTree.Add(SIG_FILE.Data());//determina i file
163  netevent signal(&sigTree,nIFO);
164  int sig_entries2 = signal.GetEntries();
165  cout << "sig entries2 : " << sig_entries2 << endl;
166 
167  TChain bgTree("waveburst");
168  bgTree.Add(BG_FILE.Data());
169  netevent background(&bgTree,nIFO);
170  int bg_entries2 = background.GetEntries();
171  cout << "bg entries2 : " << bg_entries2 << endl;
172 
173  cout<<"b: "<<b<<endl;
174  cout<<"s: "<<s<<endl;
175 
176  // add leaf
177  Float_t x[nINP];
178  for (int jj=0; jj<nINP;jj++) x[jj]=0.;
179  char ilabel[nINP][16];
180 
181  //costruzione una foglia per ogni input
182  for(int i=0;i<nINP;i++) {
183  sprintf(ilabel[i],"x%i",i+1);
184  NNTree->SetBranchAddress(ilabel[i], &x[i]);
185  }
186 
187 char ofile[1024];
188 sprintf(ofile,"outfile/%s.root",NOMEtot);
189 //TFile*f=new TFile(ofile.Data(),"RECREATE");
190 TFile*f=new TFile(ofile,"RECREATE");
191  TTree* NNTree2=new TTree("Parameters","Parameters");
192 NNTree2->SetDirectory(f);
193  double out=0.;
194  double NNcc=0.;
195  double NNrho=0.;
196  NNTree2->Branch("ANNout",&out,"ANNout/D");
197  NNTree2->Branch("cc",&NNcc,"cc/D");
198  NNTree2->Branch("rho",&NNrho,"rho/D");
199  char NNfilename[1024];
200  NNTree2->Branch("NNname",&NNfilename,"NNname/C");
201  sprintf(NNfilename,"%s",NN_FILE.Data());
202  char Testf[1024];
203  NNTree2->Branch("TestFile",&Testf,"TestFile/C");
204  sprintf(Testf,"%s",TEST_FILE.Data());
205  int nTestS;
206  NNTree2->Branch("#TestSig",&nTestS,"#TestSig/I");
207 // nTestS=TS;
208  cout<<"nTestS: "<<nTestS<<" TS: "<<TS<<endl;
209  cout<<"dopo def tree"<<endl;
210 /* TChain TreeS("waveburst");//cerca il Tree "waveburst nei file
211  TreeS.Add(TEST_FILE.Data());//determina i file
212  netevent entryS(&TreeS,nIFO);
213  int entriesTot=0;
214  entriesTot = entryS.GetEntries();
215  cout << "Sig entries : " << entriesTot << endl;
216  std::vector<double>* frameS = new vector<double>;
217  entryS.fChain->SetBranchAddress("nnframe", &frameS);//fa esattamente come ha fatto per le altre variabili nella macro netevent.cc
218 
219  //cout<<"cc"<<(double)entryS.netcc[1]<<endl;
220 */ //cout<<"rho"<<(double)entryS.rho[0]<<endl;
221 
222 int scount=0;
223  if(uf==0) nTestS=TS;
224  else {
225 for(int n=s;n<s+TS;n++) {
226  if (uf!=0&&n>=NNs&&n<=(NNs+NNnTS)) continue;
227  scount=scount+1;
228  }
229 }
230 
231  double params[nINP];
232  int sig_05=0;
233  for (int i=0; i<nINP; i++) params[i]=0.;
234  //for(int n=0;n<entryS.GetEntries();n++) {
235  for(int n=s;n<s+TS;n++) {
236  if (uf!=0&&n>=NNs&&n<=(NNs+NNnTS)) continue;
237  NNTree->GetEntry(n);
238  signal.GetEntry(n);
239  NNcc=(double)signal.netcc[1];
240  NNrho=(double)signal.rho[0];
241  for (int i=0; i<nINP; i++){
242  //x[i]=(*frameS)[i];
243  params[i]=x[i];
244  }
245  double output=mlp->Evaluate(0,params);
246  out=output;
247  if (out<0.6) sig_05=sig_05+1;
248  cout<<"rho: "<<signal.rho[0]<<" cc "<<signal.netcc[1]<<" out: "<<out<<endl;
249  NNTree2->Fill();
250  }
251 
252 cout<<"riempito sig"<<endl;
253 int bg_05=0;
254  //for(int n=sig_entries+b;n<sig_entries+b+TB;n++) {
255  for(int n=b;n<b+TB;n++) {
256  if (uf!=0&&n>=NNb&&n<=(NNb+NNnTB)) continue;
257  NNTree->GetEntry(n);
258  cout<<"n: "<<n<<"Bg index"<<(n-sig_entries)<<endl;
259  background.GetEntry(n-sig_entries);
260  NNcc=(double)background.netcc[1];
261  NNrho=(double)background.rho[0];
262  for (int i=0; i<nINP; i++){
263  //x[i]=(*frameS)[i];
264  params[i]=x[i];
265  }
266 
267  double output=mlp->Evaluate(0,params);
268  out=output;
269  if(out>0.6) bg_05=bg_05+1;
270  cout<<"rho: "<<background.rho[0]<<" cc "<<background.netcc[1]<<" out: "<<out<<endl;
271  NNTree2->Fill();
272  }
273 cout<<"riempito bg"<<endl;
274 int Realentries=0;
275 Realentries=NNTree2->GetEntries();
276 NNTree2->Write();
277 f->Close();
278 cout<<"chiuso file"<<endl;
279 
280 //graph(ofile.Data());
281 graph(ofile);
282 cout<<"dopo richiamo funzione"<<endl;
283 Annth(ofile);
284 Plots(ofile);
285 cout<<" Sig with out<06: "<<sig_05<<endl;
286 cout<<" Bg with out>06: "<<bg_05<<endl;
287 cout<<" Realentries "<<Realentries<<endl;
288 
289 }
290 
291 void graph(TString ifile){
292  TString name(ifile);
293  name.ReplaceAll("outfile/","");
294  TFile* fTEST =TFile::Open(ifile.Data());
295  TTree* NNTree2=(TTree*)fTEST->Get("Parameters");
296  double out;
297  double cc;
298  double rho;
299  int nSi;
300  NNTree2->SetBranchAddress("ANNout",&out);
301  NNTree2->SetBranchAddress("cc",&cc);
302  NNTree2->SetBranchAddress("rho",&rho);
303  NNTree2->SetBranchAddress("#TestSig",&nSi);
304  NNTree2->GetEntry(0);
305  int const nSig=nSi;
306  cout<<"nSig: "<<nSig<<" nSi: "<<nSi<<endl;
307  int const ncurve=nANN*nCC;
308  cout<<"dentro funzione dopodef"<<endl;
309 //LOG(NBg)vs RHO------------------------------------------
310  double* rhoSig[ncurve];
311  for (int i=0;i<ncurve;i++) rhoSig[i]=new double[nSig];
312  cout<<"dopo def rhoSig"<<endl;
313  //double rhoSig[20][10];
314  int NSig[ncurve];
315 // int nBg=0;
316  int const nBg=NNTree2->GetEntries()-nSig;
317  for (int i=0;i<ncurve;i++) {
318  NSig[i]=0;
319  for (int j=0;j<nSig;j++) rhoSig[i][j]=0.;
320  }
321  cout<<"dopo def rhoSig"<<endl;
322  double* rhoBg[ncurve];
323  for (int i=0;i<ncurve;i++) rhoBg[i]=new double[nBg];
324  //double rhoBg[20][10];
325  int NBg[ncurve];
326  for (int i=0;i<ncurve;i++) {
327  NBg[i]=0;
328  for (int j=0;j<nBg;j++) rhoBg[i][j]=0.;
329  }
330  cout<<"dopo def rhoBg"<<endl;
331  double ccTh[nCC];
332  for (int i=0;i<nCC;i++) ccTh[i]=0.;
333  double NNTh[nANN];
334  for (int i=0;i<nANN;i++) NNTh[i]=0.;
335  double deltacc=0.;
336 // double deltaANN=0.;
337  deltacc=0.2/nCC;
338  //deltaANN=0.6/(nANN-1);
339 
340  cout<<NNTree2->GetEntries()<<endl;
341  for(int n=0;n<NNTree2->GetEntries();n++){
342  cout<<n<<endl;
343  NNTree2->GetEntry(n);
344  cout<<"rho "<<rho<<" cc "<<cc<<" out "<<out<<endl;
345  for(int i=0; i<nCC;i++){
346  ccTh[i]=0.5+i*deltacc;
347  if(cc<ccTh[i]) continue;
348  for(int j=0; j<nANN;j++){
349  if(j==0) NNTh[j]=-1000.;
350  //else NNTh[j]=0.5+(j-1)*deltaANN;
351  //else NNTh[j]=0.1+(j-1)*deltaANN;
352  else NNTh[j]=0.+(j)*deltaANN;
353  //else NNTh[j]=0.+(j-1)*deltaANN/1000.;
354  //else NNTh[j]=1.;
355  if(out<NNTh[j]) continue;
356  int ni=0;
357  if(n>nSig) {
358  NBg[i*nANN+j]= NBg[i*nANN+j]+1;
359  while(rhoBg[i*nANN+j][ni]!=0)ni=ni+1;
360  rhoBg[i*nANN+j][ni]=rho;
361  // cout<<"rho: "<<rho<<" colonna "<<i*nANN+j<<" riga "<<ni<<endl;
362 // cout<<" soglia_cc "<<ccTh[i]<<" soglia_ANN "<<NNTh[j]<<endl;
363  }
364  else {
365  NSig[i*nANN+j]= NSig[i*nANN+j]+1;
366  while(rhoSig[i*nANN+j][ni]!=0)ni=ni+1;
367  rhoSig[i*nANN+j][ni]=rho;
368  // cout<<"rho: "<<rho<<" colonna "<<i*nANN+j<<" riga "<<ni<<endl;
369 // cout<<" soglia_cc "<<ccTh[i]<<" soglia_ANN "<<NNTh[j]<<endl;
370  }
371  // }
372  }
373  }
374  }
375  cout<<"dopo riempimento variabili"<<endl;
376  int* indexS[ncurve];
377  for (int i=0;i<ncurve;i++) indexS[i]=new int[nSig];
378 
379  TGraph * gS[ncurve];
380  for (int y=0;y<ncurve;y++) {
381  int igS=0;
382  int igS_p=0;
383  /* int indexS[nSig];
384  double rhoS[nSig];
385  for(int i=0;i<nSig;i++) {
386  rhoS[i]=0.;
387  indexS[i]=0;
388  }*/
389  // ig=1;
390  // for(int i=0;i<nSig;i++) rhoS[i]=rhoSig[y][i];
391  gS[y]=new TGraph();
392 // gS[y]->SetMarkerStyle(7);
393  TMath::Sort(nSig,rhoSig[y],indexS[y],false);
394 // cout<<"dopo Sort "<<y<<endl;
395  for (int k=0;k<nSig;k++) {
396  int ii=indexS[y][k];
397  int yy=0;
398  if (k>0){
399 // cout<<"k "<<k<<endl;
400  int ij=indexS[y][k-1];
401  //if(rhoSig[y][ii]!=0&&rhoSig[y][ii]!=rhoSig[y][ij]) {
402  if(rhoSig[y][ii]!=0) {
403  yy=NSig[y]-igS;
404  //gS[y]->SetPoint(igS,rhoSig[y][ii],yy);
405  if(rhoSig[y][ii]!=rhoSig[y][ij]) gS[y]->SetPoint(igS_p++,rhoSig[y][ii],yy);
406  cout<<"igS"<<igS<<" x "<<rhoSig[y][ii]<<" y: "<<yy<<endl;
407  igS=igS+1;
408  }
409  }
410  else {
411  if(rhoSig[y][ii]!=0){
412  yy=NSig[y]-igS;
413  if(rhoSig[y][ii]>100)continue;
414  gS[y]->SetPoint(0,rhoSig[y][ii],yy);
415  igS=igS+1;
416 // cout<<" x "<<rhoSig[y][ii]<<" y: "<<yy<<endl;
417  }
418 
419  }
420  }
421  }
422  // cout<<"dopo inserimento puntiiS"<<endl;
423 // cout<<ncurve<<endl;
424  int* indexB[ncurve];
425  for (int i=0;i<ncurve;i++) indexB[i]=new int[nBg];
426 
427  TGraph * gB[ncurve];
428  for (int y=0;y<ncurve;y++) {
429  int igB=0;
430  int igB_p=0;
431  gB[y]=new TGraph();
432  TMath::Sort(nBg,rhoBg[y],indexB[y],false);
433 // cout<<"dopo Sort "<<y<<endl;
434  for (int k=0;k<nBg;k++) {
435  int ii=indexB[y][k];
436  int yy=0;
437  if (k>0){
438  int ij=indexB[y][k-1];
439  //if(rhoBg[y][ii]!=0&&rhoBg[y][ii]!=rhoBg[y][ij]) {
440  if(rhoBg[y][ii]!=0) {
441  yy=NBg[y]-igB;
442  //gB[y]->SetPoint(igB,rhoBg[y][ii],yy);
443  if(rhoBg[y][ii]!=rhoBg[y][ij]) gB[y]->SetPoint(igB_p++,rhoBg[y][ii],yy);
444  igB=igB+1;
445 // cout<<"igB"<<igB<<" x "<<rhoBg[y][ii]<<" y: "<<yy<<endl;
446  }
447  }
448  else {
449  if(rhoBg[y][ii]!=0) {
450  yy=NBg[y]-igB;
451  gB[y]->SetPoint(0,rhoBg[y][ii],yy);
452  igB=igB+1;
453 // cout<<"igB"<<igB<<" x "<<rhoBg[y][ii]<<" y: "<<yy<<endl;
454  }
455  }
456  }
457  }
458  cout<<"dopo inserimento puntiB"<<endl;
459  //gB[1]->SetMarkerStyle(7);
460  //gB[1]->SetPoint(1,1,2);
461  TCanvas* cS=new TCanvas("Efficiency_vs_rho","Efficiency_vs_rho",0,0,1200,700);
462  cS->Divide(2,2);
463  cS->cd(1)->SetLogy();
464  TMultiGraph* mg1=new TMultiGraph();
465 // gS[0]->SetMarkerStyle(7);
466  gS[0]->SetMarkerColor(2);
467  gS[0]->SetLineColor(2);
468  mg1->SetTitle("cc=0.5;rho;#Events");
469  if(gS[1]->GetN()!=0) mg1->Add(gS[1]);
470  if(gS[0]->GetN()!=0) mg1->Add(gS[0]);
471 // gS[0]->Draw("apl");
472  for (int h=2;h<nANN;h++){
473  gS[h]->SetMarkerColor(3);
474  gS[h]->SetLineColor(3);
475 // gS[h]->SetMarkerStyle(h+1);
476  if(gS[h]->GetN()!=0) mg1->Add(gS[h]);
477  // gS[h]->Draw("apl,same");
478  }
479  mg1->Draw("apl");
480  cS->cd(2)->SetLogy();
481  TMultiGraph* mg2=new TMultiGraph();
482 // gS[nANN]->SetMarkerStyle(7);
483  gS[nANN]->SetMarkerColor(2);
484  gS[nANN]->SetLineColor(2);
485  mg2->SetTitle("cc=0.55;rho;#Events");
486  if(gS[nANN]->GetN()!=0) mg2->Add(gS[nANN]);
487  for (int h=1;h<nANN;h++){
488  gS[nANN+h]->SetMarkerColor(3);
489  gS[nANN+h]->SetLineColor(3);
490 // gS[nANN+h]->SetMarkerStyle(h+1);
491  if(gS[nANN+h]->GetN()!=0) mg2->Add(gS[nANN+h]);
492  // gS[nANN+h]->Draw("apl,same");
493  }
494  mg2->Draw("apl");
495  cS->cd(3)->SetLogy();
496  TMultiGraph* mg3=new TMultiGraph();
497 // gS[nANN*2]->SetMarkerStyle(7);
498  gS[nANN*2]->SetMarkerColor(2);
499  gS[nANN*2]->SetLineColor(2);
500  mg3->SetTitle("cc=0.6;rho;#Events");
501  if(gS[nANN*2+1]->GetN()!=0) mg3->Add(gS[nANN*2+1]);
502  if(gS[nANN*2]->GetN()!=0) mg3->Add(gS[nANN*2]);
503  for (int h=2;h<nANN;h++){
504  gS[2*nANN+h]->SetMarkerColor(3);
505  gS[2*nANN+h]->SetLineColor(3);
506 // gS[2*nANN+h]->SetMarkerStyle(h+1);
507  if(gS[2*nANN+h]->GetN()!=0) mg3->Add(gS[2*nANN+h]);
508  // gS[2*nANN+h]->Draw("apl,same");
509  }
510  mg3->Draw("apl");
511  cS->cd(4)->SetLogy();
512  TMultiGraph* mg4=new TMultiGraph();
513 // gS[nANN*3]->SetMarkerStyle(7);
514  gS[nANN*3]->SetMarkerColor(2);
515  gS[nANN*3]->SetLineColor(2);
516  mg4->SetTitle("cc=0.65;rho;#Events");
517  if(gS[nANN*3]->GetN()!=0) mg4->Add(gS[nANN*3]);
518  // gS[nANN*3]->Draw("apl");
519  for (int h=1;h<nANN;h++){
520  gS[3*nANN+h]->SetMarkerColor(3);
521  gS[3*nANN+h]->SetLineColor(3);
522 // gS[3*nANN+h]->SetMarkerStyle(h+1);
523  if(gS[3*nANN+h]->GetN()!=0) mg4->Add(gS[3*nANN+h]);
524  // gS[3*nANN+h]->Draw("apl,same");
525  }
526  mg4->Draw("apl");
527  cout<<"nuovo canv"<<endl;
528  TCanvas* cB=new TCanvas("Number_vs_rho","Number_vs_rho",0,0,1200,700);
529  cB->Divide(2,2);
530  cB->cd(1)->SetLogy();
531  TMultiGraph* mg1B=new TMultiGraph();
532 // gB[0]->SetMarkerStyle(7);
533  gB[0]->SetMarkerColor(2);
534  gB[0]->SetLineColor(2);
535  mg1B->SetTitle("cc=0.5;rho;#Events");
536  if(gB[0]->GetN()!=0) mg1B->Add(gB[0]);
537  for (int h=1;h<nANN;h++){
538  gB[h]->SetMarkerColor(3);
539  gB[h]->SetLineColor(3);
540  // gB[h]>SetMarkerStyle(h+1);
541  if(gB[h]->GetN()!=0) mg1B->Add(gB[h]);
542  //gB[h]->Draw("apl,same");
543  }
544  mg1B->Draw("apl");
545  cB->cd(2)->SetLogy();
546  TMultiGraph* mg2B=new TMultiGraph();
547  //gB[nANN]->SetMarkerStyle(7);
548  gB[nANN]->SetMarkerColor(2);
549  gB[nANN]->SetLineColor(2);
550  mg2B->SetTitle("cc=0.55;rho;#Events");
551  if(gB[nANN]->GetN()!=0) mg2B->Add(gB[nANN]);
552  for (int h=1;h<nANN;h++){
553  gB[nANN+h]->SetMarkerColor(3);
554  gB[nANN+h]->SetLineColor(3);
555  //gB[nANN+h]>SetMarkerStyle(h+1);
556  if(gB[nANN+h]->GetN()!=0) mg2B->Add(gB[nANN+h]);
557  //gB[h]->Draw("apl,same");
558  }
559  mg2B->Draw("apl");
560  cB->cd(3)->SetLogy();
561  TMultiGraph* mg3B=new TMultiGraph();
562 // gB[2*nANN]->SetMarkerStyle(7);
563  gB[2*nANN]->SetMarkerColor(2);
564  gB[2*nANN]->SetLineColor(2);
565  mg3B->SetTitle("cc=0.6;rho;#Events");
566  if(gB[2*nANN]->GetN()!=0) mg3B->Add(gB[2*nANN]);
567  for (int h=1;h<nANN;h++){
568  gB[2*nANN+h]->SetMarkerColor(3);
569  gB[2*nANN+h]->SetLineColor(3);
570 // gB[2*nANN+h]>SetMarkerStyle(h+1);
571  if(gB[2*nANN+h]->GetN()!=0) mg3B->Add(gB[2*nANN+h]);
572  //gB[h]->Draw("apl,same");
573  }
574  mg3B->Draw("apl");
575  cB->cd(4)->SetLogy();
576  TMultiGraph* mg4B=new TMultiGraph();
577 // gB[3*nANN]->SetMarkerStyle(7);
578  gB[3*nANN]->SetMarkerColor(2);
579  gB[3*nANN]->SetLineColor(2);
580  mg4B->SetTitle("cc=0.65;rho;#Events");
581  if(gB[3*nANN]->GetN()!=0) mg4B->Add(gB[3*nANN]);
582  for (int h=1;h<nANN;h++){
583  gB[3*nANN+h]->SetMarkerColor(3);
584  gB[3*nANN+h]->SetLineColor(3);
585 // gB[3*nANN+h]>SetMarkerStyle(h+1);
586  if(gB[3*nANN+h]->GetN()!=0) mg4B->Add(gB[3*nANN+h]);
587  //gB[h]->Draw("apl,same");
588  }
589  mg4B->Draw("apl");
590 
591 // cout<<"dopo Draw()"<<endl;
592  TString CnameS(name);
593  TString CnameB(name);
594  TString CnameSroot(name);
595  TString CnameBroot(name);
596  char CnameS2[1024];
597  char CnameB2[1024];
598  char CnameS2root[1024];
599  char CnameB2root[1024];
600  CnameS.ReplaceAll(".root",".png");
601  CnameB.ReplaceAll(".root",".png");
602  sprintf(CnameS2,"logN_rho/logN_rho_S_dANN%1.2f_%s",deltaANN,CnameS.Data());
603  sprintf(CnameB2,"logN_rho/logN_rho_B_dANN%1.2f_%s",deltaANN,CnameB.Data());
604  sprintf(CnameS2root,"logN_rho/logN_rho_S_dANN%1.2f_%s",deltaANN,CnameSroot.Data());
605  sprintf(CnameB2root,"logN_rho/logN_rho_B_dANN%1.2f_%s",deltaANN,CnameBroot.Data());
606  cS->SaveAs(CnameS2);
607  cB->SaveAs(CnameB2);
608  cS->Print(CnameS2root);
609  cB->Print(CnameB2root);
610 // cout<<"fine"<<endl;
611 
612 }
613 
614 
615 void Annth(TString ifile){
616  TString name(ifile);
617  name.ReplaceAll("outfile/","");
618  TFile* fTEST =TFile::Open(ifile.Data());
619  TTree* NNTree2=(TTree*)fTEST->Get("Parameters");
620  double out;
621  double cc;
622  double rho;
623  int nSi;
624  NNTree2->SetBranchAddress("ANNout",&out);
625  NNTree2->SetBranchAddress("cc",&cc);
626  NNTree2->SetBranchAddress("rho",&rho);
627  NNTree2->SetBranchAddress("#TestSig",&nSi);
628  NNTree2->GetEntry(0);
629  int const nSig=nSi;
630 // cout<<"nSig: "<<nSig<<" nSi: "<<nSi<<endl;
631  int const ncurve2=nRHO*nCC;
632 
633 
634  double* ANNSig[ncurve2];
635  for (int i=0;i<ncurve2;i++) ANNSig[i]=new double[nSig];
636  //double rhoSig[20][10];
637  int NSig[ncurve2];
638 // int nBg=0;
639  int const nBg=NNTree2->GetEntries()-nSig;
640  for (int i=0;i<ncurve2;i++) {
641  NSig[i]=0;
642  for (int j=0;j<nSig;j++) ANNSig[i][j]=0.;
643  }
644  double* ANNBg[ncurve2];
645  for (int i=0;i<ncurve2;i++) ANNBg[i]=new double[nBg];
646 
647  //double rhoBg[20][10];
648  int NBg[ncurve2];
649  for (int i=0;i<ncurve2;i++) {
650  NBg[i]=0;
651  for (int j=0;j<nBg;j++) ANNBg[i][j]=0.;
652  }
653  double ccTh[nCC];
654  for (int i=0;i<nCC;i++) ccTh[i]=0.;
655  double rhoTh[nRHO];
656  for (int i=0;i<nRHO;i++) rhoTh[i]=0.;
657  double deltacc=0.;
658  double deltarho=0.;
659  deltacc=0.2/nCC;
660  deltarho=1./(nRHO);
661 
662 
663  for(int n=0;n<NNTree2->GetEntries();n++){
664  NNTree2->GetEntry(n);
665  cout<<"rho "<<rho<<" cc "<<cc<<" out "<<out<<endl;
666  for(int i=0; i<nCC;i++){
667  ccTh[i]=0.5+i*deltacc;
668  if(cc<ccTh[i]) continue;
669  for(int j=0; j<nRHO;j++){
670  rhoTh[j]=5+j*deltarho;
671  if(rho<rhoTh[j]) continue;
672  int ni=0;
673  if(n>nSig) {
674  NBg[i*nRHO+j]= NBg[i*nRHO+j]+1;
675  while(ANNBg[i*nRHO+j][ni]!=0)ni=ni+1;
676  ANNBg[i*nRHO+j][ni]=out;
677  // cout<<" soglia_cc "<<ccTh[i]<<" soglia_RHO "<<rhoTh[j]<<endl;
678  }
679  else {
680  NSig[i*nRHO+j]= NSig[i*nRHO+j]+1;
681  while(ANNSig[i*nRHO+j][ni]!=0)ni=ni+1;
682  ANNSig[i*nRHO+j][ni]=out;
683  // cout<<" soglia_cc "<<ccTh[i]<<" soglia_RHO "<<rhoTh[j]<<endl;
684  }
685  }
686  }
687  }
688 
689  int* indexS[ncurve2];
690  for (int i=0;i<ncurve2;i++) indexS[i]=new int[nSig];
691 
692  TGraph * gS[ncurve2];
693  for (int y=0;y<ncurve2;y++) {
694  int igS=0;
695  int igS_p=0;
696  gS[y]=new TGraph();
697  TMath::Sort(nSig,ANNSig[y],indexS[y],false);
698  for (int k=0;k<nSig;k++) {
699  int ii=indexS[y][k];
700  int yy=0;
701  if (k>0){
702  //cout<<"k "<<k<<endl;
703  int ij=indexS[y][k-1];
704  //if(ANNSig[y][ii]!=0&&ANNSig[y][ii]!=ANNSig[y][ij]) {
705  if(ANNSig[y][ii]!=0) {
706  yy=NSig[y]-igS;
707  //gS[y]->SetPoint(igS,ANNSig[y][ii],yy);
708  if(ANNSig[y][ii]!=ANNSig[y][ij]) gS[y]->SetPoint(igS_p++,ANNSig[y][ii],yy);
709  // cout<<"igS"<<igS<<" x "<<ANNSig[y][ii]<<" y: "<<yy<<endl;
710  igS=igS+1;
711 
712  }
713  }
714  else {
715  if(ANNSig[y][ii]!=0){
716  yy=NSig[y]-igS;
717  gS[y]->SetPoint(0,ANNSig[y][ii],yy);
718  igS=igS+1;
719  // cout<<" x "<<ANNSig[y][ii]<<" y: "<<yy<<endl;
720  }
721 
722  }
723  }
724  }
725 
726 
727 
728  int* indexB[ncurve2];
729  for (int i=0;i<ncurve2;i++) indexB[i]=new int[nBg];
730 
731  TGraph * gB[ncurve2];
732  for (int y=0;y<ncurve2;y++) {
733  int igB=0;
734  int igB_p=0;
735  gB[y]=new TGraph();
736  TMath::Sort(nBg,ANNBg[y],indexB[y],false);
737  // cout<<"dopo Sort "<<y<<endl;
738  for (int k=0;k<nBg;k++) {
739  int ii=indexB[y][k];
740  int yy=0;
741  if (k>0){
742  int ij=indexB[y][k-1];
743  //if(ANNBg[y][ii]!=0&&ANNBg[y][ii]!=ANNBg[y][ij]) {
744  if(ANNBg[y][ii]!=0) {
745  yy=NBg[y]-igB;
746  //gB[y]->SetPoint(igB,ANNBg[y][ii],yy);
747  if(ANNBg[y][ii]!=ANNBg[y][ij]) gB[y]->SetPoint(igB_p++,ANNBg[y][ii],yy);
748  igB=igB+1;
749  // cout<<"igB"<<igB<<" x "<<ANNBg[y][ii]<<" y: "<<yy<<endl;
750  }
751  }
752  else {
753  if(ANNBg[y][ii]!=0) {
754  yy=NBg[y]-igB;
755  gB[y]->SetPoint(0,ANNBg[y][ii],yy);
756  igB=igB+1;
757  // cout<<"igB"<<igB<<" x "<<ANNBg[y][ii]<<" y: "<<yy<<endl;
758  }
759  }
760  }
761  }
762 
763 
764  TCanvas* cS=new TCanvas("Efficiency_vs_ANN","Efficiency_vs_ANN",0,0,1200,700);
765  cS->Divide(2,2);
766  //cS->cd(1)->SetLogy();
767  cS->cd(1);
768  TMultiGraph* mg1=new TMultiGraph();
769  mg1->SetTitle("cc=0.5;ANN;#Events");
770  for (int h=0;h<nRHO;h++){
771  gS[h]->SetLineColor(4);
772  if(gS[h]->GetN()!=0) mg1->Add(gS[h]);
773  }
774  mg1->Draw("al");
775  cS->cd(2);
776  // cS->cd(2)->SetLogy();
777  TMultiGraph* mg2=new TMultiGraph();
778  mg2->SetTitle("cc=0.55;ANN;#Events");
779  for (int h=0;h<nRHO;h++){
780  gS[nRHO+h]->SetLineColor(4);
781  if(gS[nRHO+h]->GetN()!=0) mg2->Add(gS[nRHO+h]);
782  }
783  mg2->Draw("al");
784 
785  //cS->cd(3)->SetLogy();
786  cS->cd(3);
787  TMultiGraph* mg3=new TMultiGraph();
788  mg3->SetTitle("cc=0.6;ANN;#Events");
789  for (int h=0;h<nRHO;h++){
790  gS[2*nRHO+h]->SetLineColor(4);
791  if(gS[2*nRHO+h]->GetN()!=0) mg3->Add(gS[2*nRHO+h]);
792  }
793  mg3->Draw("al");
794  //cS->cd(4)->SetLogy();
795  cS->cd(4);
796  TMultiGraph* mg4=new TMultiGraph();
797  mg4->SetTitle("cc=0.65;ANN;#Events");
798  for (int h=0;h<nRHO;h++){
799  gS[3*nRHO+h]->SetLineColor(4);
800  if(gS[3*nRHO+h]->GetN()!=0) mg4->Add(gS[3*nRHO+h]);
801  }
802  mg4->Draw("al");
803 
804  TCanvas* cB=new TCanvas("Number_vs_ANN","Number_vs_ANN",0,0,1200,700);
805  cB->Divide(2,2);
806  cB->cd(1)->SetLogy();
807  TMultiGraph* mg1B=new TMultiGraph();
808  mg1B->SetTitle("cc=0.5;ANN;#Events");
809  for (int h=0;h<nRHO;h++){
810  gB[h]->SetLineColor(4);
811  if(gB[h]->GetN()!=0) mg1B->Add(gB[h]);
812  }
813  mg1B->Draw("al");
814  cB->cd(2)->SetLogy();
815  TMultiGraph* mg2B=new TMultiGraph();
816  mg2B->SetTitle("cc=0.55;ANN;#Events");
817  for (int h=0;h<nRHO;h++){
818  gB[nRHO+h]->SetLineColor(4);
819  if(gB[nRHO+h]->GetN()!=0) mg2B->Add(gB[nRHO+h]);
820  }
821  mg2B->Draw("al");
822  cB->cd(3)->SetLogy();
823  TMultiGraph* mg3B=new TMultiGraph();
824  mg3B->SetTitle("cc=0.6;ANN;#Events");
825  for (int h=0;h<nRHO;h++){
826  gB[2*nRHO+h]->SetLineColor(4);
827  if(gB[2*nRHO+h]->GetN()!=0) mg3B->Add(gB[2*nRHO+h]);
828  }
829  mg3B->Draw("al");
830  cB->cd(4)->SetLogy();
831  TMultiGraph* mg4B=new TMultiGraph();
832  mg4B->SetTitle("cc=0.6;ANN;#Events");
833  for (int h=0;h<nRHO;h++){
834  gB[3*nRHO+h]->SetLineColor(4);
835  if(gB[3*nRHO+h]->GetN()!=0) mg4B->Add(gB[3*nRHO+h]);
836  }
837  mg4B->Draw("al");
838 
839 
840  //cout<<"dopo Draw()"<<endl;
841  TString CnameS(name);
842  TString CnameB(name);
843  TString CnameSroot(name);
844  TString CnameBroot(name);
845  char CnameS2[1024];
846  char CnameB2[1024];
847  char CnameS2root[1024];
848  char CnameB2root[1024];
849  CnameS.ReplaceAll(".root",".png");
850  CnameB.ReplaceAll(".root",".png");
851  sprintf(CnameS2,"ANNthres/N_ANN_S_%s",CnameS.Data());
852  sprintf(CnameB2,"ANNthres/N_ANN_B_%s",CnameB.Data());
853  sprintf(CnameS2root,"ANNthres/N_ANN_S_%s",CnameSroot.Data());
854  sprintf(CnameB2root,"ANNthres/N_ANN_B_%s",CnameBroot.Data());
855  cS->SaveAs(CnameS2);
856  cB->SaveAs(CnameB2);
857  cS->Print(CnameS2root);
858  cB->Print(CnameB2root);
859  //cout<<"fine"<<endl;
860 
861 //CARTELLA Annth
862 
863 
864 
865 }
866 
867 
868 void Plots(TString ifile){
869  TString name(ifile);
870  name.ReplaceAll("outfile/","");
871  TFile* fTEST =TFile::Open(ifile.Data());
872  TTree* NNTree2=(TTree*)fTEST->Get("Parameters");
873  double out;
874  double cc;
875  double rho;
876  int nSi;
877  NNTree2->SetBranchAddress("ANNout",&out);
878  NNTree2->SetBranchAddress("cc",&cc);
879  NNTree2->SetBranchAddress("rho",&rho);
880  NNTree2->SetBranchAddress("#TestSig",&nSi);
881  NNTree2->GetEntry(0);
882  int const nSig=nSi;
883  cout<<"nSig: "<<nSig<<" nSi: "<<nSi<<endl;
884  int const nBg=NNTree2->GetEntries()-nSig;
885 
886 /* TH2D* hglitch=new TH2D("Purezza","Purezza",NCC,minCC,maxCC,NRHO,minRHO,maxRHO);
887  TH2D* h2BgO_R=new TH2D("rho_cc_out","rho_cc_out",50,0.,50,50,-0.5,1.5);
888  TH2D* h2SigO_R=new TH2D("rho_cc_out","rho_cc_out",50,0.,50,50,-0.5,1.5);
889  TH2D* h2BgR_C=new TH2D("rho_cc_out","rho_cc_out",50,0.,1.2,50,0.,50);
890  TH2D* h2SigR_C=new TH2D("rho_cc_out","rho_cc_out",50,0.,1.2,50,0.,50);
891  TH2D* h2Bg=new TH2D("rho_cc_out","rho_cc_out",50,0.,1.2,50,-0.5,1.5);
892  TH2D* h2Sig=new TH2D("rho_cc_out","rho_cc_out",50,0.,1.2,50,-0.5,1.5);
893  h2BgO_R->SetDirectory(0);
894  h2SigO_R->SetDirectory(0);
895  h2BgR_C->SetDirectory(0);
896  h2SigR_C->SetDirectory(0);
897  h2Bg->SetDirectory(0);
898  h2Sig->SetDirectory(0);
899  hglitch->SetDirectory(0);
900 */
901 
902  TGraph * gS[3];
903  TGraph * gB[3];
904  gB[0]=new TGraph();
905  gS[0]=new TGraph();
906  gB[1]=new TGraph();
907  gS[1]=new TGraph();
908  gB[2]=new TGraph();
909  gS[2]=new TGraph();
910  double aANNS[nSig];
911  double aRHOS[nSig];
912  double aCCS[nSig];
913  for (int i=0;i<nSig;i++){
914  aCCS[i]=0.;
915  aRHOS[i]=0.;
916  aANNS[i]=0.;
917  }
918  double aANNB[nBg];
919  double aRHOB[nBg];
920  double aCCB[nBg];
921  for (int i=0;i<nBg;i++){
922  aCCB[i]=0.;
923  aRHOB[i]=0.;
924  aANNB[i]=0.;
925  }
926 // cout<<"nSig: "<<nSig<<endl;
927  for (int n = 0; n <NNTree2->GetEntries(); n++){
928  NNTree2->GetEntry(n);
929  if(n<nSig) {
930  aRHOS[n]=rho;
931  aANNS[n]=out;
932  aCCS[n]=cc;
933 /* gS[0]->SetPoint(n,cc,rho);
934  gS[1]->SetPoint(n,cc,out);
935  gS[2]->SetPoint(n,rho,out);
936  cout<<"Sig_graph1: x="<<cc<<" y: "<<rho<<endl;
937  cout<<"Sig_graph2: x="<<cc<<" y: "<<out<<endl;
938  cout<<"Sig_graph3: x="<<rho<<" y: "<<out<<endl;
939 */
940  }
941  else {
942  aRHOB[n-nSig]=rho;
943  aANNB[n-nSig]=out;
944  aCCB[n-nSig]=cc;
945 /* gB[0]->SetPoint(n-nSig,cc,rho);
946  gB[1]->SetPoint(n-nSig,cc,out);
947  gB[2]->SetPoint(n-nSig,rho,out);
948  cout<<"Bg_graph1: x="<<cc<<" y: "<<rho<<endl;
949  cout<<"Bg_graph2: x="<<cc<<" y: "<<out<<endl;
950  cout<<"Bg_graph3: x="<<rho<<" y: "<<out<<endl;
951 */
952  }
953  }
954 
955  int indexB[nBg];
956  for (int k=0;k<nBg;k++)indexB[k]=0;
957  TMath::Sort(nBg,aRHOB,indexB,false);
958  int Z[(nth+1)*nth];
959  for (int k=0;k<(nth+1)*nth;k++)Z[k]=0;
960  double cc1th[nth];
961  double ANN1th[nth+1];
962  ANN1th[0]=0.;
963  for (int k=0;k<(nth);k++){
964  ANN1th[k+1]=0.;
965  cc1th[k]=0.;
966  }
967  double dANNth=0.;
968  double dccth=0.;
969  dANNth=1./nth;
970  dccth=0.5/nth;
971  ANN1th[0]=-100.;
972  for (int k=0;k<(nth);k++){
973  ANN1th[k+1]=0.+dANNth*(k+1);
974  cc1th[k]=0.5+dccth*k;
975  }
976 
977  /*for (int a=0;a<nBg;a++){
978  if(aRHOB[a]>=RHOth){
979  for (int b=0;b<nth;b++){
980  if(aCCB[a]>=cc1th[b]){
981  for(int c=0;c<nth+1;c++){
982  if(aANNB[a]>=ANN1th[c]) Z[b*nth+c]=Z[b*nth+c]+1;
983  if(c=10) cout<<" ANN "<<ANN1th[c]<<" Zcount: "<<Z[b*nth+c]<<" aNNB "<<aANNB[a]<<" a "<<a<<" b "<<b<<endl;
984  }
985  }
986  }
987  }
988  }*/
989  for (int a=0;a<nBg;a++){
990  cout<<"tutti: a "<<a<<" rho "<<aRHOB[a]<<" aNNB "<<aANNB[a]<<" cc "<<aCCB[a]<<endl;
991  if(aRHOB[a]<RHOth)continue;
992  cout<<"dopo rho: a "<<a<<" rho "<<aRHOB[a]<<" aNNB "<<aANNB[a]<<" cc "<<aCCB[a]<<endl;
993  for (int b=0;b<nth;b++){
994  if(aCCB[a]<cc1th[b])continue;
995  cout<<"dopo cc: "<<cc1th[b]<<"a "<<a<<" rho "<<aRHOB[a]<<" aNNB "<<aANNB[a]<<" cc "<<aCCB[a]<<endl;
996  for(int c=0;c<nth+1;c++){
997  if(aANNB[a]<ANN1th[c]) continue;
998  cout<<"dopo ANN:"<<ANN1th[c]<<" a "<<a<<" rho "<<aRHOB[a]<<" aNNB "<<aANNB[a]<<" cc "<<aCCB[a]<<endl;
999  Z[b*(nth+1)+c]=Z[b*(nth+1)+c]+1;
1000  cout<<"Zindex: "<<b*(nth+1)+c<<" Z: "<<Z[b*(nth+1)+c]<<endl;
1001 
1002  //if(c=10) cout<<" ANN "<<ANN1th[c]<<" Zcount: "<<Z[b*(nth+1)+c]<<" aNNB "<<aANNB[a]<<" a "<<a<<" b "<<b<<endl;
1003  }
1004 
1005  }
1006 
1007  }
1008 
1009  //TH2D* hglitch=new TH2D("Survived_glitches(RHO_cut:5.5)","Survived_glitches(RHO_cut:5.5)",nth,0.5,1,nth+1,-1,1);
1010  double bins[nth+1];
1011  for (int d=0; d<nth+1;d++) bins[d]=0.;
1012  //bins[0]=0.;
1013  for (int d=0; d<nth+1;d++) bins[d+1]=0.+(d+1)*dANNth;
1014  TH2D* hglitch=new TH2D("Survived_glitches(RHO_cut:5.5)","Survived_glitches(RHO_cut:5.5)",nth,0.5,1,nth+1,bins);
1015 // TH2D* hglitch=new TH2D("Survived_glitches(RHO_cut:5.5)","Survived_glitches(RHO_cut:5.5)",nth,0.5,1,nth+1,0.,1.);
1016  hglitch->SetDirectory(0);
1017  //if(tot[NRHO*ii+ji]>0) hglitch->Fill(cc1[ji]+deltacc/2,rho1[ii]+deltarho/2,pur[NRHO*ii+ji]);
1018  for (int b=0;b<nth;b++) for(int c=0;c<nth+1;c++) {
1019  //hglitch->Fill(cc1th[b]+dccth/2.,ANN1th[c]+dANNth/2.,Z[(nth+1)*b+c]);
1020  hglitch->Fill(cc1th[b]+dccth/2.,bins[c]+dANNth/2.,Z[(nth+1)*b+c]);
1021  cout<<" Zindex: "<<(nth+1)*b+c<< " x: "<<cc1th[b]<<" y: "<<ANN1th[c]<<" z: "<<Z[(nth+1)*b+c]<<endl;
1022  }
1023 
1024  gB[0]=new TGraph(nBg,aCCB,aRHOB);
1025  gS[0]=new TGraph(nSig,aCCS,aRHOS);
1026  gB[1]=new TGraph(nBg,aCCB,aANNB);
1027  gS[1]=new TGraph(nSig,aCCS,aANNS);
1028  gB[2]=new TGraph(nBg,aRHOB,aANNB);
1029  gS[2]=new TGraph(nSig,aRHOS,aANNS);
1030 
1031  for(int i=0;i<3;i++){
1032  gS[i]->SetMarkerColor(2);
1033  gB[i]->SetMarkerColor(4);
1034  gS[i]->SetMarkerStyle(6);
1035  gB[i]->SetMarkerStyle(7);
1036  }
1037 
1038  TCanvas* c=new TCanvas("Plots","Plots",0,0,1200,700);
1039  c->Divide(2,2);
1040  c->cd(1);
1041  TMultiGraph* mg1=new TMultiGraph();
1042  mg1->SetTitle("cc_rho;cc;rho");
1043  if(gB[0]->GetN()!=0) mg1->Add(gB[0]);
1044  if(gS[0]->GetN()!=0) mg1->Add(gS[0]);
1045  mg1->Draw("ap");
1046  c->cd(2);
1047  TMultiGraph* mg2=new TMultiGraph();
1048  mg2->SetTitle("cc_ANNoutput;cc;ANNoutput");
1049  if(gB[1]->GetN()!=0) mg2->Add(gB[1]);
1050  if(gS[1]->GetN()!=0) mg2->Add(gS[1]);
1051  mg2->Draw("ap");
1052  c->cd(3);
1053  TMultiGraph* mg3=new TMultiGraph();
1054  mg3->SetTitle("rho_ANNoutput;rho;ANNoutput");
1055  if(gB[2]->GetN()!=0) mg3->Add(gB[2]);
1056  if(gS[2]->GetN()!=0) mg3->Add(gS[2]);
1057  mg3->Draw("ap");
1058  c->cd(4);
1059  TText* text=new TText(0.37,0.0,"no cuts on ANN");
1060  hglitch->SetStats(0);
1061  hglitch->GetXaxis()->SetTitle("cc");
1062  hglitch->GetYaxis()->SetTitle("ANNoutput");
1063  hglitch->GetZaxis()->SetTitle("count");
1064  hglitch->Draw("colz");
1065  text->Draw();
1066  gPad->SetLogz();
1067 
1068  TString Cname(name);
1069  TString Cnameroot(name);
1070  char Cname2[1024];
1071  char Cname2root[1024];
1072  Cname.ReplaceAll(".root",".png");
1073  sprintf(Cname2,"CC_RHO_ANN_Plots/CC_RHO_ANN_Plots_%s",Cname.Data());
1074  sprintf(Cname2root,"CC_RHO_ANN_Plots/CC_RHO_ANN_Plots_%s",Cnameroot.Data());
1075  c->SaveAs(Cname2);
1076  c->Print(Cname2root);
1077 
1078  TCanvas* c2=new TCanvas("Plots_Bkg_on_Sig","Plots_Bkg_on_Sig",0,0,1200,700);
1079  c2->Divide(2,2);
1080  c2->cd(1);
1081  gS[0]->Draw("ap");
1082  gB[0]->Draw("p,same");
1083 
1084  c2->cd(2);
1085  gS[1]->Draw("ap");
1086  gB[1]->Draw("p,same");
1087  c2->cd(3);
1088  gS[2]->Draw("ap");
1089  gB[2]->Draw("p,same");
1090  c2->cd(4);
1091  TText* text2=new TText(0.37,0.0,"no cuts on ANN");
1092  hglitch->SetStats(0);
1093  hglitch->GetXaxis()->SetTitle("cc");
1094  hglitch->GetYaxis()->SetTitle("ANNoutput");
1095  hglitch->GetZaxis()->SetTitle("count");
1096  hglitch->Draw("colz");
1097  text2->Draw();
1098  gPad->SetLogz();
1099 
1100  TString Cname_2(name);
1101  TString Cname_2root(name);
1102  char Cname2_2[1024];
1103  char Cname2_2root[1024];
1104  Cname_2.ReplaceAll(".root",".png");
1105  sprintf(Cname2_2,"CC_RHO_ANN_Plots/CC_RHO_ANN_Plots_BoS_%s",Cname_2.Data());
1106  sprintf(Cname2_2root,"CC_RHO_ANN_Plots/CC_RHO_ANN_Plots_BoS_%s",Cname_2root.Data());
1107  c2->SaveAs(Cname2_2);
1108  c2->Print(Cname2_2root);
1109 
1110 //CARTELLA CCi_RHO_ANN_Plots
1111 }
1112 
text
TText * text
Definition: cbc_plots.C:717
rho
double rho
Definition: cwb_report_prod_2.C:149
graph
void graph(TString ifile)
Definition: Test0_dANN_var.C:291
netevent::rho
Float_t * rho
biased null statistics
Definition: netevent.hh:93
cwb_online.f
tuple f
Definition: cwb_online.py:91
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
TString
TString("c")
Definition: cwb_report_skymap.C:111
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
std
STL namespace.
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
Plots
void Plots(TString ifile)
Definition: Test0_dANN_var.C:868
TB
CWB::Toolbox TB
Definition: ComputeSNR.C:5
out
ofstream out
Definition: cwb_merge.C:196
nth
#define nth
Definition: Test0_dANN_var.C:31
NDIM
#define NDIM
Definition: ClusterToFrameNew01_dtmin_approx_cent.C:31
h
wavearray< double > h
Definition: Regression_H1.C:25
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
nIFO
#define nIFO
Definition: Test0_dANN_var.C:24
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
ifile
TFile * ifile
Definition: cwb_dump_config.C:23
yy
wavearray< double > yy
Definition: TestFrame5.C:12
s
s s
Definition: cwb_net.C:137
nRHO
#define nRHO
Definition: Test0_dANN_var.C:25
Annth
void Annth(TString ifile)
Definition: Test0_dANN_var.C:615
deltarho
#define deltarho
Definition: Average_FAout_moreGraphs_3.C:32
nCC
#define nCC
Definition: Test0_dANN_var.C:28
sprintf
sprintf(tfres,"(1/%g)x(%g) (sec)x(Hz)", 2 *df, df)
c2
TCanvas * c2
Definition: slag.C:207
netevent::GetEntries
Int_t GetEntries()
Definition: netevent.cc:381
RHOth
#define RHOth
Definition: Test0_dANN_var.C:30
Test0_dANN_var
void Test0_dANN_var(TString NN_FILE, TString TEST_FILE, int TS=0, int TB=0, int s=0, int b=0, int uf=0)
Definition: Test0_dANN_var.C:37
y
wavearray< double > y
Definition: Test10.C:31
nANN
#define nANN
Definition: Test0_dANN_var.C:26
deltacc
#define deltacc
Definition: Average_FAout_moreGraphs_3.C:31
deltaANN
#define deltaANN
Definition: Test0_dANN_var.C:27
exit
exit(0)