1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
16 //----------------------------------------------------------------------------
17 // Implementation of the BtoJPSItoEle reconstruction and analysis class
18 // Note: the two decay tracks are labelled: 0 (positive track)
20 // An example of usage can be found in the macro AliBtoJPSItoEleTest.C
21 // Origin: G.E. Bruno giuseppe.bruno@ba.infn.it
22 // based on Class for charm golden channel (D0->Kpi)
23 //----------------------------------------------------------------------------
29 #include <TParticle.h>
33 #include "AliRunLoader.h"
34 #include "AliVertexerTracks.h"
35 #include "AliESDVertex.h"
37 #include "AliBtoJPSItoEle.h"
38 #include "AliBtoJPSItoEleAnalysis.h"
53 ClassImp(AliBtoJPSItoEleAnalysis)
55 //----------------------------------------------------------------------------
56 AliBtoJPSItoEleAnalysis::AliBtoJPSItoEleAnalysis() {
57 // Default constructor
66 fVertexOnTheFly = kFALSE;
69 fOnlyPrimaryJpsi = kFALSE;
71 //----------------------------------------------------------------------------
72 AliBtoJPSItoEleAnalysis::~AliBtoJPSItoEleAnalysis() {}
73 //----------------------------------------------------------------------------
74 void AliBtoJPSItoEleAnalysis::ApplySelection(const Char_t *inName,const Char_t *outName) const {
75 // select candidates that pass fBCuts and write them to a new file
77 TFile *inFile = TFile::Open(inName);
79 TTree *treeBin=(TTree*)inFile->Get("TreeB");
80 AliBtoJPSItoEleAnalysis *inAnalysis = (AliBtoJPSItoEleAnalysis*)inFile->Get("BtoJPSItoEleAnalysis");
81 printf("+++\n+++ I N P U T S T A T U S:\n+++\n");
82 inAnalysis->PrintStatus();
85 AliBtoJPSItoEle *d = 0;
86 treeBin->SetBranchAddress("BtoJPSItoEle",&d);
87 Int_t entries = (Int_t)treeBin->GetEntries();
89 printf("+++\n+++ Number of B candidates in input tree: %d\n+++\n",entries);
91 TTree *treeBout = new TTree("TreeB","Tree with selected B candidates");
92 treeBout->Branch("BtoJPSItoEle","AliBtoJPSItoEle",&d,200000,0);
98 for(Int_t i=0; i<entries; i++) {
99 // get event from tree
100 treeBin->GetEvent(i);
102 if(fSim && fOnlySignal && !d->IsSignal()) continue;
104 // check if candidate passes selection (as B or Bbar)
105 if(d->Select(fBCuts,okB)) {
112 AliBtoJPSItoEleAnalysis *outAnalysis = (AliBtoJPSItoEleAnalysis*)inAnalysis->Clone("BtoJPSItoEleAnalysis");
113 outAnalysis->SetBCuts(fBCuts);
114 printf("------------------------------------------\n");
115 printf("+++\n+++ O U T P U T S T A T U S:\n+++\n");
116 outAnalysis->PrintStatus();
118 printf("+++\n+++ Number of B mesons in output tree: %d\n+++\n",nSel);
120 TFile* outFile = new TFile(outName,"recreate");
122 outAnalysis->Write();
127 //----------------------------------------------------------------------------
128 Double_t AliBtoJPSItoEleAnalysis::CalculateTOFmass(Double_t mom,Double_t length,
129 Double_t time) const {
130 // calculated the mass from momentum, track length from vertex to TOF
131 // and time measured by the TOF
132 if(length==0.) return -1000.;
133 Double_t a = time*time/length/length;
135 a = TMath::Sqrt(a-1.);
137 a = -TMath::Sqrt(1.-a);
142 //----------------------------------------------------------------------------
143 void AliBtoJPSItoEleAnalysis::FindCandidates(Int_t evFirst,Int_t evLast,
144 const Char_t *outName) {
145 // Find candidates and calculate parameters
148 TString esdName="AliESDs.root";
149 if(gSystem->AccessPathName(esdName.Data(),kFileExists)) {
150 printf("AliBtoJPSItoEleAnalysis::FindCandidatesESD(): No ESDs file found!\n");
154 TString outName1=outName;
155 TString outName2="nTotEvents.dat";
157 Int_t nTotEv=0,nBrec=0,nBrec1ev=0;
159 Double_t v2[3],mom[6],d0[2];
160 Int_t iTrkP,iTrkN,trkEntries;
161 Int_t nTrksP=0,nTrksN=0;
165 AliESDtrack *postrack = 0;
166 AliESDtrack *negtrack = 0;
168 // create the AliVertexerTracks object
169 // (it will be used only if fVertexOnTheFly=kTrue)
170 AliVertexerTracks *vertexer1 = new AliVertexerTracks;
171 if(fVertexOnTheFly) {
172 // open the mean vertex
173 TFile *invtx = new TFile("AliESDVertexMean.root");
174 AliESDVertex *initVertex = (AliESDVertex*)invtx->Get("vtxmean");
176 vertexer1->SetVtxStart(initVertex);
182 // create tree for reconstructed decayes
183 AliBtoJPSItoEle *ioBtoJPSItoEle=0;
184 TTree *treeB = new TTree("TreeB","Tree with candidates");
185 treeB->Branch("BtoJPSItoEle","AliBtoJPSItoEle",&ioBtoJPSItoEle,200000,0);
187 // open file with tracks
188 TFile *esdFile = TFile::Open(esdName.Data());
189 AliESD* event = new AliESD;
190 TTree* tree = (TTree*) esdFile->Get("esdTree");
192 Error("FindCandidatesESD", "no ESD tree found");
195 tree->SetBranchAddress("ESD",&event);
197 // loop on events in file
198 for(Int_t iEvent = evFirst; iEvent < tree->GetEntries(); iEvent++) {
199 if(iEvent > evLast) break;
200 tree->GetEvent(iEvent);
201 Int_t ev = (Int_t)event->GetEventNumberInFile();
202 printf("--- Finding B -> JPSI -> e+ e- in event %d\n",ev);
203 // count the total number of events
206 trkEntries = (Int_t)event->GetNumberOfTracks();
207 printf(" Number of tracks: %d\n",trkEntries);
208 if(trkEntries<2) continue;
210 // retrieve primary vertex from file
211 if(!event->GetPrimaryVertex()) {
212 printf(" No vertex\n");
215 event->GetPrimaryVertex()->GetXYZ(fV1);
217 // call function which applies sigle-track selection and
218 // separetes positives and negatives
219 TObjArray trksP(trkEntries/2);
220 Int_t *trkEntryP = new Int_t[trkEntries];
221 TObjArray trksN(trkEntries/2);
222 Int_t *trkEntryN = new Int_t[trkEntries];
223 TTree *trkTree = new TTree();
224 SelectTracks(event,trksP,trkEntryP,nTrksP,
225 trksN,trkEntryN,nTrksN);
227 printf(" pos. tracks: %d neg .tracks: %d\n",nTrksP,nTrksN);
232 // LOOP ON POSITIVE TRACKS
233 for(iTrkP=0; iTrkP<nTrksP; iTrkP++) {
234 if(iTrkP%1==0) printf(" Processing positive track number %d of %d\n",iTrkP,nTrksP);
236 // get track from track array
237 postrack = (AliESDtrack*)trksP.UncheckedAt(iTrkP);
238 trkNum[0] = trkEntryP[iTrkP];
240 // LOOP ON NEGATIVE TRACKS
241 for(iTrkN=0; iTrkN<nTrksN; iTrkN++) {
243 // get track from tracks array
244 negtrack = (AliESDtrack*)trksN.UncheckedAt(iTrkN);
245 trkNum[1] = trkEntryN[iTrkN];
248 // ----------- DCA MINIMIZATION ------------------
250 // find the DCA and propagate the tracks to the DCA
251 Double_t b=event->GetMagneticField();
252 AliESDtrack nt(*negtrack), pt(*postrack);
253 dca = nt.PropagateToDCA(&pt,b);
255 // define the AliESDv0 object
256 AliESDv0 vertex2(nt,trkNum[0],pt,trkNum[1]);
258 // get position of the secondary vertex
259 vertex2.GetXYZ(v2[0],v2[1],v2[2]);
260 vertex2.GetPPxPyPz(mom[0],mom[1],mom[2]);
261 vertex2.GetNPxPyPz(mom[3],mom[4],mom[5]);
262 // impact parameters of the tracks w.r.t. the primary vertex
263 // d0[0] = 10000.*pt.GetD(fV1[0],fV1[1],b);
264 // d0[1] = -10000.*nt.GetD(fV1[0],fV1[1],b);
267 // no vertexing if DeltaMass > fMassCut
268 if(fVertexOnTheFly) {
270 if(SelectInvMass(mom)) {
271 // primary vertex from *other* tracks in the event
272 skipped[0] = trkEntryP[iTrkP];
273 skipped[1] = trkEntryN[iTrkN];
274 vertexer1->SetSkipTracks(2,skipped);
275 AliESDVertex *vertex1onfly =
276 (AliESDVertex*)vertexer1->FindPrimaryVertex(event);
277 if(vertex1onfly->GetNContributors()>0) goodVtx1 = kTRUE;
278 vertex1onfly->GetXYZ(fV1);
279 //vertex1onfly->PrintStatus();
284 // impact parameters of the tracks w.r.t. the primary vertex
285 d0[0] = 10000.*pt.GetD(fV1[0],fV1[1],b);
286 d0[1] = -10000.*nt.GetD(fV1[0],fV1[1],b);
288 // create the object AliBtoJPSItoEle
289 AliBtoJPSItoEle theB(ev,trkNum,fV1,v2,dca,mom,d0);
291 if(goodVtx1 && theB.Select(fBCuts,okB)) {
292 // get PID info from ESD
293 AliESDtrack *t0 = (AliESDtrack*)event->GetTrack(trkNum[0]);
295 t0->GetESDpid(esdpid0);
296 if(t0->GetStatus()&AliESDtrack::kTOFpid) {
297 tofmass[0] = CalculateTOFmass(t0->GetP(),
298 t0->GetIntegratedLength(),
303 AliESDtrack *t1 = (AliESDtrack*)event->GetTrack(trkNum[1]);
305 t1->GetESDpid(esdpid1);
306 if(t1->GetStatus()&AliESDtrack::kTOFpid) {
307 tofmass[1] = CalculateTOFmass(t1->GetP(),
308 t1->GetIntegratedLength(),
314 theB.SetPIDresponse(esdpid0,esdpid1);
315 theB.SetTOFmasses(tofmass);
318 ioBtoJPSItoEle=&theB;
326 } // loop on negative tracks
328 } // loop on positive tracks
334 printf(" Number of B candidates: %d\n",nBrec1ev);
335 } // loop on events in file
338 printf("\n+++\n+++ Total number of events: %d\n+++\n",nTotEv);
339 printf("\n+++\n+++ Total number of B candidates: %d\n+++\n",nBrec);
345 // create a copy of this class to be written to output file
346 AliBtoJPSItoEleAnalysis *copy = (AliBtoJPSItoEleAnalysis*)this->Clone("BtoJPSItoEleAnalysis");
348 // add PDG codes to decay tracks in found candidates (in simulation mode)
349 // and store tree in the output file
351 TFile *outroot = new TFile(outName1.Data(),"recreate");
357 printf(" Now adding information from simulation (PDG codes) ...\n");
358 TTree *treeBsim = new TTree("TreeB","Tree with B candidates");
359 SimulationInfo(treeB,treeBsim);
361 TFile *outroot = new TFile(outName1.Data(),"recreate");
368 // write to a file the total number of events
369 FILE *outdat = fopen(outName2.Data(),"w");
370 fprintf(outdat,"%d\n",nTotEv);
375 //-----------------------------------------------------------------------------
376 void AliBtoJPSItoEleAnalysis::PrintStatus() const {
377 // Print parameters being used
379 printf("Preselections:\n");
380 printf(" fPtCut = %f GeV\n",fPtCut);
381 printf(" fd0Cut = %f micron\n",fd0Cut);
382 printf(" fMassCut = %f GeV\n",fMassCut);
383 printf(" fPidCut > %f \n",fPidCut);
384 if(fVertexOnTheFly) printf("Primary vertex on the fly\n");
386 printf("Simulation mode\n");
387 if(fOnlySignal && !(fOnlyPrimaryJpsi)) printf(" Only signal goes to file\n");
388 if(fOnlyPrimaryJpsi && !(fOnlySignal)) printf(" Only primary Jpsi go to file\n");
389 if(fOnlyPrimaryJpsi && fOnlySignal) printf(" Both signal and primary Jpsi go to file\n");
391 printf("Cuts on candidates:\n");
392 printf(" |M-MJPsi| [GeV] < %f\n",fBCuts[0]);
393 printf(" dca [micron] < %f\n",fBCuts[1]);
394 printf(" cosThetaStar < %f\n",fBCuts[2]);
395 printf(" pTP [GeV] > %f\n",fBCuts[3]);
396 printf(" pTN [GeV] > %f\n",fBCuts[4]);
397 printf(" |d0P| [micron] < %f\n",fBCuts[5]);
398 printf(" |d0N| [micron] < %f\n",fBCuts[6]);
399 printf(" d0d0 [micron^2] < %f\n",fBCuts[7]);
400 printf(" cosThetaPoint > %f\n",fBCuts[8]);
404 //-----------------------------------------------------------------------------
405 Bool_t AliBtoJPSItoEleAnalysis::SelectInvMass(const Double_t p[6]) const {
406 // Apply preselection in the invariant mass of the pair
408 Double_t mJPsi = 3.096916;
409 Double_t mel = 0.00510998902;
412 Double_t mom2[2],momTot2;
414 mom2[0] = p[0]*p[0] + p[1]*p[1] + p[2]*p[2];
415 mom2[1] = p[3]*p[3] + p[4]*p[4] + p[5]*p[5];
417 momTot2 = (p[0]+p[3])*(p[0]+p[3])+
418 (p[1]+p[4])*(p[1]+p[4])+
419 (p[2]+p[5])*(p[2]+p[5]);
422 energy[1] = TMath::Sqrt(mel*mel+mom2[1]);
423 energy[0] = TMath::Sqrt(mel*mel+mom2[0]);
425 Double_t minvJPsi = TMath::Sqrt((energy[0]+energy[1])*(energy[0]+energy[1])-momTot2);
427 if(TMath::Abs(minvJPsi-mJPsi) < fMassCut) return kTRUE;
430 //-----------------------------------------------------------------------------
431 void AliBtoJPSItoEleAnalysis::SelectTracks(AliESD *event,
432 TObjArray &trksP,Int_t *trkEntryP,Int_t &nTrksP,
433 TObjArray &trksN,Int_t *trkEntryN,Int_t &nTrksN) const {
434 // Create two TObjArrays with positive and negative tracks and
435 // apply single-track preselection
439 Int_t entr = event->GetNumberOfTracks();
441 // transfer ITS tracks from ESD to arrays and to a tree
442 for(Int_t i=0; i<entr; i++) {
444 AliESDtrack *esdtrack = event->GetTrack(i);
445 UInt_t status = esdtrack->GetStatus();
447 if(!(status&AliESDtrack::kITSin)) continue;
449 // single track selection
450 if(!SingleTrkCuts(*esdtrack,event->GetMagneticField())) continue;
452 if(esdtrack->GetSign()<0) { // negative track
453 trksN.AddLast(esdtrack);
454 trkEntryN[nTrksN] = i;
456 } else { // positive track
457 trksP.AddLast(esdtrack);
458 trkEntryP[nTrksP] = i;
462 } // loop on ESD tracks
466 //-----------------------------------------------------------------------------
467 void AliBtoJPSItoEleAnalysis::SetBCuts(Double_t cut0,Double_t cut1,
468 Double_t cut2,Double_t cut3,Double_t cut4,
469 Double_t cut5,Double_t cut6,
470 Double_t cut7,Double_t cut8) {
471 // Set the cuts for B selection
484 //-----------------------------------------------------------------------------
485 void AliBtoJPSItoEleAnalysis::SetBCuts(const Double_t cuts[9]) {
486 // Set the cuts for B selection
488 for(Int_t i=0; i<9; i++) fBCuts[i] = cuts[i];
492 //-----------------------------------------------------------------------------
494 AliBtoJPSItoEleAnalysis::SingleTrkCuts(const AliESDtrack& trk, Double_t b) const {
495 // Check if track passes some kinematical cuts
496 // Magnetic field "b" (kG)
498 if(TMath::Abs(1./trk.GetParameter()[4]) < fPtCut)
500 if(TMath::Abs(10000.*trk.GetD(fV1[0],fV1[1],b)) < fd0Cut)
502 //select only tracks with the "combined PID"
503 UInt_t status = trk.GetStatus();
504 if ((status&AliESDtrack::kESDpid)==0) return kTRUE;
507 if(r[0] < fPidCut) return kFALSE;
511 //----------------------------------------------------------------------------
512 void AliBtoJPSItoEleAnalysis::MakeTracksRefFile(AliRun *gAlice,
513 Int_t evFirst,Int_t evLast) const {
514 // Create a file with simulation info for the reconstructed tracks
516 TFile *outFile = TFile::Open("BTracksRefFile.root","recreate");
517 TFile *esdFile = TFile::Open("AliESDs.root");
519 AliMC *mc = gAlice->GetMCApp();
527 AliESD* event = new AliESD;
528 TTree* tree = (TTree*) esdFile->Get("esdTree");
529 tree->SetBranchAddress("ESD",&event);
530 // loop on events in file
531 for(Int_t iEvent=evFirst; iEvent<tree->GetEntries(); iEvent++) {
532 if(iEvent>evLast) break;
533 tree->GetEvent(iEvent);
534 Int_t ev = (Int_t)event->GetEventNumberInFile();
536 gAlice->GetEvent(ev);
538 Int_t nentr=(Int_t)event->GetNumberOfTracks();
540 // Tree for true track parameters
542 sprintf(ttname,"Tree_Ref_%d",ev);
543 TTree *reftree = new TTree(ttname,"Tree with true track params");
544 reftree->Branch("rectracks",&reftrk,"lab/I:pdg:mumlab:mumpdg:gmumlab:gmumpdg:mumprongs:Vx/F:Vy:Vz:Px:Py:Pz");
545 // reftree->Branch("rectracks",&reftrk,"lab/I:pdg:mumlab:mumpdg:Vx/F:Vy:Vz:Px:Py:Pz");
547 for(Int_t i=0; i<nentr; i++) {
548 AliESDtrack *esdtrack = (AliESDtrack*)event->GetTrack(i);
549 label = TMath::Abs(esdtrack->GetLabel());
551 part = (TParticle*)mc->Particle(label);
553 reftrk.pdg = part->GetPdgCode();
554 reftrk.mumlab = part->GetFirstMother();
555 if(part->GetFirstMother()>=0) {
556 mumpart = (TParticle*)gAlice->GetMCApp()->Particle(part->GetFirstMother());
557 reftrk.mumpdg = mumpart->GetPdgCode();
558 reftrk.mumprongs = mumpart->GetNDaughters();
559 reftrk.gmumlab = mumpart->GetFirstMother();
560 if(mumpart->GetFirstMother()>=0) {
561 gmumpart = (TParticle*)gAlice->GetMCApp()->Particle(mumpart->GetFirstMother());
562 reftrk.gmumpdg = gmumpart->GetPdgCode();
568 reftrk.gmumlab=part->GetFirstMother(); // If it hasn't any mother, then it has neither Gmother!
569 // reftrk.gmumlab=-1; // If it hasn't any mother, then it has neither Gmother!
571 reftrk.Vx = part->Vx();
572 reftrk.Vy = part->Vy();
573 reftrk.Vz = part->Vz();
574 reftrk.Px = part->Px();
575 reftrk.Py = part->Py();
576 reftrk.Pz = part->Pz();
593 //-----------------------------------------------------------------------------
594 void AliBtoJPSItoEleAnalysis::SimulationInfo(TTree *treeBin,TTree *treeBout) const {
595 // add pdg codes to candidate decay tracks (for sim)
597 TString refFileName("BTracksRefFile.root");
598 if(fSim && gSystem->AccessPathName(refFileName.Data(),kFileExists)) {
599 printf("AliBtoJPSItoEleAnalysis::SimulationInfo: no reference file found!\n");
602 TFile *refFile = TFile::Open(refFileName.Data());
604 Char_t refTreeName[100];
606 Int_t pdg[2],mumpdg[2],mumlab[2],gmumpdg[2],gmumlab[2];
609 // read-in reference tree for event 0 (the only event for Pb-Pb)
610 sprintf(refTreeName,"Tree_Ref_%d",0);
611 TTree *refTree0 = (TTree*)refFile->Get(refTreeName);
612 refTree0->SetBranchAddress("rectracks",&reftrk);
614 AliBtoJPSItoEle *theB = 0;
615 treeBin->SetBranchAddress("BtoJPSItoEle",&theB);
616 treeBout->Branch("BtoJPSItoEle","AliBtoJPSItoEle",&theB,200000,0);
618 Int_t entries = (Int_t)treeBin->GetEntries();
620 for(Int_t i=0; i<entries; i++) {
621 if(i%100==0) printf(" done %d candidates of %d\n",i,entries);
623 treeBin->GetEvent(i);
624 event = theB->EventNo();
626 if(event==0) { // always true for Pb-Pb (avoid to read-in tree every time)
627 refTree0->GetEvent(theB->GetTrkNum(0));
629 mumpdg[0] = reftrk.mumpdg;
630 mumlab[0] = reftrk.mumlab;
631 gmumpdg[0] = reftrk.gmumpdg;
632 gmumlab[0] = reftrk.gmumlab;
633 refTree0->GetEvent(theB->GetTrkNum(1));
635 mumpdg[1] = reftrk.mumpdg;
636 mumlab[1] = reftrk.mumlab;
637 gmumpdg[1] = reftrk.gmumpdg;
638 gmumlab[1] = reftrk.gmumlab;
640 sprintf(refTreeName,"Tree_Ref_%d",event);
641 TTree *refTree = (TTree*)refFile->Get(refTreeName);
642 refTree->SetBranchAddress("rectracks",&reftrk);
643 refTree->GetEvent(theB->GetTrkNum(0));
645 mumpdg[0] = reftrk.mumpdg;
646 mumlab[0] = reftrk.mumlab;
647 gmumpdg[0] = reftrk.gmumpdg;
648 gmumlab[0] = reftrk.gmumlab;
649 refTree->GetEvent(theB->GetTrkNum(1));
651 mumpdg[1] = reftrk.mumpdg;
652 mumlab[1] = reftrk.mumlab;
653 gmumpdg[1] = reftrk.gmumpdg;
654 gmumlab[1] = reftrk.gmumlab;
658 theB->SetPdgCodes(pdg);
659 theB->SetMumPdgCodes(mumpdg);
660 theB->SetGMumPdgCodes(gmumpdg);
662 if(gmumpdg[0]==gmumpdg[1] && // Both GrandMothers are of the same sign
663 (TMath::Abs(gmumpdg[0])==521 || TMath::Abs(gmumpdg[0])==511 || // GrandMother Bplus/Bminus or B0/B0bar
664 TMath::Abs(gmumpdg[0])==523 || TMath::Abs(gmumpdg[0])==513 || // B0s/B0sbar
665 TMath::Abs(gmumpdg[0])==515 || TMath::Abs(gmumpdg[0])==525 || //
666 TMath::Abs(gmumpdg[0])==531 || TMath::Abs(gmumpdg[0])==533 || //
667 TMath::Abs(gmumpdg[0])==535 || TMath::Abs(gmumpdg[0])==541 || //
668 TMath::Abs(gmumpdg[0])==543 || TMath::Abs(gmumpdg[0])==545 || //
669 TMath::Abs(gmumpdg[0])==10521 || TMath::Abs(gmumpdg[0])==10511 || // all possible
670 TMath::Abs(gmumpdg[0])==10523 || TMath::Abs(gmumpdg[0])==10513 || // B mesons
671 TMath::Abs(gmumpdg[0])==20523 || TMath::Abs(gmumpdg[0])==20513 || //
672 TMath::Abs(gmumpdg[0])==10531 || TMath::Abs(gmumpdg[0])==10533 || //
673 TMath::Abs(gmumpdg[0])==20533 || TMath::Abs(gmumpdg[0])==10541 || //
674 TMath::Abs(gmumpdg[0])==20543 || TMath::Abs(gmumpdg[0])==10543 || //
675 TMath::Abs(gmumpdg[0])==4122 || TMath::Abs(gmumpdg[0])==4222 || // All possible B baryons
676 TMath::Abs(gmumpdg[0])==4212 || TMath::Abs(gmumpdg[0])==4112 || // All possible B baryons
677 TMath::Abs(gmumpdg[0])==4224 || TMath::Abs(gmumpdg[0])==4214 || // All possible B baryons
678 TMath::Abs(gmumpdg[0])==4114 || TMath::Abs(gmumpdg[0])==4232 || // All possible B baryons
679 TMath::Abs(gmumpdg[0])==4132 || TMath::Abs(gmumpdg[0])==4322 || // All possible B baryons
680 TMath::Abs(gmumpdg[0])==4312 || TMath::Abs(gmumpdg[0])==4324 || // All possible B baryons
681 TMath::Abs(gmumpdg[0])==4314 || TMath::Abs(gmumpdg[0])==4332 || // All possible B baryons
682 TMath::Abs(gmumpdg[0])==4334 || TMath::Abs(gmumpdg[0])==4412 || // All possible B baryons
683 TMath::Abs(gmumpdg[0])==4422 || TMath::Abs(gmumpdg[0])==4414 || // All possible B baryons
684 TMath::Abs(gmumpdg[0])==4424 || TMath::Abs(gmumpdg[0])==4432 || // All possible B baryons
685 TMath::Abs(gmumpdg[0])==4434 || TMath::Abs(gmumpdg[0])==4444 // All possible B baryons
687 mumpdg[0]==443 && mumpdg[1]== 443 &&
688 mumlab[0]==mumlab[1] &&
689 reftrk.mumprongs==2 &&
690 pdg[0]==-11 && pdg[1]==11
693 else if ( // here consider the case of primary J/psi
694 gmumlab[0]<0 && gmumlab[1]<0 &&
695 mumpdg[0]==443 && mumpdg[1]== 443 &&
696 mumlab[0]==mumlab[1] &&
697 reftrk.mumprongs==2 &&
698 pdg[0]==-11 && pdg[1]==11
699 ) theB->SetJpsiPrimary();
701 // if(!fOnlySignal || theB->IsSignal()) treeBout->Fill();
703 // write it out 1) always if you have not asked for onlySignal or OnlyPrimaryJpsi (or both)
704 // or 2) if you have asked for Signal and it is Signal
705 // or 3) if you have asked for Primary Jpsi and it is a Primary Jpsi
706 if ( (!fOnlySignal && !fOnlyPrimaryJpsi) || (fOnlySignal && theB->IsSignal())
707 || (fOnlyPrimaryJpsi && theB->IsJpsiPrimary()) ) treeBout->Fill();