--- /dev/null
+//
+// This macro reads comparison files and makes plots.
+//
+
+void DrawComparison
+(const char *fileName = "match-comparison.root")
+{
+ TH1D *hg = new TH1D("hg", "Good matches / true matches", 20, 0.0, 10.0);
+ TH1D *hf = new TH1D("hf", "Fake matches / found matches", 20, 0.0, 10.0);
+
+ //
+ // Open file
+ //
+ TFile *file = TFile::Open(fileName);
+ if (!file) return;
+
+ TH1D *hgood = (TH1D*)file->Get("hgood");
+ TH1D *hfake = (TH1D*)file->Get("hfake");
+ TH1D *htrue = (TH1D*)file->Get("htrue");
+ TH1D *hfound = (TH1D*)file->Get("hfound");
+
+ gROOT->SetStyle("Plain");
+ gStyle->SetOptStat(0);
+
+ hg->Divide(hgood, htrue, 100.0, 1.0, "b");
+ hf->Divide(hfake, hfound, 100.0, 1.0, "b");
+
+ TCanvas *c = new TCanvas("c", "", 0, 0, 800, 600);
+
+ hg->SetMarkerStyle(21);
+ hf->SetMarkerStyle(25);
+
+ hg->SetXTitle("p_{T} (GeV/c)");
+
+ hg->GetXaxis()->SetRangeUser(0.0, 6.0);
+ hf->GetXaxis()->SetRangeUser(0.0, 6.0);
+
+ hg->SetMaximum(120.0);
+ hg->SetMinimum(0.0);
+ hg->Draw("PE1");
+ hf->Draw("PE1same");
+}
--- /dev/null
+//
+// This macro performs the matching operation between ESD tracks
+// and EMCAL clusters (stored as AliESDCaloCluster objects).
+// Due to the necessity to know the magnetic field, it is supposed that
+// this macro runs in the directory where event files are stored
+// (galice.root --> recovery magnetic field), and saves its output in
+// a directory specified in the argument (default = work directory)
+//
+
+void FindMatches(const char *fileOut = "matchESD.root")
+{
+ //
+ // Initialize AliRun manager.
+ //
+ if (gAlice) {
+ delete gAlice;
+ gAlice = 0;
+ }
+
+ //
+ // Initialize AliRunLoader, in order
+ // to read magnetic field from simulation.
+ //
+ AliRunLoader *rl = AliRunLoader::Open("galice.root");
+ if (!rl) return;
+ rl->LoadgAlice();
+ gAlice = rl->GetAliRun();
+ AliMagF *magf = gAlice->Field();
+ Bool_t constField = (magf->Type() == 1);
+ AliTracker::SetFieldMap(magf, constField);
+
+ //
+ // Open ESD file and recoveries TTree of ESD objects.
+ //
+ TFile *esdFile = new TFile("AliESDs.root");
+ TTree *esdTree = (TTree*)esdFile->Get("esdTree");
+ AliESD *esd = 0;
+ esdTree->SetBranchAddress("ESD", &esd);
+ Long64_t nEvents = esdTree->GetEntries();
+
+ //
+ // Set this important flag to the tracker.
+ // This example works with and AliESD file already saved
+ // after the tracking procedure done during AliReconstruction::Run().
+ // All tracks saved in such a file, are already propagated to the vertex,
+ // while EMCAL matching needs to be done using track status in the outermost
+ // available position.
+ // Then, we must use the "outer" parameters in the ESD track, and we must
+ // tell to the track to copy the outer parameters from its AliESDtrack seed.
+ //
+ AliEMCALTrack::SetUseOuterParams(kTRUE);
+
+ //
+ // Instantiate match finder, and set some cuts.
+ // Distances are in centimeters, angles in degrees.
+ //
+ AliEMCALTracker *mf = new AliEMCALTracker;
+ mf->SetCutX(5.0);
+ mf->SetCutY(5.0);
+ mf->SetCutZ(5.0);
+ mf->SetCutAlpha(-50., 50.); // --> i.e. exclude this cut
+ mf->SetCutAngle(10.);
+ mf->SetMaxDistance(5.0);
+
+ //
+ // Define how to manage energy loss correction.
+ // Actually, these corrections are exlcluded.
+ //
+ // mf->SetCorrection(0.0604557, 34.5437); // at the moment, we exclude energy loss correction
+ mf->SetTrackCorrectionMode("NONE");
+ mf->SetNumberOfSteps(0);
+ //TGeoManager::Import("misaligned_geometry.root");
+
+ //
+ // Before starting looping on files, the output file is created.
+ // It will be structurally identical to the ESD source tree, with the
+ // only difference that here the "fEMCALindex" datamember of the ESD tracks
+ // will have been set to a meaningful value.
+ //
+ TFile *outFile = TFile::Open(fileOut, "RECREATE");
+ TTree *outTree = new TTree("esdTree", "ESD with matched clusters");
+ outTree->Branch("ESD", "AliESD", &esd);
+
+ //
+ // Loop on events.
+ // The whole track matching procedure is compiled in the
+ // method "PropagateBack" which accepts and input ESD collection.
+ // This method modifies the passed object then, the same object is linked
+ // to source tree and target tree branch.
+ //
+ Int_t nTracks, nStep1, nStep2, nSaved;
+ for (Long64_t iev = 0; iev < nEvents; iev++) {
+ cout << "Finding matches in event " << iev + 1 << "/" << nEvents << endl;
+ esdTree->GetEntry(iev);
+ mf->PropagateBack(esd);
+ outTree->Fill();
+ }
+
+ //
+ // Save processed tree and close output file.
+ //
+ outFile->cd();
+ outTree->Write("esdTree", TObject::kOverwrite);
+ outFile->Close();
+}
--- /dev/null
+//
+// This macros compares the matches found by algorithm with the true matches
+// found with the "SaveTrueMatchesSimple.C" macro.
+// It saves 4 histogram, which contain the Pt distribution of:
+// - all found matches
+// - all correctly found matches
+// - all wrong (fake) found matches
+// - all true matches
+// which will then be availabel for computing efficiency and contamination.
+//
+
+class match_t
+{
+ public:
+
+ Int_t label; // GEANT label of particle
+ Int_t indexT; // index of track in ESD collection
+ Int_t indexC; // index of cluster in ESD collection
+ Double_t p[3]; // track momentum
+ Double_t v[3]; // track vertex
+};
+
+void MatchComparison()
+{
+ //
+ // Initialize AliRun manager
+ //
+ if (gAlice) {
+ delete gAlice;
+ gAlice = 0;
+ }
+
+ //
+ // Initialize run loader and load Kinematics
+ //
+ AliRunLoader *runLoader = AliRunLoader::Open("galice.root");
+ if (!runLoader) return;
+ runLoader->LoadgAlice();
+ gAlice = runLoader->GetAliRun();
+ runLoader->LoadKinematics();
+
+ //
+ // Initialize histograms with their error computation
+ //
+ TH1D *hgood = new TH1D("hgood", "Well matched tracks", 20, 0.0, 10.0);
+ TH1D *hfake = new TH1D("hfake", "Fake matched tracks", 20, 0.0, 10.0);
+ TH1D *htrue = new TH1D("htrue", "True matches", 20, 0.0, 10.0);
+ TH1D *hfound = new TH1D("hfound", "Found matches", 20, 0.0, 10.0);
+ hgood->Sumw2();
+ hfake->Sumw2();
+ htrue->Sumw2();
+ hfound->Sumw2();
+
+ //
+ // Open file containing true matches,
+ // retrieve the Tree and link to a cursor.
+ //
+ TFile *fileTrue = TFile::Open("true-matches.root");
+ match_t trueMatch;
+
+ //
+ // Open file of found matches,
+ // link the modified ESD container.
+ //
+ TFile *fileFound = TFile::Open("matchESD.root");
+ TTree *treeFound = (TTree*)fileFound->Get("esdTree");
+ AliESD *esd = 0;
+ treeFound->SetBranchAddress("ESD", &esd);
+ Long64_t nEvents = treeFound->GetEntries();
+
+ //
+ // Loop on all events
+ //
+ Int_t im, it, ic, nTrueMatches, nTracks;
+ Int_t label, trkLabel, cluLabel;
+ for (Long64_t iev = 0; iev < nEvents; iev++) {
+
+ // get true matches tree of given event
+ TTree *treeTrue = (TTree*)fileTrue->Get(Form("tm_%d", iev));
+ treeTrue->SetBranchAddress("matches", &trueMatch);
+ nTrueMatches = treeTrue->GetEntries();
+
+ // set TTree pointers to selected event
+ runLoader->GetEvent(iev);
+ treeFound->GetEntry(iev);
+ AliStack *stack = runLoader->Stack();
+ nTracks = esd->GetNumberOfTracks();
+
+ // read all true pairs
+ for (im = 0; im < nTrueMatches; im++) {
+ treeTrue->GetEntry(im);
+ AliESDtrack *track = esd->GetTrack(trueMatch.indexT);
+ label = TMath::Abs(track->GetLabel());
+ TParticle *p = stack->Particle(label);
+ htrue->Fill(p->Pt());
+ }
+
+ // compare found matches
+ for (Int_t it = 0; it < nTracks; it++) {
+ AliESDtrack *track = esd->GetTrack(it);
+ ic = track->GetEMCALcluster();
+ if (ic == -999999999) continue;
+ ic = TMath::Abs(ic);
+ AliESDCaloCluster *cl = esd->GetCaloCluster(ic);
+ if (!cl) continue;
+ trkLabel = track->GetLabel();
+ cluLabel = cl->GetPrimaryIndex();
+ if (trkLabel == cluLabel && trkLabel > 0) {
+ TParticle *p = stack->Particle(TMath::Abs(trkLabel));
+ hgood->Fill(p->Pt());
+ hfound->Fill(p->Pt());
+ }
+ else {
+ TParticle *p = stack->Particle(TMath::Abs(trkLabel));
+ hfake->Fill(p->Pt());
+ hfound->Fill(p->Pt());
+ }
+ }
+ }
+
+ cout << "True matches : " << htrue->GetEntries() << endl;
+ cout << "Found matches: " << hfound->GetEntries() << endl;
+ cout << "Good matches : " << hgood->GetEntries() << endl;
+ cout << "Fake matches : " << hfake->GetEntries() << endl;
+
+ TFile *fout = TFile::Open("match-comparison.root", "RECREATE");
+ hgood->Write();
+ hfake->Write();
+ htrue->Write();
+ hfound->Write();
+ fout->Close();
+}
--- /dev/null
+//
+// This macro generates a simple TTree containing
+// all true matches collected from one event.
+//
+// For each match it is stored:
+// - label ID of the track
+// - label ID of the EMCAL cluster
+// - momentum (3D) of the track
+// - vertex (3D) of the track
+//
+// All tracks which come from kinks are excluded
+// from this computation, and fake tracks are collected.
+// If desired, it is possible to change these settings
+// operating on the Bool_t variables listed below.
+//
+
+Bool_t rejectKinks = kTRUE; // switch to TRUE to include kinks in computation
+Bool_t rejectFakes = kFALSE; // switch to TRUE to include fake tracks in computation
+
+class match_t
+{
+ public:
+
+ Int_t label; // GEANT label of particle
+ Int_t indexT; // index of track in ESD collection
+ Int_t indexC; // index of cluster in ESD collection
+ Double_t p[3]; // track momentum
+ Double_t v[3]; // track vertex
+};
+
+//
+// Read AliESDs.root file and saves all true pairs of trak-cluster.
+//
+void SaveTrueMatchesSimple_compiled(const char *outFileName)
+{
+ //
+ // open ESD file, retrieve tree and link to branch cursor
+ //
+ TFile *srcFile = TFile::Open("AliESDs.root");
+ if (!srcFile) return;
+ TTree *srcTree = (TTree*)srcFile->Get("esdTree");
+ if (!srcTree) return;
+ AliESD *esd = 0;
+ srcTree->SetBranchAddress("ESD", &esd);
+ Long64_t nEvents = srcTree->GetEntries();
+
+ //
+ // Open output file and create output tree
+ //
+ TFile *outFile = new TFile(outFileName, "RECREATE");
+ match_t match;
+
+ //
+ // Loop on events and store true matches
+ //
+ Bool_t isKink;
+ Int_t label, count, nTracks, firstCluster, lastCluster;
+ for (Long64_t iev = 0; iev < nEvents; iev++) {
+
+ srcTree->GetEntry(iev);
+ cout << "Event " << iev + 1 << " of " << nEvents << ": " << endl;
+
+ nTracks = esd->GetNumberOfTracks();
+ firstCluster = esd->GetFirstEMCALCluster();
+ lastCluster = esd->GetFirstEMCALCluster() + esd->GetNumberOfEMCALClusters();
+ cout << "Tracks found : " << nTracks << endl;
+ cout << "EMC clusters found: " << lastCluster - firstCluster << endl;
+
+ // create new matches tree
+ TTree *outTree = new TTree(Form("tm_%d", iev), Form("True matches from event %d", iev));
+ outTree->Branch("matches", &match, "label/I:indexT/I:indexC/I:p[3]/D:v[3]/D");
+
+ // external loop on tracks
+ for (Int_t it = 0; it < nTracks; it++) {
+ AliESDtrack *esdTrack = esd->GetTrack(it);
+ // start check to reject kinks
+ if (rejectKinks) {
+ isKink = kFALSE;
+ for (Int_t i = 0; i < 3; i++) {
+ if (esdTrack->GetKinkIndex(i) > 0) isKink = kTRUE;
+ }
+ if (isKink) continue;
+ }
+ // get track GEANT label (to be checked for match)
+ label = esdTrack->GetLabel();
+ if (rejectFakes) {
+ if (label < 0) continue;
+ }
+ else {
+ label = TMath::Abs(label);
+ }
+ // store track data into candidate match variable
+ // anyway, it will be stored in the tree only
+ // if it matches a cluster
+ match.indexT = it;
+ esdTrack->GetPxPyPz(match.p);
+ esdTrack->GetXYZ(match.v);
+ // internal loop on clusters
+ // a counter counts how many true matches are
+ // generated for the same track
+ count = 0;
+ for (Int_t ic = firstCluster; ic < lastCluster; ic++) {
+ AliESDCaloCluster *cl = esd->GetCaloCluster(ic);
+ // reject pseudo-clusters & unmatched clusters
+ if (cl->GetClusterType() != AliESDCaloCluster::kClusterv1) continue;
+ if (cl->GetPrimaryIndex() != label) continue;
+ // if the method reaches this point, we
+ // have found a match to be stored
+ match.label = label;
+ match.indexC = ic;
+ outTree->Fill();
+ count++;
+ }
+ // alert for multiple matches
+ if (count > 0) {
+ cout << "Found " << count << " clusters which match track " << it << " in ESD";
+ if (count > 1) cout << " --> MULTIPLE MATCH";
+ cout << endl;
+ }
+ } // end loop on tracks
+
+ outFile->cd();
+ outTree->Write();
+ delete outTree;
+ }
+
+ outFile->Close();
+}
git://git.uio.no / u/mrichter/AliRoot.git / commitdiff