本文整理汇总了C++中TProfile::Write方法的典型用法代码示例。如果您正苦于以下问题：C++ TProfile::Write方法的具体用法？C++ TProfile::Write怎么用？C++ TProfile::Write使用的例子？那么恭喜您, 这里精选的方法代码示例或许可以为您提供帮助。您也可以进一步了解该方法所在类TProfile的用法示例。
在下文中一共展示了TProfile::Write方法的15个代码示例，这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞，您的评价将有助于系统推荐出更棒的C++代码示例。
示例1: storeSignalEfficiencies
void DileptonAnalyzer::storeSignalEfficiencies( TFile & outputFile ) {
// Store in a TProfile
TProfile *sigEffOne = new TProfile("SignalEfficOneDecay","Signal efficiencies when one exotic decays in this channel",10,0.5,10.5,-9.9e9,9.9e9);
TAxis *sigEffOneAxis = sigEffOne->GetXaxis();
TProfile *sigEffTwo = new TProfile("SignalEfficTwoDecay","Signal efficiencies when two exotic decays in this channel",10,0.5,10.5,-9.9e9,9.9e9);
TAxis *sigEffTwoAxis = sigEffTwo->GetXaxis();
// Fill set of 3 bins in each TProfile for each exotic
for ( int iExotic = 0; iExotic < nSignalParticles_; iExotic++ ) {
TString axisTitle;
sigEffOne->Fill( iExotic*nSignalParticles_+1, numEvents_oneSensitiveDecay_[iExotic] );
axisTitle.Form("Number of events coming into macro for exotic %d",signalPdgId_[iExotic]);
sigEffOneAxis->SetBinLabel( iExotic*nSignalParticles_+1, axisTitle.Data());
sigEffOne->Fill( iExotic*nSignalParticles_+2, numExoticsRECO_oneSensitiveDecay_[iExotic] );
axisTitle.Form("Number of candidates RECO'd for exotic %d",signalPdgId_[iExotic]);
sigEffOneAxis->SetBinLabel( iExotic*nSignalParticles_+2, axisTitle.Data());
sigEffOne->Fill( iExotic*nSignalParticles_+3, numExoticsCorrectRECO_oneSensitiveDecay_[iExotic] );
axisTitle.Form("Number of exotics correctly RECO'd for exotic %d",signalPdgId_[iExotic]);
sigEffOneAxis->SetBinLabel( iExotic*nSignalParticles_+3, axisTitle.Data());
sigEffTwo->Fill( iExotic*nSignalParticles_+1, numEvents_twoSensitiveDecay_[iExotic] );
axisTitle.Form("Number of events coming into macro for exotic %d",signalPdgId_[iExotic]);
sigEffTwoAxis->SetBinLabel( iExotic*nSignalParticles_+1, axisTitle.Data());
sigEffTwo->Fill( iExotic*nSignalParticles_+2, numExoticsRECO_twoSensitiveDecay_[iExotic] );
axisTitle.Form("Number of candidates RECO'd for exotic %d",signalPdgId_[iExotic]);
sigEffTwoAxis->SetBinLabel( iExotic*nSignalParticles_+2, axisTitle.Data());
sigEffTwo->Fill( iExotic*nSignalParticles_+3, numExoticsCorrectRECO_twoSensitiveDecay_[iExotic] );
axisTitle.Form("Number of exotics correctly RECO'd for exotic %d",signalPdgId_[iExotic]);
sigEffTwoAxis->SetBinLabel( iExotic*nSignalParticles_+3, axisTitle.Data());
}
outputFile.cd();
sigEffOne->Write();
sigEffTwo->Write();
}示例2: reportCutFlow
// Report cut flow findings to screen
void DileptonAnalyzer::reportCutFlow( TFile & outputFile ) {
// Store in a TProfile
TProfile *cutFlow = new TProfile("cutFlow","Cut flow for final set of analysis cuts",20,0.5,20.5,-9.9e9,9.9e9);
TAxis *cutFlowAxis = cutFlow->GetXaxis();
for ( unsigned int iCut = 0; iCut < cutNamesInOrder_.size(); iCut++ ) {
cutFlow->Fill( iCut+1, cutFlowMap_[cutNamesInOrder_[iCut]]);
cutFlowAxis->SetBinLabel( iCut+1, cutNamesInOrder_[iCut]);
}
outputFile.cd();
cutFlow->Write();
}示例3: doCoinc
//.........这里部分代码省略.........
hModulationAv->Fit(fmod);
printf("Estimates from time delay: Distance = %f +/- %f m -- Angle = %f +/- %f deg\n",fmod->GetParameter(1),fmod->GetParError(1),fmod->GetParameter(2),fmod->GetParError(2));
h->SetStats(0);
hDeltaThetaBack->Sumw2();
hDeltaPhiBack->Sumw2();
hThetaRelBack->Sumw2();
hDeltaThetaBack->Scale(0.1);
hDeltaPhiBack->Scale(0.1);
hThetaRelBack->Scale(0.1);
hAngleBack->Scale(0.1);
hAngle->Add(hAngleBack,-1);
printf("bin counting: SIGNAL = %f +/- %f\n",hDeltaPhi->Integral()-hDeltaPhiBack->Integral(),sqrt(hDeltaPhi->Integral()));
rate = (hDeltaPhi->Integral()-hDeltaPhiBack->Integral())/nsecGR*86400;
rateErr = sqrt(hDeltaPhi->Integral())/nsecGR*86400;
Float_t val,eval;
TCanvas *c1=new TCanvas();
TF1 *ff = new TF1("ff","[0]*[4]/[2]/sqrt(2*TMath::Pi())*TMath::Exp(-(x-[1])*(x-[1])*0.5/[2]/[2]) + [3]*[4]/6/[2]");
ff->SetParName(0,"signal");
ff->SetParName(1,"mean");
ff->SetParName(2,"sigma");
ff->SetParName(3,"background");
ff->SetParName(4,"bin width");
ff->SetParameter(0,42369);
ff->SetParameter(1,0);
ff->SetParLimits(2,10,maxwidth);
ff->SetParameter(2,350); // fix witdh if needed
ff->SetParameter(3,319);
ff->FixParameter(4,(tmax-tmin)/nbint); // bin width
ff->SetNpx(1000);
if(cout) cout->cd();
h->Fit(ff,"EI","",-10000,10000);
val = ff->GetParameter(2);
eval = ff->GetParError(2);
printf("significance = %f\n",ff->GetParameter(0)/sqrt(ff->GetParameter(0) + ff->GetParameter(3)));
h->Draw();
new TCanvas;
TF1 *func1 = (TF1 *) h->GetListOfFunctions()->At(0);
func1->SetLineColor(2);
h->SetLineColor(4);
TPaveText *text = new TPaveText(1500,(h->GetMinimum()+(h->GetMaximum()-h->GetMinimum())*0.6),9500,h->GetMaximum());
text->SetFillColor(0);
sprintf(title,"width = %5.1f #pm %5.1f",func1->GetParameter(2),func1->GetParError(2));
text->AddText(title);
sprintf(title,"signal (S) = %5.1f #pm %5.1f",func1->GetParameter(0),func1->GetParError(0));
text->AddText(title);
sprintf(title,"background (B) (3#sigma) = %5.1f #pm %5.1f",func1->GetParameter(3),func1->GetParError(3));
text->AddText(title);
sprintf(title,"significance (S/#sqrt{S+B}) = %5.1f",func1->GetParameter(0)/sqrt(func1->GetParameter(0)+func1->GetParameter(3)));
text->AddText(title);
text->SetFillStyle(0);
text->SetBorderSize(0);
text->Draw("SAME");
// correct nsecGR for the event rejected because of the number of satellites (event by event cut)
nsecGR *= neventsGRandSat/neventsGR;
printf("n_day = %f\nn_dayGR = %f\n",nsec*1./86400,nsecGR*1./86400);
text->AddText(Form("rate = %f #pm %f per day",func1->GetParameter(0)*86400/nsecGR,func1->GetParError(0)*86400/nsecGR));
TFile *fo = new TFile("outputCERN-01-02.root","RECREATE");
h->Write();
hDeltaTheta->Write();
hDeltaPhi->Write();
hThetaRel->Write();
hDeltaThetaBack->Write();
hDeltaPhiBack->Write();
hThetaRelBack->Write();
hAngle->Write();
hModulation->Write();
hModulation2->Write();
hModulationAv->Write();
hModulationAvCorr->Write();
hSinTheta->Write();
hSinTheta2->Write();
hnsigpeak->Write();
hRunCut[0]->Write();
hRunCut[1]->Write();
fo->Close();
return nsecGR*1./86400;
}示例4: main
//.........这里部分代码省略.........
51,0,50,0,2) ;
TH1F EEPCompareCoeffDistr_R1 ("EEPCompareCoeff_R1","EEPCompareCoeff_R1",1000,0,2) ;
TH1F EEPCompareCoeffDistr_R2 ("EEPCompareCoeff_R2","EEPCompareCoeff_R2",1000,0,2) ;
TH1F EEPCompareCoeffDistr_R3 ("EEPCompareCoeff_R3","EEPCompareCoeff_R3",1000,0,2) ;
TH1F EEPCompareCoeffDistr_R4 ("EEPCompareCoeff_R4","EEPCompareCoeff_R4",1000,0,2) ;
TH1F EEPCompareCoeffDistr_R5 ("EEPCompareCoeff_R5","EEPCompareCoeff_R5",1000,0,2) ;
// ECAL endcap +
for (int ix = 1 ; ix <= 100 ; ++ix)
for (int iy = 1 ; iy <= 100 ; ++iy)
{
int rad = static_cast<int> (sqrt ((ix - 50) * (ix - 50) +
(iy - 50) * (iy - 50))) ;
if (rad < EEradStart || rad > EEradEnd) continue ;
double phiTemp = atan2 (iy - 50, ix - 50) ;
if (phiTemp < 0) phiTemp += 2 * PI_GRECO ;
int phi = static_cast<int> ( phiTemp * 180 / PI_GRECO) ;
if (phi < EEphiStart || phi > EEphiEnd) continue ;
if (!EEDetId::validDetId (ix,iy,1)) continue ;
EEDetId det = EEDetId (ix, iy, 1, EEDetId::XYMODE) ;
double factor = *(iEEcalibMap.find (det.rawId ())) /
*(iEEscalibMap.find (det.rawId ())) ;
// std::cout << " iEEcalibMap rad = " << rad << " --> " << *(iEEcalibMap.find (det.rawId ()));
// std::cout << " iEEscalibMap --> " << *(iEEscalibMap.find (det.rawId ())) << std::endl;
EEPCompareCoeffDistr.Fill (factor) ;
EEPCompareCoeffMap.Fill (ix,iy,factor) ;
EEPCompareCoeffEtaTrend.Fill (rad,factor) ;
EEPCompareCoeffEtaProfile.Fill (rad,factor) ;
if (abs(rad) < 22) continue ;
else if (abs(rad) < 27) EEPCompareCoeffDistr_R1.Fill (factor) ;
else if (abs(rad) < 32) EEPCompareCoeffDistr_R2.Fill (factor) ;
else if (abs(rad) < 37) EEPCompareCoeffDistr_R3.Fill (factor) ;
else if (abs(rad) < 42) EEPCompareCoeffDistr_R4.Fill (factor) ;
else EEPCompareCoeffDistr_R5.Fill (factor) ;
} // ECAL endcap +
// ECAL endcap-
TH1F EEMCompareCoeffDistr ("EEMCompareCoeffDistr","EEMCompareCoeffDistr",200,0,2) ;
TH2F EEMCompareCoeffMap ("EEMCompareCoeffMap","EEMCompareCoeffMap",100,0,100,100,0,100) ;
TH2F EEMCompareCoeffEtaTrend ("EEMCompareCoeffEtaTrend",
"EEMCompareCoeffEtaTrend",
51,0,50,500,0,2) ;
TProfile EEMCompareCoeffEtaProfile ("EEMCompareCoeffEtaProfile",
"EEMCompareCoeffEtaProfile",
51,0,50,0,2) ;
// ECAL endcap -
for (int ix = 1 ; ix <= 100 ; ++ix)
for (int iy = 1 ; iy <= 100 ; ++iy)
{
int rad = static_cast<int> (sqrt ((ix - 50) * (ix - 50) +
(iy - 50) * (iy - 50))) ;
if (rad < EEradStart || rad > EEradEnd) continue ;
double phiTemp = atan2 (iy - 50, ix - 50) ;
if (phiTemp < 0) phiTemp += 2 * PI_GRECO ;
if (!EEDetId::validDetId (ix,iy,-1)) continue ;
EEDetId det = EEDetId (ix, iy, -1, EEDetId::XYMODE) ;
double factor = *(iEEcalibMap.find (det.rawId ())) *
*(iEEscalibMap.find (det.rawId ())) ;
EEMCompareCoeffDistr.Fill (factor) ;
EEMCompareCoeffMap.Fill (ix,iy,factor) ;
EEMCompareCoeffEtaTrend.Fill (rad,factor) ;
EEMCompareCoeffEtaProfile.Fill (rad,factor) ;
} // ECAL endcap -
TFile out (filename.c_str (),"recreate") ;
EEMCompareCoeffMap.Write() ;
EEMCompareCoeffDistr.Write() ;
EEMCompareCoeffEtaTrend.Write () ;
EEMCompareCoeffEtaProfile.Write () ;
EEPCompareCoeffMap.Write() ;
EEPCompareCoeffDistr.Write() ;
EEPCompareCoeffEtaTrend.Write () ;
EEPCompareCoeffEtaProfile.Write () ;
EEPCompareCoeffDistr_R1.Write () ;
EEPCompareCoeffDistr_R2.Write () ;
EEPCompareCoeffDistr_R3.Write () ;
EEPCompareCoeffDistr_R4.Write () ;
EEPCompareCoeffDistr_R5.Write () ;
out.Close() ;
//---- remove local database ----
gSystem->Exec("rm Uno.db");
gSystem->Exec("rm Due.db");
}示例5: if
//.........这里部分代码省略.........
const UShort_t nQ = maxH - 1;
correlations::ResultVector qs[nbin];
for(int ibin=0;ibin<nbin;ibin++)
qs[ibin] = correlations::ResultVector(nQ);
// --- Event loop --------------------------------------------------
Int_t nEvents = tree->GetEntries();
for (Int_t event = 0; event < nEvents; event++) {
tree->GetEntry(event);
int ntrk = M; int xbin=-1;
for(int j=0;j<nbin;j++)
if(ntrk<trkbin[j]&&ntrk>=trkbin[j+1])
xbin=j;
if(xbin<0 || xbin==nbin) continue;
tottrk[xbin]+=ntrk;
q[xbin].reset();
// printf("Event # %4u %4d particles ", event++, phis.GetSize());
for (UShort_t pa = 0; pa < M; pa++){
if(fabs(eta[pa])>etamax) continue;
if(pt[pa]<ptmin||pt[pa]>ptmax) continue; //event selection
// phis.Set(n,pPhis);
q[xbin].fill(phi[pa], 1.);
}
for (UShort_t i = 0; i < nQ; i++) {
UShort_t n = i + 2;
// printf("%s%d", i == 0 ? "" : "..", n);
timer.Reset();
timer.Start();
qs[xbin][i] += c[xbin]->calculate(n, h);
timer.Stop();
timing->Fill(n+.5, timer.RealTime());
}
// printf(" done\n");
Nevent[xbin]++;
}
file->Close();
for(int ibin=0;ibin<nbin;ibin++){
for (UShort_t i = 0; i < nQ; i++) {
// UShort_t iq = i+2;
// Double_t t = timing->GetBinContent(i+1);
// correlations::Complex rc = qs[i].eval();
// Printf("QC{%2d}: %12g + %12gi <t>: %10gs",
// iq, rc.real(), rc.imag(), t);
// if(i==0)Printf("v2{%2d}: %3g\n",2,sqrt(qs[0].eval().real()));
// if(i==2)Printf("v2{%2d}: %3g\n",4,TMath::Power(fabs(qs[2].eval().real()),1./4));
// if(i==4)Printf("v2{%2d}: %3g\n",6,TMath::Power(fabs(qs[4].eval().real()),1./6));
sumreals[ibin]->SetBinContent(i+1,qs[ibin][i]._sum.real());
sumimags[ibin]->SetBinContent(i+1,qs[ibin][i]._sum.imag());
weights[ibin]->SetBinContent(i+1,qs[ibin][i]._weights);
//reals->SetBinContent(i+1, rc.real());
//imags->SetBinContent(i+1, rc.imag());
}
}
/*
TCanvas* can = new TCanvas("C", "C");
can->SetTopMargin(0.15);
can->SetBottomMargin(0.15);
can->SetRightMargin(0.03);
can->Divide(1,3, 0, 0);
DrawInPad(can, 3, timing, true);
DrawInPad(can, 1, reals);
DrawInPad(can, 2, imags);
can->cd(0);
TLatex* ltx = new TLatex(0.5,0.995,c->name());
ltx->SetNDC(true);
ltx->SetTextAlign(23);
ltx->SetTextSize(0.04);
ltx->Draw();
can->Modified();
can->Update();
can->cd();
*/
TString out(mode);
out.ToLower();
file = TFile::Open(Form("%s/%s_%d.root",outdir.Data(),out.Data(),ifile), "RECREATE");
for(int ibin=0;ibin<nbin;ibin++){
sumimags[ibin]->Write();
sumreals[ibin]->Write();
weights[ibin]->Write();
}
Nevent.Write("Nevent");
tottrk.Write("tottrk");
timing->Write();
hs->Write();
file->Write();
file->Close();
for(int ibin=0;ibin<nbin;ibin++){
delete sumimags[ibin];
delete sumreals[ibin];
delete weights[ibin];
}
delete timing;
delete hs;
}示例6: Loop
//.........这里部分代码省略.........
Long64_t nbytes = 0, nb = 0;
for (Long64_t jentry=0; jentry<nentries;jentry++) {
Long64_t ientry = LoadTree(jentry);
if (ientry < 0) break;
nb = fChain->GetEntry(jentry); nbytes += nb;
// if (Cut(ientry) < 0) continue;
if(jentry % 1000 == 0) std::cout << "Processing entry " << jentry << std::endl;
// barrel
for(int h=0;h<std::min(nEBRecHits,10000);++h) {
double r=timeEBRecHits2[h]-timeEBRecHits[h];
if(energyEBRecHits[h]>10) {
timesCorr_EB->Fill(timeEBRecHits[h],timeEBRecHits2[h]);
timeRatioVsEta->Fill(etaEBRecHits[h],r);
timeRatioVsChi2_EB->Fill(chi2EBRecHits[h],r);
if(abs(timeEBRecHits[h])<t_window) time_0_EB->Fill(timeEBRecHits2[h]);
if(abs(timeEBRecHits[h]-25.)<t_window) time_p1_EB->Fill(timeEBRecHits2[h]);
if(abs(timeEBRecHits[h]-50.)<t_window) time_p2_EB->Fill(timeEBRecHits2[h]);
if(abs(timeEBRecHits[h]+25.)<t_window) time_m1_EB->Fill(timeEBRecHits2[h]);
if(abs(timeEBRecHits[h]+50.)<t_window) time_m2_EB->Fill(timeEBRecHits2[h]);
}
time1EnergyCorr_EB->Fill(timeEBRecHits[h],energyEBRecHits[h]);
time2EnergyCorr_EB->Fill(timeEBRecHits2[h],energyEBRecHits[h]);
timeRatioVsE_EB->Fill(energyEBRecHits[h],r);
}
// endcap
for(int h=0;h<std::min(nEERecHits,10000);++h) {
double r=timeEERecHits2[h]-timeEERecHits[h];
if(energyEERecHits[h]>10) {
timesCorr_EE->Fill(timeEERecHits[h],timeEERecHits2[h]);
timeRatioVsEta->Fill(etaEERecHits[h],r);
timeRatioVsChi2_EE->Fill(chi2EERecHits[h],r);
if(abs(timeEERecHits[h])<t_window) time_0_EE->Fill(timeEERecHits2[h]);
if(abs(timeEERecHits[h]-25.)<t_window) time_p1_EE->Fill(timeEERecHits2[h]);
if(abs(timeEERecHits[h]-50.)<t_window) time_p2_EE->Fill(timeEERecHits2[h]);
if(abs(timeEERecHits[h]+25.)<t_window) time_m1_EE->Fill(timeEERecHits2[h]);
if(abs(timeEERecHits[h]+50.)<t_window) time_m2_EE->Fill(timeEERecHits2[h]);
}
time1EnergyCorr_EE->Fill(timeEERecHits[h],energyEERecHits[h]);
time2EnergyCorr_EE->Fill(timeEERecHits2[h],energyEERecHits[h]);
timeRatioVsE_EE->Fill(energyEERecHits[h],r);
}
} // event loop
timeRatioVsEta->Write();
timesCorr_EB->Write();
time1EnergyCorr_EB->Write();
time2EnergyCorr_EB->Write();
timeRatioVsChi2_EB->Write();
timeRatioVsE_EB->Write();
timesCorr_EE->Write();
time1EnergyCorr_EE->Write();
time2EnergyCorr_EE->Write();
timeRatioVsChi2_EE->Write();
timeRatioVsE_EE->Write();
for(int i=0; i<(int) time1Dplots.size(); ++i) time1Dplots[i]->Write();
fileo->Close();
#endif
#if ANALYSIS == 2
ofstream txtfile;
txtfile.open(outputfilename,std::ios::trunc);
Long64_t nentries = fChain->GetEntries();
std::cout << "Total entries = " << nentries << std::endl;
Long64_t nbytes = 0, nb = 0;
for (Long64_t jentry=0; jentry<nentries;jentry++) {
Long64_t ientry = LoadTree(jentry);
if (ientry < 0) break;
nb = fChain->GetEntry(jentry); nbytes += nb;
if(jentry%1000==0) std::cout << "Processing entry " << jentry << "..." << std::endl;
for(int h=0; h<std::min(nEERecHits,10000); ++h) {
if(energyEERecHits[h]>5.0)
txtfile << eventNumber << "\t"
<< ixEERecHits[h] << "\t"
<< iyEERecHits[h] << "\t"
<< nTruePU[12] << "\t" // the bx = 0
<< energyEERecHits[h] << "\t"
<< timeEERecHits[h] << "\t"
<< chi2EERecHits[h] << "\t"
<< ootenergyEERecHits[h] << std::endl;
}
// if (Cut(ientry) < 0) continue;
}
txtfile.close();
#endif
}示例7: analisi_tree_1hit
//.........这里部分代码省略.........
TF1 *fstw = new TF1("fstw","gaus",-400.,400.);
h2t1_texp_deff_tw->FitSlicesY(fstw);
TH1D *h2t1_texp_deff_tw_1 = (TH1D *) gDirectory->FindObject("h2t1_texp_deff_tw_1");
TH1D *h2t1_texp_deff_tw_2 = (TH1D *) gDirectory->FindObject("h2t1_texp_deff_tw_2");
h2t1_texp_deff_tw_1->Draw("same");
h2t1_texp_deff_tw_2->Draw("same");
// TF1 *f1 = new TF1("f1","pol1",-1.25,0.5);
// //hprofx->GetYaxis()->SetRangeUser(-40.,100.);
// hprofx->Fit("f1","R");
// Double_t offset_p1,x1_p1;
// offset_p1= f1->GetParameter(0);
// x1_p1= f1->GetParameter(1);
// for(Int_t i=0;i<nentries;i++)
// {
// T->GetEntry(i);
// for(Int_t ip=0;ip<ncluster;ip++)
// tempo[ip] -= StartTime;
// if(ncluster == 1)
// {
// if(impulso_trasv>0.8 && impulso_trasv<1.2) // serve per gli exp time
// {
// //if( TMath::Abs(DeltaZ[0])<1.75)
// //{
// if(TMath::Abs(tempo[0]-exp_time_pi[0])<800./* per avere circa 3 sigmca che sia un pi*/ )
// {
//
// //Float_t posx = (DeltaX[0])* dch;
// Float_t tw1=tempo[0]- exp_time_pi[0];
//
// Float_t tw1corr=tw1-(offset_p1 + x1_p1 *DeltaX[0]);
//
// hprofxcorr->Fill(DeltaX[0],tw1corr,1);
// hprofxcorr->GetXaxis()->SetTitle("Dx1 (cm)");
// hprofxcorr->GetYaxis()->SetTitle("t1_corr=t1-t_exp_pi corr(ps)");
//
// htbest->Fill(tw1corr);
// htbest->GetXaxis()->SetTitle("t_best=t1_corr(ps)");
//
// //}
// }
// }
// }
// }
//
// TF1 *f1c = new TF1("f1c","pol1",-1.25,0.5);
// hprofxcorr->Fit("f1c","R");
TFile *fo2 = new TFile("output_ist_tree_1hit.root","RECREATE");
hch->Write();
ht1->Write();
h2resxz1->Write();
hresx1->Write();
hresz1->Write();
hresdist1->Write();
//hdeltax->Write();
//pxt->Write();
//pzt->Write();
//hdeltach->Write();
ht1_texp->Write();
ht1_texptot->Write();
//hprofx->Write();
//hprofz->Write();
//hprofxcorr->Write();
//htbest->Write();
hprofd->Write();
h2dxdzt1_texp->Write();
//h2dxdzt1_texp_dummy->Write();
h2t1_texp_deff->Write();
h2t1_texp_deff_1->Write();
//h2t1_texp_deff_2->Write();
hprofdeff->Write();
htcorrtw->Write();
h2t1_texp_deff_tw->Write();
h2t1_texp_deff_tw_1->Write();
h2t1_texp_deff_tw_2->Write();
h2t1_texp_TOT->Write();
hproft1_texp_TOT->Write();
h2t1_deff_TOT->Write();
h2t1_TOT_deff->Write();
hx->Write();
hz->Write();
fo2->Close();
if(!kCal){
printf("write calibration\n");
fcal = new TFile("calibration.root","RECREATE");
hCalX->Write();
hCalZ->Write();
fcal->Close();
}
system("say Ehi you, I have done");
}示例8: main
//.........这里部分代码省略.........
EEcalibMap.prefillMap () ;
MiscalibReaderFromXMLEcalEndcap calibEndcapreader (EEcalibMap) ;
if (!calibEndcapfile.empty ()) calibEndcapreader.parseXMLMiscalibFile (calibEndcapfile) ;
EcalIntercalibConstants* EECconstants =
new EcalIntercalibConstants (EEcalibMap.get ()) ;
EcalIntercalibConstants::EcalIntercalibConstantMap iEEcalibMap = EECconstants->getMap () ; //MF prende i vecchi coeff
//PG fill the histograms
//PG -------------------
TH1F EEPCompareCoeffDistr ("EEPCompareCoeffDistr","EEPCompareCoeffDistr",5000,0,2) ;
TH2F EEPCompareCoeffMap ("EEPCompareCoeffMap","EEPCompareCoeffMap",101,0,101,101,0,101) ;
TH2F EEPCompareCoeffEtaTrend ("EEPCompareCoeffEtaTrend",
"EEPCompareCoeffEtaTrend",
51,0,50,500,0,2) ;
TProfile EEPCompareCoeffEtaProfile ("EEPCompareCoeffEtaProfile",
"EEPCompareCoeffEtaProfile",
51,0,50,0,2) ;
TH1F EEPCompareCoeffDistr_R1 ("EEPCompareCoeff_R1","EEPCompareCoeff_R1",1000,0,2) ;
TH1F EEPCompareCoeffDistr_R2 ("EEPCompareCoeff_R2","EEPCompareCoeff_R2",1000,0,2) ;
TH1F EEPCompareCoeffDistr_R3 ("EEPCompareCoeff_R3","EEPCompareCoeff_R3",1000,0,2) ;
TH1F EEPCompareCoeffDistr_R4 ("EEPCompareCoeff_R4","EEPCompareCoeff_R4",1000,0,2) ;
TH1F EEPCompareCoeffDistr_R5 ("EEPCompareCoeff_R5","EEPCompareCoeff_R5",1000,0,2) ;
// ECAL endcap +
for (int ix = 1 ; ix <= 100 ; ++ix)
for (int iy = 1 ; iy <= 100 ; ++iy)
{
int rad = static_cast<int> (sqrt ((ix - 50) * (ix - 50) +
(iy - 50) * (iy - 50))) ;
if (rad < EEradStart || rad > EEradEnd) continue ;
double phiTemp = atan2 (iy - 50, ix - 50) ;
if (phiTemp < 0) phiTemp += 2 * PI_GRECO ;
int phi = static_cast<int> ( phiTemp * 180 / PI_GRECO) ;
if (phi < EEphiStart || phi > EEphiEnd) continue ;
if (!EEDetId::validDetId (ix,iy,1)) continue ;
EEDetId det = EEDetId (ix, iy, 1, EEDetId::XYMODE) ;
double factor = (iEEcalibMap.find (det.rawId ()))->second *
(iEEscalibMap.find (det.rawId ()))->second ;
EEPCompareCoeffDistr.Fill (factor) ;
EEPCompareCoeffMap.Fill (ix,iy,factor) ;
EEPCompareCoeffEtaTrend.Fill (rad,factor) ;
EEPCompareCoeffEtaProfile.Fill (rad,factor) ;
if (abs(rad) < 22) continue ;
else if (abs(rad) < 27) EEPCompareCoeffDistr_R1.Fill (factor) ;
else if (abs(rad) < 32) EEPCompareCoeffDistr_R2.Fill (factor) ;
else if (abs(rad) < 37) EEPCompareCoeffDistr_R3.Fill (factor) ;
else if (abs(rad) < 42) EEPCompareCoeffDistr_R4.Fill (factor) ;
else EEPCompareCoeffDistr_R5.Fill (factor) ;
} // ECAL endcap +
// ECAL endcap-
TH1F EEMCompareCoeffDistr ("EEMCompareCoeffDistr","EEMCompareCoeffDistr",200,0,2) ;
TH2F EEMCompareCoeffMap ("EEMCompareCoeffMap","EEMCompareCoeffMap",100,0,100,100,0,100) ;
TH2F EEMCompareCoeffEtaTrend ("EEMCompareCoeffEtaTrend",
"EEMCompareCoeffEtaTrend",
51,0,50,500,0,2) ;
TProfile EEMCompareCoeffEtaProfile ("EEMCompareCoeffEtaProfile",
"EEMCompareCoeffEtaProfile",
51,0,50,0,2) ;
// ECAL endcap -
for (int ix = 1 ; ix <= 100 ; ++ix)
for (int iy = 1 ; iy <= 100 ; ++iy)
{
int rad = static_cast<int> (sqrt ((ix - 50) * (ix - 50) +
(iy - 50) * (iy - 50))) ;
if (rad < EEradStart || rad > EEradEnd) continue ;
double phiTemp = atan2 (iy - 50, ix - 50) ;
if (phiTemp < 0) phiTemp += 2 * PI_GRECO ;
int phi = static_cast<int> ( phiTemp * 180 / PI_GRECO) ;
if (!EEDetId::validDetId (ix,iy,-1)) continue ;
EEDetId det = EEDetId (ix, iy, -1, EEDetId::XYMODE) ;
double factor = (iEEcalibMap.find (det.rawId ()))->second *
(iEEscalibMap.find (det.rawId ()))->second ;
EEMCompareCoeffDistr.Fill (factor) ;
EEMCompareCoeffMap.Fill (ix,iy,factor) ;
EEMCompareCoeffEtaTrend.Fill (rad,factor) ;
EEMCompareCoeffEtaProfile.Fill (rad,factor) ;
} // ECAL endcap -
TFile out (filename.c_str (),"recreate") ;
EEMCompareCoeffMap.Write() ;
EEMCompareCoeffDistr.Write() ;
EEMCompareCoeffEtaTrend.Write () ;
EEMCompareCoeffEtaProfile.Write () ;
EEPCompareCoeffMap.Write() ;
EEPCompareCoeffDistr.Write() ;
EEPCompareCoeffEtaTrend.Write () ;
EEPCompareCoeffEtaProfile.Write () ;
EEPCompareCoeffDistr_R1.Write () ;
EEPCompareCoeffDistr_R2.Write () ;
EEPCompareCoeffDistr_R3.Write () ;
EEPCompareCoeffDistr_R4.Write () ;
EEPCompareCoeffDistr_R5.Write () ;
out.Close() ;
}示例9: main
//.........这里部分代码省略.........
TVector3 vLje; vLje.SetPtEtaPhi(sortedJets[0].pt(), sortedJets[0].eta(), sortedJets[0].phi());
TVector3 vJet;
int kJrl = -1; double dJrl = -1.;
int kTrl = -1; double dTrl = -1.;
for (int j=0; j<sortedJets.size(); j++) {
double dJet = sortedJets[j].pt();
sortedJets[j].set_user_index(-1);
vJet.SetPtEtaPhi(dJet, sortedJets[j].eta(), sortedJets[j].phi());
if (TMath::Abs(vJet.DeltaPhi(vLje))>dCut) { sortedJets[j].set_user_index(1); if (dJet>dJrl) { dJrl = dJet; kJrl = j; } }
if (TMath::Abs(vJet.DeltaPhi(vLtk))>dCut) { sortedJets[j].set_user_index(2); if (dJet>dTrl) { dTrl = dJet; kTrl = j; } }
}
//=============================================================================
TVector3 v1sj, v2sj;
for (int j=0; j<sortedJets.size(); j++) {
Float_t dVar[kVar]; for (Int_t i=0; i<kVar; i++) dVar[i] = -1.;
dVar[kWgt] = 1.; dVar[kXsc] = 1.;
vJet.SetPtEtaPhi(sortedJets[j].pt(), sortedJets[j].eta(), sortedJets[j].phi());
//=============================================================================
dVar[kLje] = vLje.Pt(); if (sortedJets[j].user_index()==1) { dVar[kLjr] = ((kJrl==j) ? 1.5 : 0.5); }
dVar[kLtk] = vLtk.Pt(); if (sortedJets[j].user_index()==2) { dVar[kLtr] = ((kTrl==j) ? 1.5 : 0.5); }
//=============================================================================
dVar[kJet] = sortedJets[j].pt();
dVar[kAje] = sortedJets[j].area();
dVar[kMje] = sortedJets[j].m();
//=============================================================================
fastjet::Filter trimmer(subjDef, fastjet::SelectorPtFractionMin(0.));
fastjet::PseudoJet trimmdJet = trimmer(sortedJets[j]);
std::vector<fastjet::PseudoJet> trimmdSj = trimmdJet.pieces();
double d1sj = -1.; int k1sj = -1;
double d2sj = -1.; int k2sj = -1;
for (int i=0; i<trimmdSj.size(); i++) {
double dIsj = trimmdSj[i].pt(); if (dIsj<0.001) continue;
if (dIsj>d1sj) {
d2sj = d1sj; k2sj = k1sj;
d1sj = dIsj; k1sj = i;
} else if (dIsj>d2sj) {
d2sj = dIsj; k2sj = i;
}
}
//=============================================================================
if (d1sj>0.) {
v1sj.SetPtEtaPhi(d1sj, trimmdSj[k1sj].eta(), trimmdSj[k1sj].phi());
dVar[k1sz] = d1sj;
dVar[k1sA] = trimmdSj[k1sj].area();
dVar[k1sm] = trimmdSj[k1sj].m();
dVar[k1sr] = v1sj.DeltaR(vJet);
}
if (d2sj>0.) {
v2sj.SetPtEtaPhi(d2sj, trimmdSj[k2sj].eta(), trimmdSj[k2sj].phi());
dVar[k2sz] = d2sj;
dVar[k2sA] = trimmdSj[k2sj].area();
dVar[k2sm] = trimmdSj[k2sj].m();
dVar[k2sr] = v2sj.DeltaR(vJet);
}
if ((d1sj>0.) && (d2sj>0.)) dVar[kDsr] = v2sj.DeltaR(v1sj);
nt->Fill(dVar);
}
}
//=============================================================================
delete evt;
ascii_in >> evt;
}
//=============================================================================
file->cd(); nt->Write(); file->Close();
//=============================================================================
TString sXsec = sFile; sXsec.ReplaceAll("out", "xsecs");
file = TFile::Open(Form("%s/xsecs/%s.root",sPath.Data(),sXsec.Data()), "READ");
TH1D *hPtHat = (TH1D*)file->Get("hPtHat"); hPtHat->SetDirectory(0);
TH1D *hWeightSum = (TH1D*)file->Get("hWeightSum"); hWeightSum->SetDirectory(0);
TProfile *hSigmaGen = (TProfile*)file->Get("hSigmaGen"); hSigmaGen->SetDirectory(0);
file->Close();
//=============================================================================
sFile.ReplaceAll("out", "wgt");
file = TFile::Open(Form("%s.root",sFile.Data()), "NEW");
hPtHat->Write();
hWeightSum->Write();
hSigmaGen->Write();
file->Close();
//=============================================================================
cout << "DONE" << endl;
//=============================================================================
return 0;
}示例10: main
//.........这里部分代码省略.........
TFile *file = TFile::Open(Form("%s.root",sFile.Data()), "NEW");
TNtuple *nt = new TNtuple("nt", "", "fWgt:fJet:fAje:fMje:f1sj:f1sA:f1sr:f1sm:f2sj:f2sA:f2sr:f2sm:fDsr:fIsm");
//=============================================================================
HepMC::IO_GenEvent ascii_in(Form("%s/%s.hepmc",sPath.Data(),sFile.Data()), std::ios::in);
HepMC::GenEvent *evt = ascii_in.read_next_event();
while (evt) {
fjInput.resize(0);
TLorentzVector vPar;
for (HepMC::GenEvent::particle_const_iterator p=evt->particles_begin(); p!=evt->particles_end(); ++p) if ((*p)->status()==1) {
vPar.SetPtEtaPhiM((*p)->momentum().perp(), (*p)->momentum().eta(), (*p)->momentum().phi(), dMass);
if ((TMath::Abs(vPar.Eta())<dCutEtaMax)) {
fjInput.push_back(fastjet::PseudoJet(vPar.Px(), vPar.Py(), vPar.Pz(), vPar.E()));
}
}
//=============================================================================
fastjet::ClusterSequenceArea clustSeq(fjInput, jetsDef, areaDef);
std::vector<fastjet::PseudoJet> includJets = clustSeq.inclusive_jets(dJetsPtMin);
std::vector<fastjet::PseudoJet> selectJets = selectJet(includJets);
TLorentzVector vJet, v1sj, v2sj, vIsj;
for (int j=0; j<selectJets.size(); j++) {
double dJet = selectJets[j].pt();
vJet.SetPtEtaPhiM(dJet, selectJets[j].eta(), selectJets[j].phi(), selectJets[j].m());
fastjet::Filter trimmer(subjDef, fastjet::SelectorPtFractionMin(0.));
fastjet::PseudoJet trimmdJet = trimmer(selectJets[j]);
std::vector<fastjet::PseudoJet> trimmdSj = trimmdJet.pieces();
double d1sj = -1.; int k1sj = -1;
double d2sj = -1.; int k2sj = -1;
for (int i=0; i<trimmdSj.size(); i++) {
double dIsj = trimmdSj[i].pt(); if (dIsj<0.001) continue;
if (dIsj>d1sj) {
d2sj = d1sj; k2sj = k1sj;
d1sj = dIsj; k1sj = i;
} else if (dIsj>d2sj) {
d2sj = dIsj; k2sj = i;
}
}
//=============================================================================
Float_t dVar[] = { 1., dJet, selectJets[j].area(), selectJets[j].m(), d1sj, -1., -1., -1., d2sj, -1., -1., -1., -1., -1. };
if (d1sj>0.) {
v1sj.SetPtEtaPhiM(d1sj, trimmdSj[k1sj].eta(), trimmdSj[k1sj].phi(), trimmdSj[k1sj].m());
dVar[k1sr] = v1sj.DeltaR(vJet);
dVar[k1sm] = trimmdSj[k1sj].m();
dVar[k1sA] = trimmdSj[k1sj].area();
}
if (d2sj>0.) {
v2sj.SetPtEtaPhiM(d2sj, trimmdSj[k2sj].eta(), trimmdSj[k2sj].phi(), trimmdSj[k2sj].m());
dVar[k2sr] = v2sj.DeltaR(vJet);
dVar[k2sm] = trimmdSj[k2sj].m();
dVar[k2sA] = trimmdSj[k2sj].area();
}
if ((d1sj>0.) && (d2sj>0.)) {
vIsj = v1sj + v2sj;
dVar[kIsm] = vIsj.M();
dVar[kDsr] = v2sj.DeltaR(v1sj);
}
nt->Fill(dVar);
}
//=============================================================================
delete evt;
ascii_in >> evt;
}
//=============================================================================
file->cd(); nt->Write(); file->Close();
//=============================================================================
TString sXsec = sFile; sXsec.ReplaceAll("out", "xsecs");
file = TFile::Open(Form("%s/xsecs/%s.root",sPath.Data(),sXsec.Data()), "READ");
TH1D *hPtHat = (TH1D*)file->Get("hPtHat"); hPtHat->SetDirectory(0);
TH1D *hWeightSum = (TH1D*)file->Get("hWeightSum"); hWeightSum->SetDirectory(0);
TProfile *hSigmaGen = (TProfile*)file->Get("hSigmaGen"); hSigmaGen->SetDirectory(0);
file->Close();
//=============================================================================
sFile.ReplaceAll("out", "wgt");
file = TFile::Open(Form("%s.root",sFile.Data()), "NEW");
hPtHat->Write();
hWeightSum->Write();
hSigmaGen->Write();
file->Close();
//=============================================================================
cout << "DONE" << endl;
return 0;
}示例11: FitProfiles
//.........这里部分代码省略.........
{ std::cout<<" Wrong number of voltage steps. "<<std::endl; i++; continue;}
TProfile* Prof = iter->second[i];
if(!Prof)
{ std::cout<<" Profile "<<thestring.Data()<<"_"<<i<<" not found."<<std::endl; i++; continue;}
//if(Histo->GetEntries()) hNhits->Fill(Histo->Integral());
/*if(Histo->Integral()<20) //0.1
{ //std::cout<<" Not enought entries for histo "<<thestring.Data()<<std::endl;
i++; continue;}*/
detid = iter->first;
bool rmfit=false;
TF1* pol = new TF1("pol", "pol1", -3, 3);
pol->SetRange(-1,0);
Prof->Fit("pol", "qr");
double chi2 = pol->GetChisquare()/pol->GetNDF();
hChi2->Fill(chi2);
if(chi2>10) rmfit=true;
if( rmfit ||
// TIB modules
// TIB - 1.4.2.5
detid==369121605 || detid==369121606 || detid==369121614 ||
detid==369121613 || detid==369121610 || detid==369121609 ||
// TIB - 1.2.2.1
detid==369121390 || detid==369121382 || detid==369121386 ||
detid==369121385 || detid==369121389 || detid==369121381 ||
// others in TIB
detid==369121437 || detid==369142077 || detid==369121722 ||
detid==369125534 || detid==369137018 || detid==369121689 ||
detid==369121765 || detid==369137045 || detid==369169740 ||
detid==369121689 ||
// TOB modules
// TOB + 4.3.3.8
detid/10==436281512 || detid/10==436281528 || detid/10==436281508 ||
detid/10==436281524 || detid/10==436281520 || detid/10==436281516 ||
// others in TOB
detid/10==436228249 || detid/10==436232694 || detid/10==436228805 ||
detid/10==436244722 || detid/10==436245110 || detid/10==436249546 ||
detid/10==436310808 || detid/10==436312136 || detid/10==436315600 ||
// without 'sensors' option
detid==436281512 || detid==436281528 || detid==436281508 ||
detid==436281524 || detid==436281520 || detid==436281516 ||
detid==436228249 || detid==436232694 || detid==436228805 ||
detid==436244722 || detid==436245110 || detid==436249546 ||
detid==436310808 || detid==436312136 || detid==436315600 ||
// TID modules
detid==402664070 || detid==402664110 ||
// TEC modules in small scans
detid==470148196 || detid==470148200 || detid==470148204 ||
detid==470148228 || detid==470148232 || detid==470148236 ||
detid==470148240 || detid==470148261 || detid==470148262 ||
detid==470148265 || detid==470148266 || detid==470148292 ||
detid==470148296 || detid==470148300 || detid==470148304 ||
detid==470148324 || detid==470148328 || detid==470148332 ||
detid==470148336 || detid==470148340 ) {
Prof->Write();
std::cout << " Saving histo : " << thestring.Data() << std::endl;
}
if(rmfit) {nfitrm++; i++; continue;}
int subdet = ((detid>>25)&0x7);
int TECgeom=0;
if(subdet==6) TECgeom = ((detid>>5)&0x7);
// save values
detid = iter->first;
voltage = *itVolt;
index = i;
errvoltage = 2 ;
Slope = pol->GetParameter(1);
errSlope = pol->GetParError(1);
Origin = pol->GetParameter(0);
errOrigin = pol->GetParError(0);
Chi2 = chi2;
tree->Fill();
i++;
}
}
tree->Write();
//hNhits->Write();
//// If you want to store all the individual detId histograms uncomments this line !!!!
//myFile->Write();
myFile->Close();
}示例12: main
//.........这里部分代码省略.........
//=============================================================================
AliGenPythiaEventHeader *pHeadPy = (AliGenPythiaEventHeader*)pHeader->GenEventHeader();
if (!pHeadPy) continue;
hPtHat->Fill(pHeadPy->GetPtHard());
//=============================================================================
for (Int_t i=0; i<pStack->GetNtrack(); i++) if (pStack->IsPhysicalPrimary(i)) {
TParticle *pTrk = pStack->Particle(i); if (!pTrk) continue;
if (TMath::Abs(pTrk->Eta())>dCutEtaMax) { pTrk = 0; continue; }
// TParticlePDG *pPDG = pTrk->GetPDG(); if (!pPDG) { pTrk = 0; continue; }
fjInput.push_back(fastjet::PseudoJet(pTrk->Px(), pTrk->Py(), pTrk->Pz(), pTrk->P()));
// pPDG = 0;
pTrk = 0;
}
//=============================================================================
fastjet::ClusterSequenceArea clustSeq(fjInput, jetsDef, areaDef);
std::vector<fastjet::PseudoJet> includJets = clustSeq.inclusive_jets(dJetsPtMin);
// std::vector<fastjet::PseudoJet> subtedJets = bkgSubtractor(includJets);
std::vector<fastjet::PseudoJet> selectJets = selectJet(includJets);
// std::vector<fastjet::PseudoJet> sortedJets = fastjet::sorted_by_pt(selectJets);
for (int j=0; j<selectJets.size(); j++) {
double dJet = selectJets[j].pt();
hJet->Fill(dJet);
//=============================================================================
fastjet::Filter trimmer(subjDef, fastjet::SelectorPtFractionMin(0.));
fastjet::PseudoJet trimmdJet = trimmer(selectJets[j]);
std::vector<fastjet::PseudoJet> trimmdSj = trimmdJet.pieces();
double nIsj = 0.;
double d1sj = -1.; int k1sj = -1;
double d2sj = -1.; int k2sj = -1;
for (int i=0; i<trimmdSj.size(); i++) {
double dIsj = trimmdSj[i].pt(); if (dIsj<0.001) continue;
hJetIsj->Fill(dJet, dIsj);
hJetIsz->Fill(dJet, dIsj/dJet);
if (dIsj>d1sj) {
d2sj = d1sj; k2sj = k1sj;
d1sj = dIsj; k1sj = i;
} else if (dIsj>d2sj) {
d2sj = dIsj; k2sj = i;
} nIsj += 1.;
}
hJetNsj->Fill(dJet, nIsj);
if (d1sj>0.) { hJet1sj->Fill(dJet, d1sj); hJet1sz->Fill(dJet, d1sj/dJet); }
if (d2sj>0.) { hJet2sj->Fill(dJet, d2sj); hJet2sz->Fill(dJet, d2sj/dJet); }
if ((d1sj>0.) && (d2sj>0.)) {
double dsj = d1sj - d2sj;
double dsz = dsj / dJet;
hJetDsj->Fill(dJet, dsj);
hJetDsz->Fill(dJet, dsz);
}
}
//=============================================================================
pStack = 0;
pHeadPy = 0;
pHeader = 0;
}
//=============================================================================
rl->UnloadgAlice();
rl->UnloadHeader();
rl->UnloadKinematics();
rl->RemoveEventFolder();
//=============================================================================
TFile *file = TFile::Open(Form("%s/pyxsec_hists.root",sPath.Data()), "READ");
TList *lXsc = (TList*)file->Get("cFilterList");
file->Close();
TH1D *hWeightSum = (TH1D*)lXsc->FindObject("h1Trials"); hWeightSum->SetName("hWeightSum");
TProfile *hSigmaGen = (TProfile*)lXsc->FindObject("h1Xsec"); hSigmaGen->SetName("hSigmaGen");
//=============================================================================
file = TFile::Open(Form("%s.root",sFile.Data()), "NEW");
hPtHat->Write();
hWeightSum->Write();
hSigmaGen->Write();
list->Write();
file->Close();
//=============================================================================
cout << "DONE" << endl;
//=============================================================================
return 0;
}示例13: recoil_analysis
//.........这里部分代码省略.........
RooDataSet* RecoilRespData_phi3 = (RooDataSet*) RecoilRespData[i]->reduce("zphi>=3.&&zphi<4.");
RooDataSet* RecoilRespData_phi4 = (RooDataSet*) RecoilRespData[i]->reduce("zphi>=4.&&zphi<5.");
RooDataSet* RecoilRespData_phi5 = (RooDataSet*) RecoilRespData[i]->reduce("zphi>=5.");
RecoilRespData_phi0->printMultiline(cout, 1);
RecoilRespData_phi1->printMultiline(cout, 1);
RecoilRespData_phi2->printMultiline(cout, 1);
RecoilRespData_phi3->printMultiline(cout, 1);
RecoilRespData_phi4->printMultiline(cout, 1);
RecoilRespData_phi5->printMultiline(cout, 1);
cout<<"RESP bin start HIST "<<i<<endl;
genHist0 = (TH2D*) RecoilRespData_phi0->createHistogram(TString::Format("recoil_resp_hist_truephi0_bin_%d", i).Data(),
resp, Binning(250, 0., resp_bdr[i]),
YVar(dphi_norm, Binning(128, -1*dphi_bdr[i], dphi_bdr[i])));
genHist1 = (TH2D*) RecoilRespData_phi1->createHistogram(TString::Format("recoil_resp_hist_truephi1_bin_%d", i).Data(),
resp, Binning(250, 0., resp_bdr[i]),
YVar(dphi_norm, Binning(128, -1*dphi_bdr[i], dphi_bdr[i])));
genHist2 = (TH2D*) RecoilRespData_phi2->createHistogram(TString::Format("recoil_resp_hist_truephi2_bin_%d", i).Data(),
resp, Binning(250, 0., resp_bdr[i]),
YVar(dphi_norm, Binning(128, -1*dphi_bdr[i], dphi_bdr[i])));
genHist3 = (TH2D*) RecoilRespData_phi3->createHistogram(TString::Format("recoil_resp_hist_truephi3_bin_%d", i).Data(),
resp, Binning(250, 0., resp_bdr[i]),
YVar(dphi_norm, Binning(128, -1*dphi_bdr[i], dphi_bdr[i])));
genHist4 = (TH2D*) RecoilRespData_phi4->createHistogram(TString::Format("recoil_resp_hist_truephi4_bin_%d", i).Data(),
resp, Binning(250, 0., resp_bdr[i]),
YVar(dphi_norm, Binning(128, -1*dphi_bdr[i], dphi_bdr[i])));
genHist5 = (TH2D*) RecoilRespData_phi5->createHistogram(TString::Format("recoil_resp_hist_truephi5_bin_%d", i).Data(),
resp, Binning(250, 0., resp_bdr[i]),
YVar(dphi_norm, Binning(128, -1*dphi_bdr[i], dphi_bdr[i])));
cout<<"RESP bin start Write "<<i<<endl;
genHist0->Smooth(1);
genHist0->SetName(TString::Format("recoil_resp_hist_truephi0_bin_%d", i));
genHist0->SetTitle(TString::Format("recoil_resp_hist_truephi0_bin_%d", i));
double integral0 = genHist0->Integral();
genHist0->Scale(1/integral0);
genHist1->Smooth(1);
genHist1->SetName(TString::Format("recoil_resp_hist_truephi1_bin_%d", i));
genHist1->SetTitle(TString::Format("recoil_resp_hist_truephi1_bin_%d", i));
double integral1 = genHist1->Integral();
genHist1->Scale(1/integral1);
genHist2->Smooth(1);
genHist2->SetName(TString::Format("recoil_resp_hist_truephi2_bin_%d", i));
genHist2->SetTitle(TString::Format("recoil_resp_hist_truephi2_bin_%d", i));
double integral2 = genHist2->Integral();
genHist2->Scale(1/integral2);
genHist3->Smooth(1);
genHist3->SetName(TString::Format("recoil_resp_hist_truephi3_bin_%d", i));
genHist3->SetTitle(TString::Format("recoil_resp_hist_truephi3_bin_%d", i));
double integral3 = genHist3->Integral();
genHist3->Scale(1/integral3);
genHist4->Smooth(1);
genHist4->SetName(TString::Format("recoil_resp_hist_truephi4_bin_%d", i));
genHist4->SetTitle(TString::Format("recoil_resp_hist_truephi4_bin_%d", i));
double integral4 = genHist4->Integral();
genHist4->Scale(1/integral4);
genHist5->Smooth(1);示例14: extractFlowVZERO
//.........这里部分代码省略.........
hPidAll->Reset();
hPidAll->Add(hPidPrPiTPC);
hPidAll->Add(hPidPrKaTPC);
hPidAll->Add(hPidPrPrTPC);
hPidAll->Add(hPidPrElTPC);
hPidPrPiTPC->Divide(hPidAll);
hPidPrKaTPC->Divide(hPidAll);
hPidPrPrTPC->Divide(hPidAll);
hPidPrElTPC->Divide(hPidAll);
for(Int_t k=1;k <=hPi->GetNbinsX();k++){
Float_t pt = hPi->GetBinCenter(k);
Float_t xPi=hPi->GetBinContent(k)*hPidPiPi->Interpolate(pt) + hK->GetBinContent(k)*hPidPiKa->Interpolate(pt) + hPr->GetBinContent(k)*hPidPiPr->Interpolate(pt) + hEl->GetBinContent(k)*hPidPiEl->Interpolate(pt);
if(pt < 0.5)
xPi=hPi->GetBinContent(k)*hPidPiPiTPC->Interpolate(pt) + hK->GetBinContent(k)*hPidPiKaTPC->Interpolate(pt) + hPr->GetBinContent(k)*hPidPiPrTPC->Interpolate(pt) + hEl->GetBinContent(k)*hPidPiElTPC->Interpolate(pt);
Float_t xKa=hPi->GetBinContent(k)*hPidKaPi->Interpolate(pt) + hK->GetBinContent(k)*hPidKaKa->Interpolate(pt) + hPr->GetBinContent(k)*hPidKaPr->Interpolate(pt) + hEl->GetBinContent(k)*hPidKaEl->Interpolate(pt);
if(pt < 0.45)
xKa=hPi->GetBinContent(k)*hPidKaPiTPC->Interpolate(pt) + hK->GetBinContent(k)*hPidKaKaTPC->Interpolate(pt) + hPr->GetBinContent(k)*hPidKaPrTPC->Interpolate(pt) + hEl->GetBinContent(k)*hPidKaElTPC->Interpolate(pt);
Float_t xPr=hPi->GetBinContent(k)*hPidPrPi->Interpolate(pt) + hK->GetBinContent(k)*hPidPrKa->Interpolate(pt) + hPr->GetBinContent(k)*hPidPrPr->Interpolate(pt) + hEl->GetBinContent(k)*hPidPrEl->Interpolate(pt);
if(pt < 1)
xPr=hPi->GetBinContent(k)*hPidPrPiTPC->Interpolate(pt) + hK->GetBinContent(k)*hPidPrKaTPC->Interpolate(pt) + hPr->GetBinContent(k)*hPidPrPrTPC->Interpolate(pt) + hEl->GetBinContent(k)*hPidPrElTPC->Interpolate(pt);
hPi->SetBinContent(k,hPi->GetBinContent(k)*2 - xPi);
hK->SetBinContent(k,hK->GetBinContent(k)*2 - xKa);
hPr->SetBinContent(k,hPr->GetBinContent(k)*2 - xPr);
}
}
// antiproton Feed down
TProfile *pFromLambda = extractFlowVZEROsingle(icentr,11,arm,0,pTh,addbin,"pFromLambda",1,1);
TProfile *piFromK = extractFlowVZEROsingle(icentr,12,arm,0,pTh,addbin,"piFromK",1,1,1,1);
TProfile *pFromLambda2 = extractFlowVZEROsingle(icentr,11,arm,0,0.6,addbin,"pFromLambdanoPID",0,1);
TProfile *piFromK2 = extractFlowVZEROsingle(icentr,12,arm,0,0.6,addbin,"piFromKnoPID",0,1);
TProfile *piFromK3 = extractFlowVZEROsingle(icentr,12,arm,0,0.6,addbin,"piFromKnoPIDtof",1,1);
TH1D *hFeedSyst = NULL;
if(kFEEDcorr){
hFeedSyst = new TH1D(*hPr);
hFeedSyst->SetName("hFeedSyst");
hFeedSyst->Reset();
for(Int_t k=1;k <=hPr->GetNbinsX();k++){
Float_t contam = 3.23174e-01 * TMath::Exp(- 9.46743e-01 * hPr->GetBinCenter(k));
Float_t corr = contam * pFromLambda->GetBinContent(k)/(1-contam);
Float_t corrErr = contam * pFromLambda->GetBinError(k)/(1-contam);
Float_t value = hPr->GetBinContent(k)/(1-contam) - corr;
Float_t valueErr = hPr->GetBinError(k)/(1-contam);
hFeedSyst->SetBinContent(k,hPr->GetBinContent(k) - value);
hFeedSyst->SetBinContent(k,sqrt(corrErr*corrErr + valueErr*valueErr - hPr->GetBinError(k)*hPr->GetBinError(k)));
hPr->SetBinContent(k,value);
hPr->SetBinError(k,sqrt(corrErr*corrErr + valueErr*valueErr));
}
hFeedSyst->Divide(hPr);
}
// write output
snprintf(name,100,"results%03i-%03iv%i_pTh%3.1f%s.root",cMin[icentr],cMax[icentr+addbin],arm,pTh,type);
TFile *fout = new TFile(name,"RECREATE");
pAll->ProjectionX()->Write();
hPi->Write();
hK->Write();
hPr->Write();
if(isMC){
TH1D *pTmp = extractFlowVZEROsingle(icentr,0,arm,1,pTh,addbin,"allMC",1,1,-1,1)->ProjectionX();
pTmp->SetLineColor(6);
pTmp->SetMarkerColor(6);
pTmp->SetMarkerStyle(24);
pTmp->Write();
pTmp = extractFlowVZEROsingle(icentr,1,arm,1,pTh,addbin,"piMC",1,1,-1,1)->ProjectionX();
pTmp->SetLineColor(4);
pTmp->SetMarkerColor(4);
pTmp->SetMarkerStyle(24);
pTmp->Write();
pTmp = extractFlowVZEROsingle(icentr,2,arm,1,pTh,addbin,"kaMC",1,1,-1,1)->ProjectionX();
pTmp->SetLineColor(1);
pTmp->SetMarkerColor(1);
pTmp->SetMarkerStyle(25);
pTmp->Write();
pTmp = extractFlowVZEROsingle(icentr,3,arm,1,pTh,addbin,"prMC",1,1,-1,-1)->ProjectionX();
pTmp->SetLineColor(2);
pTmp->SetMarkerColor(2);
pTmp->SetMarkerStyle(26);
pTmp->Write();
}
extractFlowVZEROsingle(icentr,2,arm,0,pTh,addbin,"kProf")->Write();
extractFlowVZEROsingle(icentr,9,arm,0,pTh,addbin,"ks",0,1,1,1)->Write();
extractFlowVZEROsingle(icentr,9,arm,0,pTh,addbin,"ksMy",0,1,-1,-1)->Write();
extractFlowVZEROsingle(icentr,10,arm,0,pTh,addbin,"lambda")->Write();
pFromLambda->Write();
piFromK->Write();
pFromLambda2->Write();
piFromK2->Write();
piFromK3->Write();
if(hFeedSyst) hFeedSyst->Write();
fout->Close();
}示例15: main
//.........这里部分代码省略.........
continue;
CrystalCalibration cryCalib = *(eeCryCalibs[hashedIndex]);
int x = det.ix();
int y = det.iy();
//chiSquaredTotalHist->Fill(cryCalib.totalChi2);
//expectedStatPresHistEB->Fill(sqrt(1/expectedPresSumEB));
//expectedStatPresVsObservedMeanErrHistEB->Fill(sigmaM,sqrt(1/expectedPresSumEB));
//XXX: Filter events at default 0.5*meanE threshold
cryCalib.filterOutliers();
//numPointsErasedHist->Fill(numPointsErased);
//Write cryTimingHists
vector<TimingEvent> times = cryCalib.timingEvents;
for(vector<TimingEvent>::const_iterator timeItr = times.begin();
timeItr != times.end(); ++timeItr)
{
if(det.zside() < 0)
{
float weight = 1/((timeItr->sigmaTime)*(timeItr->sigmaTime));
cryTimingHistsEEM[x-1][y-1]->Fill(timeItr->time,weight);
}
else
{
float weight = 1/((timeItr->sigmaTime)*(timeItr->sigmaTime));
cryTimingHistsEEP[x-1][y-1]->Fill(timeItr->time,weight);
}
}
if(det.zside() < 0)
{
cryDirEEM->cd();
cryTimingHistsEEM[x-1][y-1]->Write();
}
else
{
cryDirEEP->cd();
cryTimingHistsEEP[x-1][y-1]->Write();
}
outfile->cd();
if(det.zside() < 0)
{
hitsPerCryHistEEM->SetBinContent(hashedIndex+1,cryCalib.timingEvents.size());
hitsPerCryMapEEM->Fill(x,y,cryCalib.timingEvents.size());
}
else
{
hitsPerCryHistEEP->SetBinContent(hashedIndex+1,cryCalib.timingEvents.size());
hitsPerCryMapEEP->Fill(x,y,cryCalib.timingEvents.size());
}
// Make timing calibs
double p1 = cryCalib.mean;
double p1err = cryCalib.meanE;
//cout << "cry ieta: " << ieta << " cry iphi: " << iphi << " p1: " << p1 << " p1err: " << p1err << endl;
if(cryCalib.timingEvents.size() < 10)
{
fileStreamProb << "Cry (only " << cryCalib.timingEvents.size() << " events) was calibrated to avg: " << det.zside() << ", "
<< x <<", " << y << ", hash: "
<< hashedIndex
<< "\t Calib: " << p1 << "\t Error: " << p1err << std::endl;
hashesToCalibrateToAvg.push_back(hashedIndex);
continue;
}