Logo coherent WaveBurst  
Library Reference Guide
Logo
cwb_mkcfad.C
Go to the documentation of this file.
1 /*
2 # Copyright (C) 2019 Gabriele Vedovato
3 #
4 # This program is free software: you can redistribute it and/or modify
5 # it under the terms of the GNU General Public License as published by
6 # the Free Software Foundation, either version 3 of the License, or
7 # (at your option) any later version.
8 #
9 # This program is distributed in the hope that it will be useful,
10 # but WITHOUT ANY WARRANTY; without even the implied warranty of
11 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 # GNU General Public License for more details.
13 #
14 # You should have received a copy of the GNU General Public License
15 # along with this program. If not, see <https://www.gnu.org/licenses/>.
16 */
17 
18 
19 // This macro is used to computed the false alarm probability of combined searches
20 // Uses as input parameters the files produced by the report sim when pp_fad is enabled
21 //
22 // fapList :
23 // the fapList is the name of the file which contains the list of searches
24 // for each search the format is :
25 //
26 // $SEARCH TAG COLOR STYPE
27 // - TAG : label to be used in the legend
28 // - COLOR : line color
29 // - STYPE : line style
30 //
31 // $FAD FADvsRHO PRODvsFAD
32 // - FADvsRHO : file name of FAD vs RHO
33 // - PRODvsFAD : file name of PROD vs FAD
34 //
35 // $FAR FARvsRHO OBSTYPEvsFAR
36 // - FARvsRHO : file name of FAR vs RHO
37 // - OBSTIMEvsFAR : file name of OBSTIME vs FAR
38 //
39 // Example :
40 //
41 // $SEARCH LH_S6B_LF_LF 1 1
42 // $FAR FARvsRHO_ALL.txt OBSTIMEvsFAR_ALL.txt
43 // $FAD FADvsRHO_ALL.txt PRODvsFAD_ALL.txt
44 //
45 // fapMode :
46 // the fapMode is the mode used to produce the combined searches
47 // the fapMode can contains a combination of the following options : FAD/FAR/RHO
48 // the combined curves are obtained from the minimum FAP_FAD, FAP_FAR, FAP_RHO
49 // where FAP_XXX is the False Alarm Probability computed using the ranked
50 // events according to XXX=(FAD/FAR/RHO)
51 //
52 // Examples :
53 //
54 // fadMode = "FAD"
55 // fadMode = "FAD FAR"
56 // fadMode = "FAD FAR RHO"
57 //
58 // Auxiliary parameters :
59 // rho_min : minimum rho showed in the plot (def=0, the value is computed from input files)
60 // rho_max : maximum rho showed in the plot (def=0, the value is computed from input files)
61 // refline : vertical line showed as reference in the plot (def=0, not showed)
62 // reflabel: refernce label used in the legend (def="reference")
63 // pfile : saved file name. If !="" -> save plot with name pfile (def="")
64 // if == "batch" the macro is executed in batch mode
65 //
66 // How to run (example) :
67 // root 'cwb_mkcfad.C("FAP.lst","FAD FAR RHO",5,40,7.35,"BigDog")'
68 //
69 
70 #define MAX_FAP 20
71 #define SHOW_3SIGMA // plot 3 sigma probability line
72 #define PROB_3SIGMA (1./370.398) // percentage outside 3 sigma gaussian probability
73 
74 // NOTE : FAX = FAD or FAR
75 
76 //#define XAXIS_IS_IFAX // uncomment to use as xaxis IFAX instead of RHO
77 
78 double GetFAP(double rho, int n, int nFAP, TString* FAXvsRHO, TString* PRODvsFAX);
79 
80 void GetFAPvsRHO(int n, vector<double>& x, vector<double>& y,
81  int nFAP, TString* FAXvsRHO, TString* PRODvsFAX);
82 
83 void GetFAPvsRHO(vector<double>& x, vector<double>& y,
84  int nFAP, TString* FARvsRHO, TString* OBSTIMEvsFAR);
85 
86 TGraph* GetGraph(vector<double> x, vector<double> y,
87  TString ptitle, TString xtitle, TString ytitle);
88 
89 double gRHO_MIN = 0;
90 double gRHO_MAX = 0;
91 
92 double gFAP_MIN = 1.79769e+308;
93 double gFAP_MAX = -1.79769e+308;
94 
95 void cwb_mkcfad(TString fapList, TString fapMode, double rho_min=0, double rho_max=0,
96  double refline=0, TString reflabel="reference", TString pfile="") {
97 
98  bool batch=false;
99  if(pfile=="batch") {gROOT->SetBatch(true);pfile="";batch=true;}
100 
101  if(pfile!="" && !pfile.EndsWith(".png")) {
102  cout << endl;
103  cout << "cwb_mkcfad : Error in pfile name : " << pfile << endl;
104  cout << " file name must ends with .png" << endl << endl;
105  exit(1);
106  }
107 
108  gRHO_MIN = rho_min;
109  gRHO_MAX = rho_max;
110 
111  TString FARvsRHO[MAX_FAP];
112  TString OBSTIMEvsFAR[MAX_FAP];
113  TString FADvsRHO[MAX_FAP];
114  TString PRODvsFAD[MAX_FAP];
115  TString FAP_LABEL[MAX_FAP];
116  int FAP_COLOR[MAX_FAP];
117  int FAP_STYLE[MAX_FAP];
118  unsigned int FAP_TYPE[MAX_FAP];
119  for(int i=0;i<MAX_FAP;i++) FAP_TYPE[i]=0;
120 
121  // Open list
122  ifstream in;
123  in.open(fapList.Data(),ios::in);
124  if (!in.good()) {cout << "Error Opening File : " << fapList.Data() << endl;exit(1);}
125 
126  int size=0;
127  char str[1024];
128  int fpos=0;
129  while(true) {
130  in.getline(str,1024);
131  if (!in.good()) break;
132  if(str[0] != '#') size++;
133  }
134  cout << "size " << size << endl;
135  in.clear(ios::goodbit);
136  in.seekg(0, ios::beg);
137  if (size==0) {cout << "Error : File " << fapList.Data() << " is empty" << endl;exit(1);}
138 
139  char iFARvsRHO[1024];
140  char iOBSTIMEvsFAR[1024];
141  char iFADvsRHO[1024];
142  char iPRODvsFAD[1024];
143  char iFAP_LABEL[1024];
144  char iFAP_TAG[1024];
145  int iFAP_COLOR;
146  int iFAP_STYLE;
147  int nFAP=0;
148  char line[1024];
149  while(true) {
150  in.getline(line,1024);
151  if (in.eof()) break;
152  std::stringstream linestream(line);
153  if(TString(line).BeginsWith("#")) continue;
154  if(TString(line)=="") continue;
155  if(TString(line).BeginsWith("$SEARCH")) {
156  if(FAP_TYPE[nFAP]>1) nFAP++;
157  linestream >> iFAP_TAG >> iFAP_LABEL >> iFAP_COLOR >> iFAP_STYLE;
158  FAP_LABEL[nFAP]=iFAP_LABEL;
159  FAP_COLOR[nFAP]=iFAP_COLOR;
160  FAP_STYLE[nFAP]=iFAP_STYLE;
161  FAP_TYPE[nFAP]=1;
162  }
163  if(TString(line).BeginsWith("$FAR")) {
164  linestream >> iFAP_TAG >> iFARvsRHO >> iOBSTIMEvsFAR;
165  FARvsRHO[nFAP]=iFARvsRHO;
166  OBSTIMEvsFAR[nFAP]=iOBSTIMEvsFAR;
167  FAP_TYPE[nFAP]|=2;
168  }
169  if(TString(line).BeginsWith("$FAD")) {
170  linestream >> iFAP_TAG >> iFADvsRHO >> iPRODvsFAD;
171  FADvsRHO[nFAP]=iFADvsRHO;
172  PRODvsFAD[nFAP]=iPRODvsFAD;
173  FAP_TYPE[nFAP]|=4;
174  }
175  }
176  if(FAP_TYPE[nFAP]>1) nFAP++;
177  in.close();
178 
179  for(int i=0;i<nFAP;i++) {
180  cout << "SEARCH : " << FAP_LABEL[i] << " " << FAP_COLOR[i] << " " << FAP_STYLE[i] << endl;
181  if(FAP_TYPE[i]&2) cout << FARvsRHO[i] << " " << OBSTIMEvsFAR[i] << endl;
182  if(FAP_TYPE[i]&4) cout << FADvsRHO[i] << " " << PRODvsFAD[i] << endl;
183  cout << endl;
184  for(int j=i;j<nFAP;j++) if(FAP_TYPE[i]!=FAP_TYPE[j]) {
185  cout << "cwb_mkcfap : Error - searches not consistent !!!" << endl;
186  cout << " Missing FAR or FAD declaration" << endl << endl;
187  exit(1);
188  }
189  if((!fapMode.Contains("RHO"))&&(!fapMode.Contains("FAR"))&&(!fapMode.Contains("FAD"))) {
190  cout << "cwb_mkcfap : Error - fapMode not contains valid declarations (FAD/FAR/RHO)" << endl;
191  exit(1);
192  }
193  if((!FAP_TYPE[i]&2)&&(fapMode.Contains("RHO"))) {
194  cout << "cwb_mkcfap : Error - fapMode=RHO needs FAR declarations" << endl;
195  exit(1);
196  }
197  if((!FAP_TYPE[i]&2)&&(fapMode.Contains("FAR"))) {
198  cout << "cwb_mkcfap : Error - fapMode=FAR needs FAR declarations" << endl;
199  exit(1);
200  }
201  if((!FAP_TYPE[i]&2)&&(fapMode.Contains("FAD"))) {
202  cout << "cwb_mkcfap : Error - fapMode=FAD needs FAD declarations" << endl;
203  exit(1);
204  }
205  }
206 
207  if((fapMode.Contains("RHO"))&&(!fapMode.Contains("FAR"))&&(!fapMode.Contains("FAD"))) {
208  for(int i=0;i<nFAP;i++) FAP_COLOR[i]=2; // red
209  }
210 
211  // get FAP for RHO=refline
212  if(refline>0) {
213  cout << endl << "---------------------------------------------------------------------" << endl;
214  for(int i=0;i<nFAP;i++) {
215  double FAP;
216  if(FAP_TYPE[i]&2) FAP = GetFAP(refline, i, nFAP, FARvsRHO, OBSTIMEvsFAR);
217  if(FAP_TYPE[i]&4) FAP = GetFAP(refline, i, nFAP, FADvsRHO, PRODvsFAD);
218  cout << "Search : " << FAP_LABEL[i] << "\t -> FAP at RHO = " << refline << " is " << FAP << endl;
219  }
220  cout << "---------------------------------------------------------------------" << endl << endl;
221  if(batch) exit(0);
222  }
223 
224  // create canvas
225  TCanvas* gCANVAS;
226  gCANVAS = new TCanvas("canvas", "canvas",16,30,825,546);
227  gCANVAS->Range(-19.4801,-9.25,-17.4775,83.25);
228  gCANVAS->SetBorderSize(2);
229  gCANVAS->SetFrameFillColor(0);
230  gCANVAS->SetGridx();
231  gCANVAS->SetGridy();
232 //#ifdef XAXIS_IS_IFAD
233  gCANVAS->SetLogx();
234 //#endif
235  gCANVAS->SetLogy();
236  gStyle->SetOptFit(kTRUE);
237 
238  // create graphs
239  double rho_max=0;
240  TGraph* gr[MAX_FAP];
241  vector<double> xFAP[MAX_FAP],yFAP[MAX_FAP];
242  vector<double> xFAD[MAX_FAP],yFAD[MAX_FAP];
243  vector<double> xFAR[MAX_FAP],yFAR[MAX_FAP];
244  vector<double> xRHO[MAX_FAP],yRHO[MAX_FAP];
245 
246  for(int i=0;i<nFAP;i++) {
247  if(fapMode.Contains("FAD")) {
248  // get FAP from events ranked by FAD
249  GetFAPvsRHO(i, xFAD[i], yFAD[i], nFAP, FADvsRHO, PRODvsFAD);
250  int xsize = xFAP[i].size();
251  if(xsize) {for(int j=0;j<xsize;j++) if(yFAD[i][j]<yFAP[i][j]) yFAP[i][j]=yFAD[i][j];}
252  else {xFAP[i]=xFAD[i];yFAP[i]=yFAD[i];}
253  }
254  if(fapMode.Contains("FAR")) {
255  // get FAP from events ranked by FAR
256  GetFAPvsRHO(i, xFAR[i], yFAR[i], nFAP, FARvsRHO, OBSTIMEvsFAR);
257  // select minimum between events ranked by FAR and by FAD
258  int xsize = xFAP[i].size();
259  if(xsize) {for(int j=0;j<xsize;j++) if(yFAR[i][j]<yFAP[i][j]) yFAP[i][j]=yFAR[i][j];}
260  else {xFAP[i]=xFAR[i];yFAP[i]=yFAR[i];}
261  }
262  if(fapMode.Contains("RHO")) {
263  // get FAP from events ranked by RHO
264  GetFAPvsRHO(xRHO[i], yRHO[i], nFAP, FARvsRHO, OBSTIMEvsFAR);
265  // select minimum between events ranked by RHO and by FAR/FAD
266  int xsize = xFAP[i].size();
267  if(xsize) {for(int j=0;j<xsize;j++) if(yRHO[i][j]<yFAP[i][j]) yFAP[i][j]=yRHO[i][j];}
268  else {xFAP[i]=xRHO[i];yFAP[i]=yRHO[i];}
269  }
270  }
271 
272  for(int i=0;i<nFAP;i++) {
273  int nMODE=0;
274  TString subtitle="";
275  if(fapMode.Contains("FAD")) {subtitle+=" FAD ";nMODE++;}
276  if(fapMode.Contains("FAR")) {subtitle+=" FAR ";nMODE++;}
277  if(fapMode.Contains("RHO")) {subtitle+=" RHO ";nMODE++;}
278  TString title = "";
279  if(nMODE==1) {
280  title = TString::Format("False Alarm Probability - ranked by %s",subtitle.Data());
281  } else {
282  title = TString::Format("False Alarm Probability - minimum FAP ranked by (%s)",subtitle.Data());
283  }
284 #ifdef XAXIS_IS_IFAD
285  gr[i] = GetGraph(xFAP[i],yFAP[i],title,"IFAD","FAP");
286 #else
287  gr[i] = GetGraph(xFAP[i],yFAP[i],title,"rho","FAP");
288 #endif
289  gr[i]->SetLineColor(FAP_COLOR[i]);
290  gr[i]->SetLineStyle(FAP_STYLE[i]);
291  if(xFAP[i][xFAP[i].size()-1] > rho_max) rho_max=xFAP[i][xFAP[i].size()-1];
292  }
293  if((gRHO_MAX>0)&&(rho_max>gRHO_MAX)) rho_max = gRHO_MAX;
294 
295  TMultiGraph* mg = new TMultiGraph();
296  mg->SetTitle(title);
297  for(int i=0;i<nFAP;i++) mg->Add(gr[i]);
298 
299  // add reference line
300  if(refline>0) {
301  gr[nFAP] = new TGraph;
302  gr[nFAP]->SetPoint(0,refline,0);
303  gr[nFAP]->SetPoint(1,refline,1);
304  gr[nFAP]->SetLineColor(kGreen);
305  gr[nFAP]->SetLineWidth(2);
306  mg->Add(gr[nFAP]);
307  }
308 
309  // add 3sigma reference line
310 #ifdef SHOW_3SIGMA
311  gr[nFAP+1] = new TGraph;
312  gr[nFAP+1]->SetPoint(0,gRHO_MIN,PROB_3SIGMA);
313  gr[nFAP+1]->SetPoint(1,rho_max,PROB_3SIGMA);
314  gr[nFAP+1]->SetLineColor(kBlue);
315  gr[nFAP+1]->SetLineWidth(2);
316  gr[nFAP+1]->SetLineStyle(10);
317  mg->Add(gr[nFAP+1]);
318 #endif
319 
320  mg->Paint("APL");
321 
322  mg->GetHistogram()->GetXaxis()->SetLabelSize(0.04);
323  mg->GetHistogram()->GetYaxis()->SetLabelSize(0.04);
324  mg->GetHistogram()->GetXaxis()->SetTitleSize(0.04);
325  mg->GetHistogram()->GetYaxis()->SetTitleSize(0.04);
326  mg->GetHistogram()->GetXaxis()->SetLabelFont(42);
327  mg->GetHistogram()->GetXaxis()->SetRangeUser(gRHO_MIN,gRHO_MAX);
328  mg->GetHistogram()->GetYaxis()->SetLabelFont(42);
329  mg->GetHistogram()->GetYaxis()->SetLabelOffset(0.01);
330  mg->GetHistogram()->GetYaxis()->SetTitleOffset(1.5);
331  mg->GetHistogram()->GetYaxis()->SetRangeUser(0.5*gFAP_MIN,1.5*gFAP_MAX);
332 
333  mg->GetXaxis()->SetLabelFont(42);
334  mg->GetYaxis()->SetLabelFont(42);
335  mg->GetXaxis()->SetTitleFont(42);
336  mg->GetYaxis()->SetTitleFont(42);
337  mg->GetXaxis()->SetTitleOffset(1.20);
338  mg->GetYaxis()->SetTitleOffset(1.05);
339  mg->GetXaxis()->SetTitleSize(0.04);
340  mg->GetYaxis()->SetTitleSize(0.04);
341 #ifdef XAXIS_IS_IFAD
342  mg->GetXaxis()->SetTitle("IFAD");
343 #else
344  mg->GetXaxis()->SetTitle("rho");
345 #endif
346  mg->GetYaxis()->SetTitle("FAP");
347 
348  mg->GetXaxis()->SetMoreLogLabels();
349 
350  mg->Draw("ALP");
351 
352  // draw the legend
353  TLegend *leg;
354  double hleg = 0.84-nFAP*0.03;
355  leg = new TLegend(0.7369062,hleg,0.9914738,0.9265385,NULL,"brNDC");
356  leg->SetBorderSize(1);
357  leg->SetTextAlign(22);
358  leg->SetTextFont(12);
359  leg->SetLineColor(1);
360  leg->SetLineStyle(1);
361  leg->SetLineWidth(1);
362  leg->SetFillColor(0);
363  leg->SetFillStyle(1001);
364  leg->SetTextSize(0.03);
365  leg->SetLineColor(kBlack);
366  leg->SetFillColor(kWhite);
367 
368  for(int i=0;i<nFAP;i++) leg->AddEntry(gr[i],FAP_LABEL[i].Data(),"lp");
369  if(refline>0) leg->AddEntry(gr[nFAP],reflabel.Data(),"lp");
370 #ifdef SHOW_3SIGMA
371  leg->AddEntry(gr[nFAP+1],"3 sigma","lp");
372 #endif
373  leg->Draw();
374 
375  // save plot
376  if(pfile!="") {
377  TString gfile=pfile;
378  gfile.ReplaceAll(".png",".gif");
379  gCANVAS->Print(gfile);
380  char cmd[1024];
381  sprintf(cmd,"convert %s %s",gfile.Data(),pfile.Data());
382  cout << cmd << endl;
383  gSystem->Exec(cmd);
384  sprintf(cmd,"rm %s",gfile.Data());
385  gSystem->Exec(cmd);
386  exit(0);
387  }
388 }
389 
390 double GetFAP(double rho, int n, int nFAP, TString* FAXvsRHO, TString* PRODvsFAX) {
391 
392  double RHO = rho;
393  double FAX = CWB::Toolbox::GetStepFunction("Y",FAXvsRHO[n],RHO);
394  double PROD = 0;
395  for(int j=0;j<nFAP;j++) PROD += CWB::Toolbox::GetStepFunction("Y",PRODvsFAX[j],FAX);
396  double MU = FAX*PROD;
397  double FAP = 1.-exp(-MU);
398 
399  return FAP;
400 }
401 
402 void GetFAPvsRHO(int n, vector<double>& x, vector<double>& y, int nFAP, TString* FAXvsRHO, TString* PRODvsFAX) {
403 
404  x.clear();y.clear();
405 
406  double rho_min = CWB::Toolbox::GetStepFunction("xmin",FAXvsRHO[n]);
407  double rho_max = CWB::Toolbox::GetStepFunction("xmax",FAXvsRHO[n]);
408  if((gRHO_MIN>0)&&(rho_min<gRHO_MIN)) rho_min = gRHO_MIN;
409  if((gRHO_MAX>0)&&(rho_max>gRHO_MAX)) rho_max = gRHO_MAX;
410 
411  double drho = 0.1;
412  int nrho = TMath::Nint((rho_max-rho_min)/drho)+2;
413 
414  for(int i=0;i<nrho;i++) {
415  double RHO = rho_min+i*drho;
416  double FAX = CWB::Toolbox::GetStepFunction("Y",FAXvsRHO[n],RHO);
417  double PROD = 0;
418  for(int j=0;j<nFAP;j++) PROD += CWB::Toolbox::GetStepFunction("Y",PRODvsFAX[j],FAX);
419  double MU = FAX*PROD;
420  double FAP = 1.-exp(-MU);
421 
422 #ifdef XAXIS_IS_IFAX
423  x.push_back(1./FAX);
424 #else
425  x.push_back(RHO);
426 #endif
427  y.push_back(FAP);
428 
429  if((FAP!=0)&&(FAP<gFAP_MIN)) gFAP_MIN = FAP; // skip last y[i]=0 to avoid TGraph log issue
430  if(FAP>gFAP_MAX) gFAP_MAX = FAP;
431  }
432 }
433 
434 void GetFAPvsRHO(vector<double>& x, vector<double>& y, int nFAP, TString* FARvsRHO, TString* OBSTIMEvsFAR) {
435 
436  x.clear();y.clear();
437 
438  double rho_min = 1e10;
439  double rho_max = 0;
440  for(int j=0;j<nFAP;j++) {
441  double _rho_min = CWB::Toolbox::GetStepFunction("xmin",FARvsRHO[j]);
442  double _rho_max = CWB::Toolbox::GetStepFunction("xmax",FARvsRHO[j]);
443  if(_rho_min<rho_min) rho_min=_rho_min;
444  if(_rho_max>rho_max) rho_max=_rho_max;
445  }
446  if((gRHO_MIN>0)&&(rho_min<gRHO_MIN)) rho_min = gRHO_MIN;
447  if((gRHO_MAX>0)&&(rho_max>gRHO_MAX)) rho_max = gRHO_MAX;
448 
449  double drho = 0.1;
450  int nrho = TMath::Nint((rho_max-rho_min)/drho)+2;
451 
452  double MIN_FAR[MAX_FAP];
453  for(int j=0;j<MAX_FAP;j++) MIN_FAR[j]=1e10;
454  for(int i=0;i<nrho;i++) {
455  double RHO = rho_min+i*drho;
456  double MU = 0;
457  for(int j=0;j<nFAP;j++) {
458  double FAR = CWB::Toolbox::GetStepFunction("Y",FARvsRHO[j],RHO);
459  double OBSTIME = CWB::Toolbox::GetStepFunction("Y",OBSTIMEvsFAR[j],FAR);
460  if(FAR>0) if(FAR<MIN_FAR[j]) MIN_FAR[j]=FAR; // FAR min hold
461  MU += MIN_FAR[j]*OBSTIME;
462  }
463  double FAP = 1.-exp(-MU);
464 
465  x.push_back(RHO);
466  y.push_back(FAP);
467 
468  if((FAP!=0)&&(FAP<gFAP_MIN)) gFAP_MIN = FAP; // skip last y[i]=0 to avoid TGraph log issue
469  if(FAP>gFAP_MAX) gFAP_MAX = FAP;
470  }
471 }
472 
473 TGraph*
474 GetGraph(vector<double> x, vector<double> y, TString ptitle, TString xtitle, TString ytitle) {
475 
476  int size = x.size();
477 
478  TGraph* gr = new TGraph;
479 
480  int cnt=0;
481 
482  for(int i=1;i<size;i++) {
483  if(y[i]==0) continue; // skip y[i]=0 to avoid TGraph log issue
484  double dx = (x[i]-x[i-1])/10000.;
485  gr->SetPoint(cnt++,x[i]-dx,y[i-1]);
486  gr->SetPoint(cnt++,x[i]+dx,y[i]);
487  }
488 
489  gr->GetHistogram()->GetXaxis()->SetTitle(xtitle.Data());
490  gr->GetHistogram()->GetYaxis()->SetTitle(ytitle.Data());
491  gr->GetHistogram()->GetYaxis()->SetTitleOffset(1.4);
492  // add more log labels when y range is small
493  double max=0;double min=1e80;
494  for(int i=0;i<size;i++) {
495  if(y[i]==0) continue; // skip y[i]=0 to avoid TGraph log issue
496  if(y[i]>max) max=y[i]; if(y[i]<min && y[i]!=0) min=y[i];
497  }
498  if(max/min<10) {
499  gr->GetHistogram()->GetYaxis()->SetMoreLogLabels();
500  }
501  gr->SetTitle(ptitle);
502  gr->SetLineColor(kRed);
503  gr->SetLineWidth(2);
504 
505  return gr;
506 }
TGraph * GetGraph(vector< double > x, vector< double > y, TString ptitle, TString xtitle, TString ytitle)
Definition: cwb_mkcfad.C:474
double rho
vector< double > FAR
Definition: cwb_mkfad.C:77
char xtitle[1024]
double min(double x, double y)
Definition: eBBH.cc:31
#define PROB_3SIGMA
Definition: cwb_mkcfad.C:72
int n
Definition: cwb_net.C:28
double rho_min
TString("c")
double GetFAP(double rho, int n, int nFAP, TString *FAXvsRHO, TString *PRODvsFAX)
Definition: cwb_mkcfad.C:390
double gRHO_MIN
Definition: cwb_mkcfad.C:89
cout<< endl;cout<< "ts size = "<< ts.size()<< " ts rate = "<< ts.rate()<< endl;tf.Forward(ts, wdm);int levels=tf.getLevel();cout<< "tf size = "<< tf.size()<< endl;double dF=tf.resolution();double dT=1./(2 *dF);cout<< "rate(hz) : "<< RATE<< "\ layers : "<< nLAYERS<< "\ dF(hz) : "<< dF<< "\ dT(ms) : "<< dT *1000.<< endl;int itime=TIME_PIXEL_INDEX;int ifreq=FREQ_PIXEL_INDEX;int index=(levels+1) *itime+ifreq;double time=itime *dT;double freq=(ifreq >0) ? ifreq *dF :dF/4;cout<< endl;cout<< "PIXEL TIME = "<< time<< " sec "<< endl;cout<< "PIXEL FREQ = "<< freq<< " Hz "<< endl;cout<< endl;wavearray< double > x
vector< double > RHO
Definition: cwb_mkfad.C:73
void cwb_mkcfad(TString fapList, TString fapMode, double rho_min=0, double rho_max=0, double refline=0, TString reflabel="reference", TString pfile="")
Definition: cwb_mkcfad.C:95
Long_t size
int j
Definition: cwb_net.C:28
#define MAX_FAP
Definition: cwb_mkcfad.C:70
i drho i
char str[1024]
TGraph * gr
double FAP
double rho_max
double gFAP_MIN
Definition: cwb_mkcfad.C:92
char title[256]
Definition: SSeriesExample.C:1
char subtitle[1024]
double gRHO_MAX
Definition: cwb_mkcfad.C:90
ifstream in
static int GetStepFunction(TString fName, vector< double > &x, vector< double > &y, vector< double > &ex=DEFAULT_DOUBLE_VECTOR, vector< double > &ey=DEFAULT_DOUBLE_VECTOR)
Definition: Toolbox.cc:6429
TCanvas * gCANVAS
Definition: cwb_mkfad.C:61
char cmd[1024]
int cnt
sprintf(tfres,"(1/%g)x(%g) (sec)x(Hz)", 2 *df, df)
TString gfile
void GetFAPvsRHO(int n, vector< double > &x, vector< double > &y, int nFAP, TString *FAXvsRHO, TString *PRODvsFAX)
Definition: cwb_mkcfad.C:402
char line[1024]
wavearray< double > y
Definition: Test10.C:31
TMultiGraph * mg
double gFAP_MAX
Definition: cwb_mkcfad.C:93
exit(0)