make_cm_pullterm_pdf.C
Go to the documentation of this file.
1 #include "CAFAna/Analysis/CovMxSurface.h"
3 #include <fstream>
4 
6  //std::cout << "SystToName: " << syst << std::endl;
7  if(syst=="stats_only") return "Stats Only";
8  if(syst=="all_covmx_systs") return "All CM Systs";
9  if(syst=="all_pull_terms") return "All Pull Terms";
10  std::string det="ND";// = (syst.find("_neardet_")!=std::string::npos) ? "ND" : "FD";
11  std::string pol="FHC";// = (syst.find("_fhc_")!=std::string::npos) ? "FHC" : "RHC";
12  std::string comb, sy;
13  if(syst.find("single_covmx")!=std::string::npos) comb = "CM: ";
14  if(syst.find("single_pull_term")!=std::string::npos) comb = "Pull: ";
15  sy = syst.substr(syst.find("/")+1);
16  //return Form("%s%s %s",comb.c_str(),det.c_str(),sy.c_str());
17  return Form("%s %s",comb.c_str(),sy.c_str());
18 }
19 
21  if(label=="single_syst") return "Single Syst";
22  if(label=="swap_syst") return "Swap Pull";
23  return "";
24 }
25 
27  std::string systpath, std::string parameters,
28  std::string pdfpath){
29  /*
30  directories look like:
31  th24vsdm41_nus_fhc_neardet_covmxsyst=AbsCalib_noprofdelta24
32  th24vsdm41_nus_fhc_neardet_pullsyst=AbsCalib_noprofdelta24
33  th24vsdm41_nus_fhc_neardet_covmxsysts_noprofdelta24
34  th24vsdm41_nus_fhc_neardet_covmxsysts_pullsyst=AbsCalib_noprofdelta24
35  */
36 
37  std::string fromfile = path +"/"+ space +"/"+ systpath + "/merged.root";
38 
39  std::string fromdir = space + parameters + "_noprofdelta24";
40 
41  std::cout << " From File: " << fromfile << std::endl;
42  std::cout << " From Dir: " << fromdir << std::endl;
43 
44  std::string fittype = systpath;
45  std::string syst;
46  if(systpath.find("/")!=std::string::npos){
47  fittype = systpath.substr(0,systpath.find("/"));
48  syst = systpath.substr(systpath.find("/")+1);
49  }
50 
51  TFile * f = new TFile(fromfile.c_str());
52  TCanvas * canv = new TCanvas("","");
53  canv->SetLogx(true);
54  canv->SetLogy(true);
55  gStyle->SetPalette(kBird);
56  TH2D * h2;
57  //f->cd(Form("%s/Fit vars",fromdir.c_str())); f->ls();
58  for(auto& var : {"dmsq32","ssth34","th23"} ){
59  f->GetObject((fromdir+"/Fit vars/"+var).c_str(),h2);
60  if(!h2){
61  std::cout << "Did not find " << fromdir+"/Fit vars/"+var << std::endl;
62  } else {
63  h2->Draw("colz");
64  canv->Print((pdfpath+"/"+space+"/"+fittype+"/fitvars/"+syst+"_fit-"+var+".pdf").c_str());
65  canv->Clear();
66  }
67  }
68  if(parameters.find("pullsyst=")!=std::string::npos){
69  //f->cd(Form("%s/Syst pulls",fromdir.c_str())); //f->ls();
70  f->GetObject((fromdir+"/Syst pulls/nus_fhc_neardet_"+syst).c_str(),h2);
71  h2->Draw("colz");
72  canv->Print((pdfpath+"/"+space+"/"+fittype+"/systpull/"+syst+"_systpull.pdf").c_str());
73  canv->Clear();
74  }
75 
76  return ana::LoadFromFile<ana::CovMxSurface>( fromfile.c_str(),
77  fromdir.c_str()
78  ).release();
79 }
80 
81 
83  std::vector< ana::CovMxSurface* > surfs,
85  std::cout << "ContourComparisonPlot: " << space+"_"+label+"_"+syst << std::endl;
86  TCanvas * c = new TCanvas((space+"_"+label+"_"+syst).c_str(),"c");
87  std::vector<double> up = {0.68,0.95};//,0.997};
88 
89  TLegend* legend = new TLegend(0.15, 0.2, 0.5, 0.45);
90  legend->SetFillStyle(0);
91  legend->SetHeader(Header(label).c_str(),"C");
92  int col_interval = (int) gStyle->GetNumberOfColors() / surfs.size();
93 
94  // Each surface in the combo
95  for(unsigned int d=0; d<surfs.size(); d++){
96  ana::CovMxSurface*surf = surfs[d];
97  if(d==0){
98  c->cd();
99  TH2D* h = (TH2D*)surf->GetSurface()->Clone();
100  h->SetTitle("");
101  h->GetYaxis()->SetTitle("#Deltam^{2}_{41} (eV^{2})");
102  h->GetXaxis()->SetTitle("sin^{2}(#theta_{24})");
103  h->GetXaxis()->SetRangeUser(1e-4,1);
104  h->Draw("axis");
105  }
106 
107  // set colors
108  int col = gStyle->GetColorPalette(d*col_interval);
109  if (d==0) col = kBlack;
110  if (d==1) col = kRed-4;
111  if (d==2) col = kBlue-4;
112 
113  std::vector< std::vector<TGraph*> > graphs = surf->GetContours(up);
114  for( size_t sig_i=0; sig_i<graphs.size(); sig_i++ ){
115  for(size_t j=0; j<graphs[sig_i].size(); j++){
116  TGraph* gg = (TGraph*)graphs[sig_i][j]->Clone();
117  if(!gg){ std::cout << "TGraph not found, j = " << j << std::endl; }
118  gg->SetLineColor(col);
119  gg->SetLineStyle(sig_i+1);
120  if(j==0){
121  legend->AddEntry(gg,
122  Form("%s, %i#sigma",
123  SystToName(syst).c_str(),
124  (int)sig_i+1),"l");
125  }
126  c->cd();
127  gg->Draw("c same");
128  }
129  }
130  //delete surf;
131  } // all surfaces in the combo
132 
133  legend->Draw("same");
134  c->SetLogx(true);
135  c->SetLogy(true);
136 
137  return c;
138 }
139 
140 
142  if( -1 == std::system(( "mkdir -p " + path ).c_str() ) ){
143  std::cout << "Error creating directory: " << path << std::endl; exit(1); }
144 }
145 
146 
148  //gStyle->SetPalette(kRainBow);
149  gStyle->SetPalette(kBird);
150 
151  std::string path = "/pnfs/nova/scratch/users/talion/nus/"+project;
152  std::string pdfpath="/nova/app/users/talion/pdf/nus/"+project;
153 
154  std::map<std::string, std::map<std::string, std::string> > f;
155  std::vector< std::string > spaces = {
156  "th24vsdm41_nus_fhc_neardet"
157  //"th24vsdm41_fhc_neardet_fardet",
158  //"th24vsth34_dm41=0.5_fhc_neardet_fardet",
159  //"th24vsth34_dm41=5_fhc_neardet_fardet"
160  };
161 
162  std::ifstream systlist("systs.txt");
163  std::vector< std::string > all_systs;
164  std::string sys;
165  while( std::getline(systlist, sys) ){
166  all_systs.push_back(sys);
167  }
168 
169 
170  // output directories for saving contour
171  TFile * outf = new TFile(("contours_"+project+".root").c_str(),"recreate");
172  std::vector< TDirectory * > spacedir;
173  for( size_t i=0; i<spaces.size(); i++ ){
174  outf->cd();
175  spacedir.push_back(outf->mkdir(spaces[i].c_str()));
176  }
177 
178 
179  for(unsigned int a=0; a<spaces.size(); a++){
180  std::string space=spaces[a];
181  std::string single_syst_dir = pdfpath+"/"+space+"/single_syst";
182  std::string swap_syst_dir = pdfpath+"/"+space+"/swap_syst";
183  MakeDir(single_syst_dir);
184  MakeDir(swap_syst_dir);
185  for( auto& dir : {"fitvars"} ){
186  MakeDir(pdfpath+"/"+space+"/single_pull_term/"+dir);
187  MakeDir(pdfpath+"/"+space+"/single_covmx_syst/"+dir);
188  MakeDir(pdfpath+"/"+space+"/stats_only/"+dir);
189  MakeDir(pdfpath+"/"+space+"/all_covmx_systs/"+dir);
190  MakeDir(pdfpath+"/"+space+"/swap_single_pull_term/"+dir);
191  }
192  MakeDir(pdfpath+"/"+space+"/single_pull_term/systpull");
193  MakeDir(pdfpath+"/"+space+"/swap_single_pull_term/systpull");
194 
195  // stats_only
196  // all_covmx_systs
197  // all_pull_terms
198 
199  std::cout << all_systs.size() << " systs in space: " << space << std::endl;
200  unsigned int n_c=0;
201  for(unsigned int s_i=0; s_i < all_systs.size(); s_i++){
202  std::string syst = all_systs[s_i];
203 
204  // Single Syst
205  std::vector< ana::CovMxSurface* > surfs;
206  surfs.push_back( GetSurf(path,space,"stats_only","",pdfpath) );
207  surfs.push_back( GetSurf(path,space,
208  "single_pull_term/"+syst,
209  "_pullsyst="+syst,
210  pdfpath) );
211  surfs.push_back( GetSurf(path,space,
212  "single_covmx_syst/"+syst,
213  "_covmxsyst="+syst,
214  pdfpath) );
215  TCanvas* cc = ContourComparisonPlot( space, surfs, "single_syst", syst );
216  // pdfpath / space / combo / combo_syst_cont.pdf
217  cc->Print((single_syst_dir+"/single_syst_"+ syst + "_cont.pdf").c_str() );
218  //spacedir[a]->WriteTObject(cc, ("single_syst_"+syst).c_str());
219  cc->Clear();
220 
221  // All Syst + Pull Swap
222  surfs.clear();
223  surfs.push_back( GetSurf(path,space,"all_covmx_systs","_covmxsysts",pdfpath) );
224  surfs.push_back( GetSurf(path,space,
225  "swap_single_pull_term/"+syst,
226  "_covmxsysts_pullsyst="+syst,
227  pdfpath) );
228  cc = ContourComparisonPlot( space, surfs, "swap_syst", syst );
229  cc->Print((swap_syst_dir+"/swap_syst_"+ syst + "_cont.pdf").c_str() );
230  //spacedir[a]->WriteTObject(cc, ("swap_syst_"+syst).c_str());
231  cc->Clear();
232 
233 
234  } // all systs
235  } // all spaces
236  outf->Close();
237 }
system("rm -rf microbeam.root")
TCanvas * canv
const char * label
std::string Header(std::string label)
std::vector< std::vector< TGraph * > > GetContours(std::vector< double > up)
Get contours from surface.
Int_t col[ntarg]
Definition: Style.C:29
const double a
Float_t d
Definition: plot.C:236
TFile * outf
Definition: testXsec.C:51
const double j
Definition: BetheBloch.cxx:29
OStream cout
Definition: OStream.cxx:6
TH1F * h2
Definition: plot.C:45
void make_cm_pullterm_pdf(std::string project)
const std::string path
Definition: plot_BEN.C:43
ana::CovMxSurface * GetSurf(std::string path, std::string space, std::string systpath, std::string parameters, std::string pdfpath)
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
TDirectory * dir
Definition: macro.C:5
void cc()
Definition: test_ana.C:28
TCanvas * ContourComparisonPlot(std::string space, std::vector< ana::CovMxSurface * > surfs, std::string label, std::string syst)
TH2D * GetSurface()
Definition: CovMxSurface.h:40
void MakeDir(std::string path)
exit(0)
std::string SystToName(std::string syst)
Float_t e
Definition: plot.C:35
std::vector< std::string > all_systs
surf
Definition: demo5.py:35