Fork me on GitHub

Changeset df408d2 in git for examples/Example6.C


Ignore:
Timestamp:
Mar 2, 2021, 3:14:09 PM (4 years ago)
Author:
Franco BEDESCHI <bed@…>
Branches:
master
Children:
92b8d11
Parents:
d3165fa
Message:

Fixed compatibility with ROOT5/minor changes to TrkUtil

File:
1 edited

Legend:

Unmodified
Added
Removed
  • examples/Example6.C

    rd3165fa rdf408d2  
    5151        // Track parameter pulls
    5252        TH1* histDpull = new TH1F("h_Dpull", "Pull impact parameter", 100, -10., 10.);
    53         TH1* histP0pull = new TH1F("h_P0pull", "Pull phi_{0}", 100, -10., 10.);
     53        TH1* histP0pull = new TH1F("h_P0pull", "Pull #phi_{0}", 100, -10., 10.);
    5454        TH1* histCpull = new TH1F("h_Cpull", "Pull half curvature", 100, -10., 10.);
    55         TH1* histZ0pull = new TH1F("h_Z0pull", "Pull Z0", 100, -10., 10.);
     55        TH1* histZ0pull = new TH1F("h_Z0pull", "Pull Z_{0}", 100, -10., 10.);
    5656        TH1* histCtpull = new TH1F("h_Ctpull", "Pull cotg(#theta)", 100, -10., 10.);
    5757        //
     
    6464        // CEntral track over min Pt
    6565        TH1* histPtgenC = new TH1F("h_PtgenC", "Generated Pt - Central", 500, 0., 50.);         // pt for central tracks;
    66         TH1* histPtobsC = new TH1F("h_PtobsC", "Reconstructed Pt - Central", 500, 0., 50.);     // pt for central
     66        TH1* histPtobsC = new TH1F("h_PtobsC", "Reconstructed Pt - Central", 500, 0., 50.);     // pt for central 
    6767        TH1* histCtgenH = new TH1F("h_CtgenH", "Generateded Cotangent", 100, -10., 10.);                // cot(theta) for high pt tracks;
    6868        TH1* histCtobsH = new TH1F("h_CtobsH", "Reconstructed Cotangent", 100, -10., 10.);      // cot(theta) for high pt tracks
     
    7171        TH1* histAccPt = new TH1F("h_AccPt", "Pt acceptance", 500, 0., 50.);
    7272        // Tracks over min pt
    73         TH1* histAccCtgH = new TH1F("h_AccCtgH", "Cotangent acceptance (high pt)", 500, -10., 10.);
     73        TH1* histAccCtgH = new TH1F("h_AccCtgH", "Cotangent acceptance (high pt)", 100, -10., 10.);
    7474        // tracks above min angle
    7575        TH1* histAccPtC = new TH1F("h_AccPtC", "Pt acceptance (central tracks)", 500, 0., 50.);
    7676
    7777        //
    78         // Get magnetic field
     78        // Get magnetifc field
    7979        Double_t Bz = 2.0;
    8080
     
    106106                                histZ0obs->Fill(obsZ0);
    107107                                histCtobs->Fill(obsCtg);
    108 
     108                                //
     109                                // Get associated generated particle
    109110                                GenParticle* gp = (GenParticle*)trk->Particle.GetObject();
    110 
     111                                //
    111112                                // Position of origin in meters
    112113                                Double_t x = 1.0e-3 * gp->X;
     
    139140                                Double_t pullZ0  = (obsZ0  - genZ0)  / trk->ErrorDZ;
    140141                                Double_t pullCtg = (obsCtg - genCtg) / trk->ErrorCtgTheta;
    141 
    142                                 //
    143                                 histDpull->Fill(pullD0);
     142                                //
     143                                histDpull ->Fill(pullD0);
    144144                                histP0pull->Fill(pullPhi);
    145                                 histCpull->Fill(pullC);
     145                                histCpull ->Fill(pullC);
    146146                                histZ0pull->Fill(pullZ0);
    147147                                histCtpull->Fill(pullCtg);
     
    149149                                // Acceptance plots
    150150                                Double_t obsPt = trk->PT;
     151                                histPtobs->Fill(obsPt);
    151152                                if (TMath::Abs(genCtg) < maxCtAcc)histPtobsC->Fill(obsPt); // pt for central tracks
    152                                 if (genPt > minPtAcc)histCtobsH->Fill(obsCtg);                          // cot(theta) for high pt tracks
     153                                if (obsPt > minPtAcc)histCtobsH->Fill(obsCtg);                          // cot(theta) for high pt tracks
    153154                        }               // End loop on tracks
    154155
     
    161162                                        GenParticle* gpart = (GenParticle*)branchGenPart->At(it);
    162163                                        //
    163                                         // Plot charged particle parameters
     164                                        // Plot charged particle parameters 
    164165                                        // Only final state particles (Status = 1)
    165166                                        if (gpart->Status == 1 && TMath::Abs(gpart->Charge) > 0) {
    166167                                                //
    167                                                 // Position of origin
    168                                                 Double_t x = gpart->X;
    169                                                 Double_t y = gpart->Y;
    170                                                 Double_t z = gpart->Z;
     168                                                // Position of origin in meters
     169                                                Double_t x = 1.0e-3 * gpart->X;
     170                                                Double_t y = 1.0e-3 * gpart->Y;
     171                                                Double_t z = 1.0e-3 * gpart->Z;
    171172                                                TVector3 xv(x, y, z);
    172173                                                //
     
    179180                                                // Original parameters
    180181                                                Double_t Q = gpart->Charge;
    181                                                 TVectorD genPar = TrkUtil::XPtoPar(xv, pv, Q, Bz);      // Get parameters from position and momentum
     182                                                TVectorD gPar = TrkUtil::XPtoPar(xv, pv, Q, Bz);        // Get parameters from position and momentum
     183                                                TVectorD genPar = TrkUtil::ParToMm(gPar);                       // Convert to mm
    182184                                                // Unpack parameters
    183                                                 Double_t genD0 = genPar(0);
     185                                                Double_t genD0  = genPar(0);
    184186                                                Double_t genPhi = genPar(1);
    185                                                 Double_t genC = genPar(2);
    186                                                 Double_t genZ0 = genPar(3);
     187                                                Double_t genC   = genPar(2);
     188                                                Double_t genZ0  = genPar(3);
    187189                                                Double_t genCtg = genPar(4);
     190                                                Double_t genPt  = pv.Pt();
    188191                                                // Fill histograms
    189                                                 histDgen->Fill(genD0);
     192                                                histDgen ->Fill(genD0);
    190193                                                histP0gen->Fill(genPhi);
    191                                                 histCgen->Fill(genC);
     194                                                histCgen ->Fill(genC);
    192195                                                histZ0gen->Fill(genZ0);
    193196                                                histCtgen->Fill(genCtg);
    194197                                                //
    195198                                                // Acceptance plots
    196                                                 if (TMath::Abs(genCtg) < maxCtAcc)histPtgenC->Fill(gpart->PT);
    197                                                 if (gpart->PT > minPtAcc)histCtgenH->Fill(genCtg);
     199                                                histPtgen->Fill(genPt);
     200                                                if (TMath::Abs(genCtg) < maxCtAcc)histPtgenC->Fill(genPt);
     201                                                if (genPt > minPtAcc)histCtgenH->Fill(genCtg);
    198202                                        }
    199203                                }
     
    222226        TCanvas* Cnv2 = new TCanvas("Cnv2", "Delphes observed track pulls", 150, 150, 900, 500);
    223227        Cnv2->Divide(3, 2);
    224         Cnv2->cd(1); histDpull->Draw();
    225         Cnv2->cd(2); histP0pull->Draw();
    226         Cnv2->cd(3); histCpull->Draw();
    227         Cnv2->cd(4); histZ0pull->Draw();
    228         Cnv2->cd(5); histCtpull->Draw();
     228        Cnv2->cd(1); gPad->SetLogy(1);
     229        histDpull->Draw();
     230        Cnv2->cd(2); gPad->SetLogy(1);
     231        histP0pull->Draw();
     232        Cnv2->cd(3); gPad->SetLogy(1);
     233        histCpull->Draw();
     234        Cnv2->cd(4); gPad->SetLogy(1);
     235        histZ0pull->Draw();
     236        Cnv2->cd(5); gPad->SetLogy(1);
     237        histCtpull->Draw();
    229238        //
    230239        // Acceptance plots
    231240        //
    232         TCanvas* Cnv4 = new TCanvas("Cnv4", "Delphes acceptance", 350, 350, 500, 500);
     241        TCanvas* Cnv4 = new TCanvas("Cnv4", "Delphes acceptance", 200, 200, 500, 500);
    233242        Cnv4->Divide(2, 2);
    234243        Cnv4->cd(1);
    235244        histCtgen->Draw();
     245        histCtobs->SetLineColor(kRed);
    236246        histCtobs->Draw("SAME");
    237247        Cnv4->cd(2);
     
    246256        }
    247257        histAccCtg->Draw();
    248         Cnv4->cd(3);
     258        Cnv4->cd(3); gPad->SetLogy(1);
    249259        histPtgen->Draw();
     260        histPtobs->SetLineColor(kRed);
    250261        histPtobs->Draw("SAME");
    251262        Cnv4->cd(4);
     
    261272        histAccPt->Draw();
    262273        //
    263         TCanvas* Cnv5 = new TCanvas("Cnv5", "Delphes acceptance (constrained)", 400, 400, 500, 500);
     274        TCanvas* Cnv5 = new TCanvas("Cnv5", "Delphes acceptance (constrained)", 250, 250, 500, 500);
    264275        Cnv5->Divide(2, 2);
    265276        Cnv5->cd(1);
    266277        histCtgenH->Draw();
     278        histCtobsH->SetLineColor(kRed);
    267279        histCtobsH->Draw("SAME");
    268280        Cnv5->cd(2);
     
    277289        }
    278290        histAccCtgH->Draw();
    279         Cnv5->cd(3);
    280         histPtgenC->Draw();
     291        Cnv5->cd(3); gPad->SetLogy(1);
     292        histPtgenC->Draw();
     293        histPtobsC->SetLineColor(kRed);
    281294        histPtobsC->Draw("SAME");
    282295        Cnv5->cd(4);
Note: See TracChangeset for help on using the changeset viewer.