Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-01-18 10:18:20

0001 
0002 void pfrich(const char *dfname, const char *cfname = 0)
0003 {
0004   auto fcfg  = new TFile(cfname ? cfname : dfname);
0005   auto geometry = dynamic_cast<CherenkovDetectorCollection*>(fcfg->Get("CherenkovDetectorCollection"));
0006   auto fdata = new TFile(dfname);
0007   TTree *t = dynamic_cast<TTree*>(fdata->Get("t")); 
0008   auto event = new CherenkovEvent();
0009   t->SetBranchAddress("e", &event);
0010 
0011   int nEvents = t->GetEntries();
0012 
0013   auto hxy = new TH2D("hxy", "", 650, -650., 650., 650, -650.0, 650.);
0014 
0015   for(unsigned ev=0; ev<nEvents; ev++) {
0016     t->GetEntry(ev);
0017 
0018     for(auto particle: event->ChargedParticles()) {
0019 
0020       for(auto rhistory: particle->GetRadiatorHistory()) {
0021     auto history  = particle->GetHistory (rhistory);
0022 
0023     for(auto photon: history->Photons()) {
0024       if (!photon->WasDetected() ) continue;
0025 
0026       TVector3 phx = photon->GetDetectionPosition();
0027       hxy->Fill(phx.X(), phx.Y());
0028     } //for photon
0029       } //for rhistory
0030     } //for particle
0031   } //for ev
0032 
0033   gStyle->SetOptStat(0);
0034   auto cv = new TCanvas("cv", "", 1000, 1000);
0035   hxy->GetXaxis()->SetTitle("Sensor plane X, [mm]");
0036   hxy->GetYaxis()->SetTitle("Sensor plane Y, [mm]");
0037   hxy->Draw("COL");
0038 } // pfrich()