StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
plPi0Solo.C
1 TFile *fd=0;
2 TCanvas *c;
3 class TPad;
4 TPad *pad=0;
5 
6 //==========================
7 //==========================
8 
9 void plPi0Solo( int page=1,TString fName="R7098001" ) {
10  fName="sumP";
11  int pr=0;
12 
13  TString outDir="./out1/";
14 
15  gROOT->Reset();
16  gStyle->SetPalette(1,0);
17 
18  // Connect the input file and get the histogram
19  TString name2=outDir+fName+".hist.root";
20  printf("opening '%s'\n",name2.Data());
21  TFile *fd=new TFile(name2);
22  if(!fd->IsOpen()) {
23  printf("NOT found file=%s, STOP\n",hName.Data());
24  return ;
25  }
26  assert(fd->IsOpen());
27  // fd->ls();
28 
29  float xLo=.05, xHi=.22;
30  // float xLo=.4, xHi=.7;
31 
32 
33  TString ctit="pi0-"; ctit+=page;
34  if(page==3 || page==6)
35  c=new TCanvas(ctit,ctit,900,950);
36  else
37  c=new TCanvas(ctit,ctit,900,650);
38 
39  c->Range(0,0,1,1);
40  TPad *pad0 = new TPad("pad0", "apd0",0.0,0.95,1.,1.);
41  pad0->Draw();
42  pad0->cd();
43 
44  TPaveText *pt = new TPaveText(0,0.,1,1,"br");
45  pt->Draw();
46  TDatime dt;
47  TString txt2=fName+", page="+page+", ";
48  txt2+=dt.AsString();
49  pt->AddText(txt2);
50  txt2="--";
51  pt->AddText(txt2);
52 
53  c->cd();
54  pad = new TPad("pad1", "apd1",0.0,0.0,1,.95);
55  pad->Draw();
56 
57  gStyle->SetOptStat(11111 );
58 
59  switch(page) {
60  //..........................................
61  //..........................................
62  case 1: {
63  pad->Divide(1,1);
64  pad->cd(1);
65  printYield(xLo,xHi);
66  plotMix("invm", xLo,xHi, fd);
67  } break;
68  //..........................................
69  //..........................................
70  case 2: {
71  char *name[]={"ctbSum","tE","cE","cSh","sE","cN"};
72  int log[]={ 1,1,1,0,1,0};
73  int i;
74  pad->Divide(4,2);
75  for(i=0;i<6;i++){
76  pad->cd(i+1);
77  TH1* h= (TH1*)fd->Get(name[i]); assert(h);
78  h->Draw();
79  if(log[i]) gPad->SetLogy();
80  }
81 
82  pad->cd(7);
83  plotMix("invm", xLo,xHi, fd);
84  } break;
85  //..........................................
86  //..........................................
87  case 3: {
88  char *name[]={"xyL","xyH","0eta","eta12","0xyH","invm","0ener","0phi","0Z"};
89  int i;
90  pad->Divide(3,3);
91  for(i=0;i<9;i++){
92  pad->cd(i+1);
93  TH1* h= (TH1*)fd->Get(name[i]); assert(h);
94  if(i<6 && i%3<2 ) { h->Draw("colz"); continue; }
95  plotMix(name[i], xLo,xHi, fd);
96  }
97  } break;
98 
99  //..........................................
100  //..........................................
101  case 4: {
102  pad->Divide(2,3);
103  int i;
104  char *namex[]={"0ener","0eta","0phi","0pt","0Z","0Ang"};
105  for(i=0;i<6;i++){
106  pad->cd(i+1);
107  plotMix(namex[i], xLo,xHi, fd);
108  }
109  } break;
110 
111  //..........................................
112  //..........................................
113  case 5: {
114  pad->Divide(1,2);
115  pad->cd(1);
116  TH1* h= (TH1*)fd->Get("ytw"); assert(h);
117  h->Draw();
118 
119  pad->cd(2);
120  plotMix("0ytw", xLo,xHi, fd);
121 
122  } break;
123 
124  //..........................................
125  //..........................................
126  case 6: {
127  gStyle->SetOptStat(1);
128  TGraphErrors* gr1=new TGraphErrors;
129  gr1->SetMarkerStyle(21);
130  gr1->SetName(ctit+"eM"); //
131  gr1->SetTitle("Pi0 mass vs. Eta; eta bin; inv mass (MeV);");
132 
133  TGraphErrors* gr2=new TGraphErrors;
134  gr2->SetMarkerStyle(20);
135  gr2->SetName(ctit+"eS"); //
136  gr2->SetTitle("Pi0 mass width vs. Eta; eta bin; sigma (MeV);");
137 
138  TGraphErrors* gr3=new TGraphErrors;
139  gr3->SetMarkerStyle(19);
140  gr3->SetName(ctit+"eN"); //
141  gr3->SetTitle("Pi0 Yield vs. Eta; eta bin; yield, err=nBckg/10;");
142 
143  pad->Divide(3,4);
144  int i;
145  for(i=0;i<12;i++) {
146  char t1[100];
147  sprintf(t1,"invm%02d",i+1);
148  pad->cd(1+i);
149  plotMix(t1, xLo,xHi, fd,gr1,gr2,gr3);
150  }
151  printYield(xLo,xHi);
152 
153  c2=new TCanvas(ctit+"e",ctit+"e",400,700);
154  c2->Divide(1,3);
155  c2->cd(1); gr1->Draw("AP");
156  gr1->SetMinimum(100);gr1->SetMaximum(180);
157  TLine *ln=new TLine(0,135,13,135);
158  ln->SetLineColor(kBlue); ln->Draw();
159  c2->cd(2); gr2->Draw("AP");
160  gr2->SetMinimum(20);gr2->SetMaximum(60);
161  c2->cd(3); gr3->Draw("AP");
162  } break;
163  //..........................................
164  //..........................................
165  default: printf("page %d not defined\n",page);
166  pr=0;
167  }
168 
169  if(pr) {
170  TString outF=fName+"_"+page+".ps";
171  //TString outF=fName+".ps";
172 
173  c->Print(outF);
174  outF.ReplaceAll(".ps",".gif");
175  c->Print(outF);
176  return;
177  }
178 
179  return;
180 }
181 
182 
183 //---------------------------------------------
184 //---------------------------------------------
185 //---------------------------------------------
186  void printYield(float xLo, float xHi) {
187  ax=invm->GetXaxis();
188  int bLo=ax->FindBin(xLo);
189  int bHi=ax->FindBin(xHi);
190  int bin1=ax->FindBin(0.6);
191  int bin2=ax->FindBin(1.0);
192 
193  float nReal=invm->Integral(bin1,bin2);
194  float nMix=Xinvm->Integral(bin1,bin2);
195  float fac=-1;
196  if(nMix>0) fac=nReal/nMix;
197  printf("nReal=%d nMix=%d fac=%f\n",nReal,nMix,fac);
198 
199 
200  float s1=invm->Integral(bLo,bHi);
201  float s2=Xinvm->Integral(bLo,bHi);
202  float x=s1-s2*fac;
203  float ex=sqrt(s1+s2*fac*fac);
204  printf("s1=%f s2=%f\n",s1,s2);
205  float s2b=-1;
206  if(s2>0) s2b=x/s2/fac; // signal to background
207  //error calculated assuming s1 & s2 are independent
208  float es2b=-999;
209  if(s1>0 && s2>0 ) es2b=s2b*sqrt(1/s1+1/s2/fac/fac);
210  printf(" Nrec=%.0f Nmix=%.0f for invm=[%.2f, %.2f] , Npi0=%.0f + / - %.0f \n s2b=%.2f + / - %.3f \n",s1,s2,xLo,xHi, x,ex,s2b,es2b);
211  printf("#1 Nrec=%.0f, Nmix=%.0f mixFac=%.3f<br> Npi0=%.0f + / - %.0f <br> s2b=%.2f + / - %.3f \n",s1,s2, fac,x,ex,s2b,es2b);
212  }
213 
214 
215 //---------------------------------------------
216 //---------------------------------------------
217 //---------------------------------------------
218 void plotMix(TString hname="m1",float xLo, float xHi,TFile *fd,
219  TGraphErrors* grM=0, TGraphErrors* grS=0, TGraphErrors* grN=0) {
220  hx= (TH1*)fd->Get("invm"); assert(hx);
221  int bin1=hx->FindBin(0.6);
222  int bin2=hx->FindBin(1.0);
223  float nReal=hx->Integral(bin1,bin2);
224  hx= (TH1*)fd->Get("Xinvm"); assert(hx);
225  float nMix=hx->Integral(bin1,bin2);
226  float fac=nReal/nMix;
227 
228  printf("total: nReal=%d nMix=%d fac=%f \n",nReal,nMix,fac);
229 
230  int i=1;
231  TH1* hr= (TH1*)fd->Get(hname); assert(hr);
232  TH1* hm= (TH1*)fd->Get("X"+hname); assert(hm);
233  hm->Sumw2();
234  hm->Scale(fac);
235 
236  hm->SetLineColor(kGreen);
237  // hm->SetLineStyle(2);
238  TH1F* hD=(TH1F*)hr->Clone();
239  hD->Sumw2();
240  hD->SetTitle(hname+" Real-Mixed");
241  hD->Add(hm,-1);
242  hD->SetLineColor(kBlue);
243  hD->Draw();
244  hD->SetMinimum(-0.2 *hr->GetMaximum());
245  hD->SetMaximum(1.1*hr->GetMaximum());
246  hr->Draw("same");
247  hm->Draw("same");
248 
249  // gStyle->SetOptStat(11);
250  mzer=new TLine(-10.,0,10,0.);
251  mzer->Draw();
252 
253  // printf("total: nReal=%d nMix=%d fac=%f name='%s'\n",nReal,nMix,fac,hname.Data());
254 
255  if(strstr(hname.Data(),"invm")==0) return;
256 
257  hD->SetTitle(" Real-Mixed; inv mass GeV");
258 
259  hD->Fit("gaus","R","same",xLo,xHi);
260  hD->Draw("same");
261  mLo=new TLine(xLo,0,xLo,3500);
262  mLo->SetLineColor(kMagenta);
263  mHi=new TLine(xHi,0,xHi,3500);
264  mHi->SetLineColor(kMagenta);
265  mLo->Draw();
266  mHi->Draw();
267 
268  f=hD->GetFunction("gaus");
269  f->SetLineColor(kRed);
270  f->SetLineWidth(1.);
271 
272  double *par=f->GetParameters();
273  double *epar=f->GetParErrors();
274  float pi0m=par[1]*1000; // in MeV
275  float epi0m=epar[1]*1000; // in MeV
276  float pi0s=par[2]*1000; // in MeV
277  float epi0s=epar[2]*1000; // in MeV
278 
279  printf("#2 mass/MeV=%.0f + / - %.0f , <br>sig/MeV %.0f + / - %.0f\n",pi0m,epi0m,pi0s,epi0s);
280 
281 
282  if(grM==0) return;
283  int etaB=atoi(hname.Data()+4);
284 
285  int n=grM->GetN();
286  grM->SetPoint(n,etaB,pi0m);
287  grM->SetPointError(n,0.5,epi0m);
288 
289  n=grS->GetN();
290  grS->SetPoint(n,etaB,pi0s);
291  grS->SetPointError(n,0.5,epi0s);
292 
293  ax=invm->GetXaxis();
294  int bLo=ax->FindBin(xLo);
295  int bHi=ax->FindBin(xHi);
296  float nPi0=hD->Integral(bLo,bHi);
297  float nBckg=hr->Integral(bLo,bHi) - nPi0;
298  printf("etaB=%d nPi0=%.1f nBckg=%.1f\n",etaB,nPi0,nBckg);
299  n=grN->GetN();
300  grN->SetPoint(n,etaB,nPi0);
301  grN->SetPointError(n,0.5,0. );
302  return;
303  // TH1F* hR=(TH1F*)hr->Clone();
304 }
305 
Generated by   doxygen 1.8.5