MUON/MUONcombi.C

   1 void MUONcombi (Int_t evNumber=0)
   2 {
   3 // Dynamically link some shared libs
   4    if (gClassTable->GetID("AliRun") < 0) {
   5       gROOT->LoadMacro("loadlibs.C");
   6       loadlibs();
   7    } else {
   8       delete gAlice;
   9       gAlice = 0;
  10    }
  11 //
  12 // Connect the Root Galice file containing Geometry, Kine and Hits
  13    TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject("galice.root");
  14    if (!file) file = new TFile("galice.root");
  15
  16 // Get AliRun object from file or create it if not on file
  17    if (!gAlice) {
  18       gAlice = (AliRun*)file->Get("gAlice");
  19       if (gAlice) printf("AliRun object found on file\n");
  20       if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
  21    }
  22    TClonesArray *PartArray = gAlice->Particles();
  23    TParticle *Part;
  24
  25
  26 // Import the Kine and Hits Trees for the event evNumber in the file
  27    Int_t nparticles = gAlice->GetEvent(evNumber);
  28    if (nparticles <= 0) return;
  29 //
  30 //
  31    TH1F *dmass = new TH1F("dmass","Dimuon-Mass Distribution"
  32                           ,25,0.,5.);
  33    TH1F *dmassc = new TH1F("dmassc","Dimuon-Mass Distribution"
  34                            ,50,0.,10.);
  35    TH1F *dmassd = new TH1F("dmassd","Dimuon-Mass Distribution"
  36                            ,50,0.,10.);
  37
  38    TH1F *pt     = new TH1F("pt","pT-single"
  39                            ,50,0.,10.);
  40 //
  41 //  Generator Loop
  42 //
  43
  44    AliGenCocktailEntry *Entry, *e1, *e2;
  45    AliGenCocktail*  Cocktail = (AliGenCocktail*) gAlice->Generator();
  46 //                   Single Generator
  47    for (Entry=Cocktail->FirstGenerator();
  48         Entry;
  49         Entry=Cocktail->NextGenerator()
  50        ) {
  51        Entry->PrintInfo();
  52      }
  53 //                   Pairs of Generators
  54 //
  55 // Initialize Combinator
  56    AliDimuCombinator* Combinator = new AliDimuCombinator();
  57    Combinator->SetEtaCut(-5, 5);
  58    Combinator->SetPtMin(0.0);
  59
  60    Int_t i=0;
  61 //
  62 // Single Muon  Loop
  63 //
  64    Combinator->ResetRange();
  65
  66
  67
  68    for(Muon=Combinator->FirstMuon();
  69        Muon;
  70        Muon=Combinator->NextMuon()) {
  71 //
  72        Int_t chfirst = Muon->GetFirstDaughter();
  73        Int_t chlast  = Muon->GetLastDaughter();
  74        Int_t parent  = Muon->GetFirstMother();
  75        Float_t ptm   = Muon->Pt();
  76        Float_t eta   = Muon->Eta();
  77
  78        printf("\n Particle %d Parent %d first child %d last child %d",
  79               i,parent, chfirst, chlast);
  80        printf("\n Particle pt, eta: %f , %f ", ptm, eta);
  81        i++;
  82    }
  83 //
  84 // Di-Muon Loop
  85
  86    Float_t pt1,pt2;
  87    TParticle* Muon1;
  88    TParticle* Muon2;
  89
  90    Combinator->ResetRange();
  91 /*
  92    for (Combinator->FirstMuonPairSelected(Muon1,Muon2);
  93         (Muon1 && Muon2);
  94         Combinator->NextMuonPairSelected(Muon1,Muon2))
  95    {
  96        pt1 = Muon1->Pt();
  97        pt2 = Muon2->Pt();
  98        Float_t mass = Combinator->Mass(Muon1, Muon2);
  99        Float_t wgt  = Combinator->Weight(Muon1, Muon2);
 100        pt->Fill(pt1, wgt);
 101        pt->Fill(pt2, wgt);
 102        Float_t smeared_mass = mass;
 103        Combinator->SmearGauss(0.05*mass, smeared_mass);
 104        if (Combinator->Correlated(Muon1, Muon2)) {
 105            dmassc->Fill(mass, wgt);
 106        } else {
 107            dmass->Fill(mass, wgt);
 108        }
 109     }
 110 */
 111 //
 112 //   Dimuon Loop controlled by Generator Loop
 113 //
 114    for (Cocktail->FirstGeneratorPair(e1,e2);
 115         (e1&&e2);
 116         Cocktail->NextGeneratorPair(e1,e2)
 117        ){
 118        printf("\n ----------------------------------------------------");
 119        e1->PrintInfo();
 120        e2->PrintInfo();
 121        printf("\n ----------------------------------------------------");
 122        Combinator->SetFirstRange (e1->GetFirst(), e1->GetLast());
 123        Combinator->SetSecondRange(e2->GetFirst(), e2->GetLast());
 124        Combinator->SetRate(e1->Rate(), e2->Rate());
 125        for (Combinator->FirstMuonPairSelected(Muon1,Muon2);
 126             (Muon1 && Muon2);
 127             Combinator->NextMuonPairSelected(Muon1,Muon2))
 128        {
 129            pt1 = Muon1->Pt();
 130            pt2 = Muon2->Pt();
 131            Float_t mass = Combinator->Mass(Muon1, Muon2);
 132            Float_t wgt  = Combinator->Weight(Muon1, Muon2);
 133            pt->Fill(pt1, wgt);
 134            pt->Fill(pt2, wgt);
 135            Float_t smeared_mass = mass;
 136            Combinator->SmearGauss(0.05*mass, smeared_mass);
 137            Float_t DecayWeight =
 138                Combinator->DecayProbability(Muon2)*
 139                Combinator->DecayProbability(Muon2);
 140 //         if (TMath::Min(pt1,pt2) > -0.5*TMath::Max(pt1,pt2)+2.) {
 141                if (Combinator->Correlated(Muon1, Muon2)) {
 142                    dmassc->Fill(smeared_mass, wgt);
 143                }
 144                    dmass->Fill(smeared_mass, wgt);
 145                    dmassd->Fill(smeared_mass, wgt*DecayWeight);
 146 //         }
 147        } // Dimuon Loop
 148    }// Generator Loop
 149 //
 150 //Create a canvas, set the view range, show histograms
 151 //
 152    TCanvas *c1 = new TCanvas("c1","Dimuon Plots",400,10,600,700);
 153    TPad* pad11 = new TPad("pad11"," ",0.01,0.51,0.49,0.99);
 154    TPad* pad12 = new TPad("pad12"," ",0.51,0.51,0.99,0.99);
 155    TPad* pad13 = new TPad("pad13"," ",0.01,0.01,0.49,0.49);
 156    TPad* pad14 = new TPad("pad14"," ",0.51,0.01,0.99,0.49);
 157    pad11->SetFillColor(11);
 158    pad12->SetFillColor(11);
 159    pad13->SetFillColor(11);
 160    pad14->SetFillColor(11);
 161    pad11->Draw();
 162    pad12->Draw();
 163    pad13->Draw();
 164    pad14->Draw();
 165
 166    pad11->cd();
 167    pt->SetFillColor(42);
 168    pt->SetXTitle("pT (GeV)");
 169    pt->Draw();
 170
 171    pad12->cd();
 172    pad12->SetLogy(0);
 173    dmass->SetFillColor(42);
 174    dmass->SetXTitle("m (GeV)");
 175    dmass->Draw();
 176
 177    pad13->cd();
 178    pad13->SetLogy(0);
 179    dmassc->SetFillColor(42);
 180    dmassc->SetXTitle("m (GeV)");
 181    dmassc->Draw();
 182
 183    pad14->cd();
 184    pad14->SetLogy(0);
 185    dmassd->SetFillColor(42);
 186    dmassd->SetXTitle("m (GeV)");
 187    dmassd->Draw();
 188
 189
 190 }
 191
 192
 193
 194
 195
 196
 197
 198
 199
 200