|
| 1 | +R__LOAD_LIBRARY(libDelphes) |
| 2 | +R__LOAD_LIBRARY(libDelphesO2) |
| 3 | + |
| 4 | +double tof_radius = 100.; // [cm] |
| 5 | +double tof_length = 200.; // [cm] |
| 6 | +double tof_sigmat = 0.02; // [ns] |
| 7 | + |
| 8 | +void |
| 9 | +ftof(const char *inputFile = "delphes.root", |
| 10 | + const char *outputFile = "ftof.root") |
| 11 | +{ |
| 12 | + |
| 13 | + // Create chain of root trees |
| 14 | + TChain chain("Delphes"); |
| 15 | + chain.Add(inputFile); |
| 16 | + |
| 17 | + // Create object of class ExRootTreeReader |
| 18 | + auto treeReader = new ExRootTreeReader(&chain); |
| 19 | + auto numberOfEntries = treeReader->GetEntries(); |
| 20 | + |
| 21 | + // Get pointers to branches used in this analysis |
| 22 | + auto events = treeReader->UseBranch("Event"); |
| 23 | + auto tracks = treeReader->UseBranch("Track"); |
| 24 | + auto particles = treeReader->UseBranch("Particle"); |
| 25 | + |
| 26 | + // TOF layer |
| 27 | + o2::delphes::TOFLayer toflayer; |
| 28 | + toflayer.setup(tof_radius, tof_length, tof_sigmat); |
| 29 | + toflayer.setType(o2::delphes::TOFLayer::kForward); |
| 30 | + toflayer.setRadiusIn(10.); |
| 31 | + |
| 32 | + // smearer |
| 33 | + o2::delphes::TrackSmearer smearer; |
| 34 | + smearer.useEfficiency(true); |
| 35 | + smearer.loadTable(11, "lutCovm.el.dat"); |
| 36 | + smearer.loadTable(13, "lutCovm.mu.dat"); |
| 37 | + smearer.loadTable(211, "lutCovm.pi.dat"); |
| 38 | + smearer.loadTable(321, "lutCovm.ka.dat"); |
| 39 | + smearer.loadTable(2212, "lutCovm.pr.dat"); |
| 40 | + |
| 41 | + // logx binning |
| 42 | + const Int_t nbins = 80; |
| 43 | + double xmin = 1.e-2; |
| 44 | + double xmax = 1.e2; |
| 45 | + double logxmin = std::log10(xmin); |
| 46 | + double logxmax = std::log10(xmax); |
| 47 | + double binwidth = (logxmax - logxmin) / nbins; |
| 48 | + double xbins[nbins + 1]; |
| 49 | + xbins[0] = xmin; |
| 50 | + for (Int_t i = 1; i <= nbins; ++i) |
| 51 | + xbins[i] = xmin + std::pow(10., logxmin + i * binwidth); |
| 52 | + |
| 53 | + // histograms |
| 54 | + auto hTime0 = new TH1F("hTime0", ";t_{0} (ns)", 1000, -1., 1.); |
| 55 | + auto hBetaP = new TH2F("hBetaP", ";#it{p} (GeV/#it{c});#beta", nbins, xbins, 1000, 0.1, 1.1); |
| 56 | + TH2 *hNsigmaP[5]; |
| 57 | + const char *pname[5] = {"el", "mu", "pi", "ka", "pr"}; |
| 58 | + const char *plabel[5] = {"e", "#mu", "#pi", "K", "p"}; |
| 59 | + for (int i = 0; i < 5; ++i) |
| 60 | + hNsigmaP[i] = new TH2F(Form("hNsigmaP_%s", pname[i]), Form("#it{p} (GeV/#it{c});n#sigma_{%s}", plabel[i]), nbins, xbins, 200, -10., 10.); |
| 61 | + |
| 62 | + for (Int_t ientry = 0; ientry < numberOfEntries; ++ientry) { |
| 63 | + |
| 64 | + // Load selected branches with data from specified event |
| 65 | + treeReader->ReadEntry(ientry); |
| 66 | + |
| 67 | + // loop over tracks, smear and store TOF tracks |
| 68 | + std::vector<Track *> tof_tracks; |
| 69 | + for (Int_t itrack = 0; itrack < tracks->GetEntries(); ++itrack) { |
| 70 | + |
| 71 | + // get track and corresponding particle |
| 72 | + auto track = (Track *)tracks->At(itrack); |
| 73 | + auto particle = (GenParticle *)track->Particle.GetObject(); |
| 74 | + |
| 75 | + // smear track |
| 76 | + if (!smearer.smearTrack(*track)) continue; |
| 77 | + |
| 78 | + // select primaries based on 3 sigma DCA cuts |
| 79 | + if (fabs(track->D0 / track->ErrorD0) > 3.) continue; |
| 80 | + if (fabs(track->DZ / track->ErrorDZ) > 3.) continue; |
| 81 | + |
| 82 | + // check if has TOF |
| 83 | + if (!toflayer.hasTOF(*track)) continue; |
| 84 | + |
| 85 | + // push track |
| 86 | + tof_tracks.push_back(track); |
| 87 | + |
| 88 | + } |
| 89 | + |
| 90 | + // compute the event time |
| 91 | + std::array<float, 2> tzero; |
| 92 | + toflayer.eventTime(tof_tracks, tzero); |
| 93 | + hTime0->Fill(tzero[0]); |
| 94 | + |
| 95 | + // loop over tracks and do PID |
| 96 | + for (auto track : tof_tracks) { |
| 97 | + |
| 98 | + // fill beta-p |
| 99 | + auto p = track->P; |
| 100 | + auto beta = toflayer.getBeta(*track); |
| 101 | + hBetaP->Fill(p, beta); |
| 102 | + |
| 103 | + // fill nsigma |
| 104 | + std::array<float, 5> deltat, nsigma; |
| 105 | + toflayer.makePID(*track, deltat, nsigma); |
| 106 | + for (int i = 0; i < 5; ++i) |
| 107 | + hNsigmaP[i]->Fill(p, nsigma[i]); |
| 108 | + |
| 109 | + } |
| 110 | + } |
| 111 | + |
| 112 | + auto fout = TFile::Open(outputFile, "RECREATE"); |
| 113 | + hTime0->Write(); |
| 114 | + hBetaP->Write(); |
| 115 | + for (int i = 0; i < 5; ++i) |
| 116 | + hNsigmaP[i]->Write(); |
| 117 | + fout->Close(); |
| 118 | + |
| 119 | +} |
0 commit comments