StRoot: StEEmcPool/StBbcVertex/macros/plBbc.C Source File

00001 TFile *fd=0;
00002 TCanvas *c;
00003 class TPad;
00004 TPad *pad=0;
00005 
00006 
00007 //==========================
00008 //==========================
00009 
00010 void plBbc( int page=21,TString fName="aa3" ) {
00011   int pr=0;
00012 
00013   TString outDir="../outBbc/";
00014  
00015   gROOT->Reset();
00016   gStyle->SetOptStat(1111111);
00017   gStyle->SetOptFit(1111111);
00018 
00019   gStyle->SetPalette(1,0);
00020   
00021   // Connect the input file and get the  histogram
00022   TString name2=outDir+fName+".hist.root";
00023   printf("opening '%s'\n",name2.Data());
00024   fd=new TFile(name2);
00025   if(!fd->IsOpen()) {
00026     printf("NOT found file=%s, STOP\n",hName.Data());
00027     return ;
00028   }
00029   assert(fd->IsOpen());
00030   //  fd->ls();
00031    
00032   TString ctit=fName+"_"; ctit+=page;
00033   c=new TCanvas(ctit,ctit,450,350);
00034 
00035 #if 0  
00036   c->Range(0,0,1,1);
00037   TPad *pad0 = new TPad("pad0", "apd0",0.0,0.95,1.,1.);
00038   pad0->Draw();
00039   pad0->cd();
00040 
00041   TPaveText *pt = new TPaveText(0,0.,1,1,"br");
00042   pt->Draw();
00043   TDatime dt;
00044   TString txt2=fName+", page="+page+",  ";
00045   txt2+=dt.AsString();
00046   pt->AddText(txt2);
00047   txt2="--";
00048   pt->AddText(txt2);
00049 
00050   c->cd();
00051   pad = new TPad("pad1", "apd1",0.0,0.0,1,.95);
00052   pad->Draw();
00053 #endif
00054   pad=c;
00055 
00056  //fd->ls();
00057   gStyle->SetOptStat(11111     );
00058 
00059   switch(page) {
00060   //..........................................
00061   //..........................................
00062   case 0: {
00063     pad->Divide(1,1);
00064     pad->cd(1);
00065     TH1* h= (TH1*)fd->Get("vd"); assert(h);
00066     h->Fit("gaus","R","",-60,60);
00067   }  break;
00068   
00069   //..........................................
00070   //..........................................
00071   case 1: {
00072     pad->Divide(2,2);
00073     int i;
00074     char *namex[]={"vk","va","vt","vv"};
00075     for(i=0;i<4;i++){
00076       pad->cd(i+1);
00077       TH1* h= (TH1*)fd->Get(namex[i]); assert(h);
00078       h->Draw("colz");
00079     }
00080   }  break;
00081   //..........................................
00082   //..........................................
00083   case 2: {
00084     pad->Divide(2,3);
00085     int i;
00086     char *namex[]={"vk","va","vt","vE","vv","vd","vW"};
00087     for(i=0;i<7;i++){
00088       pad->cd(i+1);
00089       TH1* h= (TH1*)fd->Get(namex[i]); assert(h);
00090       if(i!=5)h->Draw("box");
00091       else h->Fit("gaus");
00092     }
00093   }  break;
00094   //..........................................
00095   //..........................................
00096  case 3: {
00097    pad->Divide(1,1);
00098    pad->cd(1);
00099    TString tt="vv";
00100    TH2F*h2= (TH2F*)gROOT->FindObject(tt);
00101    assert(h2);
00102    h2->FitSlicesY();
00103    TH1D*hd= (TH1D*)gROOT->FindObject(tt+"_1");
00104 
00105    h2->Draw("colz");
00106    hd->SetMarkerStyle(20);
00107    char *cFunc="pol1";
00108    hd->Fit(cFunc,"R","same",220,280);
00109    TF1* f1=hd->GetFunction(cFunc);
00110    f1->SetLineColor(kYellow);
00111   f1->SetLineWidth(2.);
00112  } break;
00113 
00114   //..........................................
00115   //..........................................
00116   case 21: 
00117   case 22: 
00118     calT(pad,fd,page-21); // time offset
00119     break;
00120   //..........................................
00121   //..........................................
00122   case 31: 
00123   case 32: 
00124     calW(pad,fd,page-31); // time walk
00125     break;
00126   //..........................................
00127   //..........................................
00128   default: printf("page %d not defined\n",page);
00129     pr=0;
00130   }
00131 
00132   if(pr) {
00133     TString outF=fName+"_"+page+".ps";
00134     //TString outF=fName+".ps";
00135 
00136     c->Print(outF);
00137     outF.ReplaceAll(".ps",".gif");
00138     c->Print(outF);
00139     return;
00140   }
00141 
00142   return;
00143 }
00144 
00145 //==========================
00146 //==========================
00147 calT(TPad *pad, TFile *fd, int iew=0) {
00148   char cEW[2]={'E','W'};
00149   int pm=-1+2*iew;
00150   int nok=0;
00151   ln0=new TLine(0,0,0,200); ln0->SetLineColor(kGreen);
00152   pad->Divide(3,4);
00153   for(int ih=0;ih<16;ih++) {
00154     char tt1[100];
00155     int iT=ih;
00156     sprintf(tt1,"w%c%02d",cEW[iew],iT+1);
00157     pad->cd(ih+1);
00158     TH2* h2= (TH2*)fd->Get(tt1);
00159     if(h2==0) break;
00160     assert(h2);
00161     // h2->Draw("box");
00162     printf("pmt=%s neve=%.0f mean=%f\n",tt1,h2->GetEntries(),h2->GetMean(1));
00163     if(h2->GetEntries()<50) continue;
00164     TH1D *hd=h2->ProjectionY(); 
00165     hd->Fit("gaus");
00166     ln0->Draw();
00167     f=hd->GetFunction("gaus");
00168     f->SetLineColor(kRed);
00169     f->SetLineWidth(2.);
00170     nok++;
00171     sprintf(tt1,"bbc%c%02d",cEW[iew],iT+1);    
00172     printf("#2 pmt=%s dz/cm=%.1f + / - %.1f sigDZ/cm %.1f + / - %.1f\n",tt1,f->GetParameter(1),f->GetParError(1),f->GetParameter(2),f->GetParError(2));
00173     // output in TDC
00174     float fac=2;  //(1chan=2 cm)
00175     printf("#1 %s %.1f + / - %.1f  N=%.0f\n",tt1,pm*f->GetParameter(1)/fac,f->GetParError(1)/fac, h2->GetEntries());
00176   }
00177     printf("#1# found %d pmt in bbc%c, sign= %.1f\n",nok,cEW[iew],pm);
00178 
00179 }
00180 
00181 //==========================
00182 //==========================
00183 calW(TPad *pad, TFile *fd, int iew=0) {// time walk
00184   char cEW[2]={'E','W'};
00185    int nok=0;
00186   int pm=-1+2*iew;
00187   pad->Divide(1,1);
00188   
00189   TH2F *hs=0;
00190   for(int ih=0;ih<16 ;ih++) {
00191     char tt1[100];
00192     int iT=ih;
00193     //  if(iT<6) continue;
00194     sprintf(tt1,"w%c%02d",cEW[iew],iT+1);
00195     pad->cd(ih+1);
00196     TH2* h2= (TH2*)fd->Get(tt1); assert(h2);
00197     // h2->Draw("box");
00198     printf("pmt=%s neve=%.0f mean=%f\n",tt1,h2->GetEntries(),h2->GetMean(1));
00199     if(h2->GetEntries()<400) continue;
00200     twFit(h2);
00201     TString tt=h2->GetName();
00202     pad->Print(tt+".gif");
00203    nok++;
00204   }
00205   printf("#1# added %d pmt in bbc%c\n",nok,cEW[iew]);
00206 }
00207 
00208 //==========================
00209 //==========================
00210 twFit( TH2* h2) {
00211 
00212   ln0=new TLine(0,0,200,0); ln0->SetLineColor(kGreen);
00213   TString tt=h2->GetName();
00214 
00215   TF1 *f1 = new TF1(tt+"fit",timeWalkFunc,-100,120,4);
00216   // f1->SetParameters(-25,50,0.01 1.5);// East
00217    f1->SetParameters(-50,40,0.01 0.1);// East #8
00218   //f1->SetParameters(-80,30,1,0.1);// West
00219   f1->SetLineColor(kRed);
00220   f1->SetLineWidth(2.);
00221     
00222   h2->FitSlicesY();
00223   TH1D*hd= (TH1D*)gROOT->FindObject(tt+"_1");
00224   hd->SetMaximum(100);
00225   hd->SetMinimum(-100);
00226   hd->SetMarkerStyle(20);
00227   h2->Draw("colz");
00228   hd->Fit(tt+"fit","R","same",4,200);
00229 
00230   ln0->Draw();
00231 
00232   //........ draw +/- 1 sig error band for projection
00233   Stat_t er1[1000];
00234   memset(er1,0,sizeof(er1));
00235 
00236   TH1D *he= (TH1D*)gROOT->FindObject(tt+"_2");
00237   TH1D *yu=(TH1D*)  hd->Clone(); yu->SetLineColor(kMagenta);
00238   yu->Add(he);  yu->SetError(er1);  
00239   yu->Draw("same");
00240   
00241   TH1D *yd=(TH1D*)  hd->Clone();
00242   yd->Add(he,-1);  yd->SetError(er1);   yd->SetLineColor(kMagenta);
00243   yd->Draw("same");
00244 
00245   double *par=f1->GetParameters();
00246   double *parE=f1->GetParErrors();
00247   int i;
00248 
00249   float fac=2.16; // 1tdc ch = 2 cm
00250   float pm= 1; // use - for E
00251   printf("#1 bbc%s %g %g %g %g \n",h2->GetName()+1,pm*par[0]/fac,pm*par[1]/fac,par[2],par[3]);    
00252 }
00253 
00254 //==========================
00255 //==========================
00256 
00257 Double_t timeWalkFunc(Double_t *x, Double_t *par){
00258   float maxAdc=80;
00259   Float_t xx =x[0];
00260   if(xx>maxAdc) xx=maxAdc; // ADC cut-off
00261   Double_t f = par[0] + par[1]*exp(-par[2]*pow(xx,par[3]));
00262   //Double_t f = par[0] + par[1]*xx;
00263   return f;
00264 }
00265 
00266 
00267