#include "commonmacro/common.h" #include "common/Name.cc" #include "commonmacro/histutil.h" void makeUA1Ratio( const char* inName= "links/P01hi.minbias.2000.hist/hianalysis_1000.hist.root", const char* psDir="ps", int cut = 1, const char* outDir="./", const char* more = "west", float extraValue = 10) { cout << "--------------------------" << endl; cout << "in name=" << inName << endl << "ps dir=" << psDir << endl << "cut=" << cut << endl; cout << "--------------------------" << endl; inRoot = new TFile(inName); if(!inRoot){ cout << "cannot find the infile" << endl; return; } gStyle->SetOptStat(0);gStyle->SetPadTickX(1);gStyle->SetPadTickY(1); TCanvas c1("c1","c1",500,400); //------------------------------------------------------- const int nbin=2; TGraphAsymmErrors* gSpec[2]; //plus+minus TGraphAsymmErrors* gSpecPM[2][2]; // plus,minus TGraphAsymmErrors* gUA1, *gRatio; // get stuff for(int ibin=0;ibinGetName() << endl; if(!gSpec[ibin]) return; for(int ipm=0;ipm<2;ipm++){ setName(name,"gSpecCorrected",ibin,sPM[ipm].Data()); cout << name << endl; gSpecPM[ibin][ipm]=(TGraphAsymmErrors*)inRoot.Get(name); if(!gSpecPM[ibin][ipm]) return; } } // ratios for(int ibin=0;ibin<2; ibin++){ cout << "setting title"<< endl; sprintf(title,"%s / UA1 (bin %d)(cut %d)", gSpec[ibin]->GetName(),ibin,cut); Divide(&c1,1,1,title,inName); cout << "making" << endl; gUA1=makeUA1(gSpec[ibin]); gRatio=makeRatio(gSpec[ibin],gUA1); gRatio->Draw("ap"); sprintf(title,"%sOverUA1",gSpec[ibin]->GetName()); gRatio->SetMaximum(1.2); gRatio->SetMinimum(.2); Print(&c1,psDir,title); for(int ipm=0;ipm<2;ipm++){ // sName=gSpec[ibin][ipm]->GetName(); // charge sign // sName.Replace(0,sName.First("."),""); //sName.Replace(sName.First("."),sName.Length(),""); sprintf(title,"%s / UA1 (bin %d)(cut %d)", sName.Data(),gSpecPM[ibin][ipm]->GetName(),ibin,cut); Divide(&c1,1,1,title,inName); //gUA1=makeUA1(gSpec[ibin][ipm]); gRatio=makeRatio(gSpecPM[ibin][ipm],gUA1); gRatio->Draw("ap"); gRatio->SetMaximum(1.2); gRatio->SetMinimum(.2); sprintf(title,"%sOverUA1",gSpecPM[ibin][ipm]->GetName()); Print(&c1,psDir,title); } } } TGraphAsymmErrors* makeUA1(const TGraphAsymmErrors* graph) { // // find the bins // const Int_t nBin = graph->GetN(); Double_t* xValues = graph->GetX(); Double_t* yValues = graph->GetY(); Double_t* errXLow = graph->GetEXlow(); Double_t* errXHigh = graph->GetEXhigh(); // // ratios // /* 1.500000 0.6146703 6.9558643E-02 0.1131642 2.000000 0.5883511 8.5044047E-03 1.4454641E-02 2.500000 0.5589531 1.7606174E-03 3.1498480E-03 3.000000 0.5283710 4.8784356E-04 9.2329737E-04 3.500000 0.5029035 1.6314255E-04 3.2440131E-04 4.000000 0.4790372 6.2226289E-05 1.2989865E-04 4.500000 0.4574749 2.6327352E-05 5.7549289E-05 5.000000 0.4378719 1.2056932E-05 2.7535298E-05 5.500000 0.4186996 5.9010908E-06 1.4093855E-05 6.000000 0.4019553 3.0456495E-06 7.5770845E-06 */ sprintf(name,"h1tmp%s",graph->GetName()); TH1D* h1 = new TH1D(name,name,9,1.5,6.0); h1->SetBinContent(1,.61467); h1->SetBinContent(2,.58835); h1->SetBinContent(3,.55895); h1->SetBinContent(4,.52837); h1->SetBinContent(5,.50290); h1->SetBinContent(6,.47903); h1->SetBinContent(7,.45747); h1->SetBinContent(8,.43787); h1->SetBinContent(9,.40195); // fixed Double_t A = 266.6; Double_t p0 = 1.895; Double_t n = 12.98; // old /* double A=286; double p0=1.8; double n=12.14; */ Double_t TAA = 26; char func[200], mean[200]; sprintf(func,"2*pi*%.0f*%.0f*(1+x/%.2f)^(-%.2f)",A,TAA,p0,n); sprintf(mean,"2*pi*x*%.0f*%.0f*(1+x/%.2f)^(-%.2f)",A,TAA,p0,n); TF1* fUA1 = new TF1("fUA1",func,lowPt,highPt); TF1* fMean = new TF1("fUA1Mean",mean,lowPt,highPt); gUA1 = new TGraphAsymmErrors; for(Int_t i=0; iIntegral(lowEdge,highEdge); Double_t denom = fUA1->Integral(lowEdge,highEdge); // // just use star's x // Double_t UA1ErrLow = fabs(lowEdge-xValues[i]); Double_t UA1ErrHigh= fabs(highEdge-xValues[i]); //Double_t y = fUA1->Eval(xValues[i]); double y = denom/(highEdge-lowEdge); // ratio //Int_t iBin = h1->GetXaxis()->FindBin(xValues[i]); //Double_t fraction = h1->GetBinContent(iBin); // cout << "ratio : " << h1->GetXaxis()->GetBinLowEdge(iBin) << " " // << h1->GetXaxis()->GetBinUpEdge(iBin) << " = " << fraction << endl; //Double_t yCorrected = fraction*y; double yCorrected = y; cout << "UA1 uncorrected " << y << ", corrected y " << yCorrected << " graph y " << yValues[i] << endl; gUA1->SetPoint(i,xValues[i],yCorrected); gUA1->SetPointError(i,UA1ErrLow,UA1ErrHigh,0,0); } gUA1->SetName("gUA1Spectrum"); gUA1->SetTitle("gUA1Spectrum"); return gUA1; } TGraphAsymmErrors* makeRatio(TGraphAsymmErrors* gStar, TGraphAsymmErrors* gUA1) { cout << "********************************" << endl; cout << gStar->GetName() << endl; const Int_t nBin = gStar->GetN(); Double_t* starY = gStar->GetY(); Double_t* starX = gStar->GetX(); Double_t* starXErrLow = gStar->GetEXlow(); Double_t* starXErrHigh = gStar->GetEXhigh(); Double_t* starYErrHigh = gStar->GetEYhigh(); Double_t* ua1Y = gUA1->GetY(); Double_t* ua1YErrHigh = gUA1->GetEYhigh(); TGraphAsymmErrors* gRatio = new TGraphAsymmErrors; for(Int_t i=0; iSetPoint(i,starX[i],ratio); gRatio->SetPointError(i,starXErrLow[i],starXErrHigh[i],yError,yError); } TString sName = gStar->GetName(); sName.Append("UA1Ratio"); gRatio->SetName(sName.Data()); gRatio->SetTitle(sName.Data()); return gRatio; }