00001 #include "commonmacro/common.h"
00002 #include "common/Name.cc"
00003 #include "commonmacro/histutil.h"
00004
00005 void compareVertex(const char* inName=
00006 "links/P01hi.minbias.2000.hist/hianalysis_1000.hist.root",
00007 const char* psDir="ps",
00008 int cut = 1,
00009 const char* outDir="./",
00010 const char* more = "embeddingfile",
00011 float extraValue = 0)
00012 {
00013
00014
00015 gSystem->Load("StHiMicroEvent");
00016 gSystem->Load("StHiMicroAnalysis");
00017
00018 TFile* inRoot = new TFile(inName);
00019 TFile* embedRoot=0;
00020 if(more){
00021 strcpy(name,more);
00022 }
00023 else if(extraValue>0){
00024
00025 TString s=inName;
00026 char* trigger = (s.Contains("central")) ? "central" :"minbias";
00027 sprintf(name,"~/wrk/assoc/links/P01hi.central.HighpT_piminus_%d.hist/assoc_cut%d.hist.root",(int)extraValue,cut);
00028 }
00029 else{
00030 cout << "Unknown embed file" << endl; return;
00031 }
00032 embedRoot = new TFile(name);
00033
00034 cout << "--------------------------" << endl;
00035 cout << "in name=" << inName << endl
00036 << "ps dir=" << psDir << endl
00037 << "cut=" << cut << endl
00038 << "embed name=" << name << endl;
00039 cout << "--------------------------" << endl;
00040
00041 if(!inRoot){
00042 cout << "cannot find the infile" << endl;
00043 return;
00044 }
00045 if(!embedRoot){
00046 cout << "cannot find " << more << endl;
00047 }
00048 TH1* h1[2]={0};
00049
00050 h1[0]=(TH1*)inRoot.Get("h1VtxZCentCut");
00051 h1[1]=(TH1*)embedRoot.Get("h1VtxZCentCut");
00052 if(!h1[0] || !h1[1]) { cout << "cannot find histo" << endl; return; }
00053
00054 int nrebin=10;
00055 h1[0]->Rebin(nrebin); h1[1]->Rebin(nrebin);
00056 TH1* ha, *hb;
00057
00058 char* type[] = {"real","embed"};
00059
00060 TCanvas c1("c1","c1",400,500);
00061
00062 float minZ=-80,maxZ=80;
00063 TText* t=new TText;
00064 gStyle->SetOptStat(0);
00065
00066 TLegend l(0.15,0.65,0.25,0.85);
00067 l.SetFillColor(kWhite); l.SetBorderSize(0);
00068
00069 int lowBin=h1[0]->GetXaxis()->FindBin(minZ);
00070 int upBin =h1[1]->GetXaxis()->FindBin(maxZ-0.0001);
00071
00072
00073 sprintf(title,"vertexZ (cut %d)",cut);
00074 Divide(&c1,1,2,title,inName);
00075 for(int i=0; i<2; i++){
00076 c1.cd(i+1); h1[i]->SetTitle(title); h1[i]->Draw();
00077 PrintMeanRms(h1[i],0.7,0.8);
00078 }
00079 Print(&c1,psDir,"vertexZ");
00080
00081
00082 sprintf(title,"vertex z scaled (cut %d)",cut);
00083 Divide(&c1,1,2,title,inName);
00084 c1.cd(1); ha=(TH1*)h1[0]->Clone(); hb=(TH1*)h1[1]->Clone();
00085 ha->SetTitle("scaled");
00086 ha->Scale(1./ha->Integral());
00087 hb->Scale(1./hb->Integral());
00088 ha->SetMarkerStyle(4); hb->SetMarkerStyle(8);
00089 ha->Draw("p"); hb->Draw("psame");
00090 l.AddEntry(ha,"real","p"); l.AddEntry(hb,"embed","p"); l.Draw();
00091
00092 c1.cd(2); ha=(TH1*)h1[0]->Clone(); hb=(TH1*)h1[1]->Clone();
00093 ha->SetTitle("weighted");
00094 cout << "min = " << ha->GetXaxis()->GetBinLowEdge(lowBin)
00095 << ",max = " << ha->GetXaxis()->GetBinUpEdge(upBin) << endl;
00096
00097 double realPeak = ha->Integral(lowBin,upBin);
00098 double realOut = ha->Integral()-realPeak;
00099 double realFrac = realPeak/realOut;
00100 double embedPeak = hb->Integral(lowBin,upBin);
00101 double embedOut = hb->Integral()-hb->Integral(lowBin,upBin);
00102 double embedFrac = embedPeak/embedOut;
00103 double weight=realFrac/embedFrac;
00104
00105 cout << "real frac="<<realFrac <<", embedFrac="<<embedFrac
00106 << ", weight = " << weight << endl;
00107 for(int i=lowBin;i<=upBin;i++){
00108 hb->SetBinContent(i,hb->GetBinContent(i)*weight);
00109
00110 }
00111 ha->Scale(1./ha->Integral()); hb->Scale(1./hb->Integral());
00112 ha->SetMarkerStyle(4); hb->SetMarkerStyle(8);
00113 ha->Draw("p"); hb->Draw("psame");
00114 sprintf(name,"weight %.2f",weight); t->DrawTextNDC(0.15,0.75,name);
00115
00116 Print(&c1,psDir,"vertexZScaled");
00117
00118 }
00119
00120
00121
00122
00123
00124
00125
00126
00127
