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>
30 #include "AliESDEvent.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 fVertexOnTheFly(kFALSE),
60 fOnlyPrimaryJpsi(kFALSE),
67 // Default constructor
72 //----------------------------------------------------------------------------
73 AliBtoJPSItoEleAnalysis::~AliBtoJPSItoEleAnalysis() {}
74 //----------------------------------------------------------------------------
75 void AliBtoJPSItoEleAnalysis::ApplySelection(const Char_t *inName,const Char_t *outName) const {
76 // select candidates that pass fBCuts and write them to a new file
78 TFile *inFile = TFile::Open(inName);
80 TTree *treeBin=(TTree*)inFile->Get("TreeB");
81 AliBtoJPSItoEleAnalysis *inAnalysis = (AliBtoJPSItoEleAnalysis*)inFile->Get("BtoJPSItoEleAnalysis");
82 printf("+++\n+++ I N P U T S T A T U S:\n+++\n");
83 inAnalysis->PrintStatus();
86 AliBtoJPSItoEle *d = 0;
87 treeBin->SetBranchAddress("BtoJPSItoEle",&d);
88 Int_t entries = (Int_t)treeBin->GetEntries();
90 printf("+++\n+++ Number of B candidates in input tree: %d\n+++\n",entries);
92 TTree *treeBout = new TTree("TreeB","Tree with selected B candidates");
93 treeBout->Branch("BtoJPSItoEle","AliBtoJPSItoEle",&d,200000,0);
99 for(Int_t i=0; i<entries; i++) {
100 // get event from tree
101 treeBin->GetEvent(i);
103 if(fSim && fOnlySignal && !d->IsSignal()) continue;
105 // check if candidate passes selection (as B or Bbar)
106 if(d->Select(fBCuts,okB)) {
113 AliBtoJPSItoEleAnalysis *outAnalysis = (AliBtoJPSItoEleAnalysis*)inAnalysis->Clone("BtoJPSItoEleAnalysis");
114 outAnalysis->SetBCuts(fBCuts);
115 printf("------------------------------------------\n");
116 printf("+++\n+++ O U T P U T S T A T U S:\n+++\n");
117 outAnalysis->PrintStatus();
119 printf("+++\n+++ Number of B mesons in output tree: %d\n+++\n",nSel);
121 TFile* outFile = new TFile(outName,"recreate");
123 outAnalysis->Write();
128 //----------------------------------------------------------------------------
129 Double_t AliBtoJPSItoEleAnalysis::CalculateTOFmass(Double_t mom,Double_t length,
130 Double_t time) const {
131 // calculated the mass from momentum, track length from vertex to TOF
132 // and time measured by the TOF
133 if(length==0.) return -1000.;
134 Double_t a = time*time/length/length;
136 a = TMath::Sqrt(a-1.);
138 a = -TMath::Sqrt(1.-a);
143 //----------------------------------------------------------------------------
144 void AliBtoJPSItoEleAnalysis::FindCandidates(Int_t evFirst,Int_t evLast,
145 const Char_t *outName) {
146 // Find candidates and calculate parameters
149 TString esdName="AliESDs.root";
150 if(gSystem->AccessPathName(esdName.Data(),kFileExists)) {
151 printf("AliBtoJPSItoEleAnalysis::FindCandidatesESD(): No ESDs file found!\n");
155 TString outName1=outName;
156 TString outName2="nTotEvents.dat";
158 Int_t nTotEv=0,nBrec=0,nBrec1ev=0;
160 Double_t v2[3],mom[6],d0[2];
161 Int_t iTrkP,iTrkN,trkEntries;
162 Int_t nTrksP=0,nTrksN=0;
166 AliESDtrack *postrack = 0;
167 AliESDtrack *negtrack = 0;
169 // create the AliVertexerTracks object
170 // (it will be used only if fVertexOnTheFly=kTrue)
171 AliVertexerTracks *vertexer1 = new AliVertexerTracks();
172 if(fVertexOnTheFly) {
173 // open the mean vertex
174 TFile *invtx = new TFile("AliESDVertexMean.root");
175 AliESDVertex *initVertex = (AliESDVertex*)invtx->Get("vtxmean");
177 vertexer1->SetVtxStart(initVertex);
183 // create tree for reconstructed decayes
184 AliBtoJPSItoEle *ioBtoJPSItoEle=0;
185 TTree *treeB = new TTree("TreeB","Tree with candidates");
186 treeB->Branch("BtoJPSItoEle","AliBtoJPSItoEle",&ioBtoJPSItoEle,200000,0);
188 // open file with tracks
189 TFile *esdFile = TFile::Open(esdName.Data());
190 AliESDEvent* event = new AliESDEvent();
191 TTree* tree = (TTree*) esdFile->Get("esdTree");
193 Error("FindCandidatesESD", "no ESD tree found");
196 event->ReadFromTree(tree);
198 // loop on events in file
199 for(Int_t iEvent = evFirst; iEvent < tree->GetEntries(); iEvent++) {
200 if(iEvent > evLast) break;
201 tree->GetEvent(iEvent);
202 Int_t ev = (Int_t)event->GetEventNumberInFile();
203 printf("--- Finding B -> JPSI -> e+ e- in event %d\n",ev);
204 // count the total number of events
207 trkEntries = (Int_t)event->GetNumberOfTracks();
208 printf(" Number of tracks: %d\n",trkEntries);
209 if(trkEntries<2) continue;
211 // retrieve primary vertex from file
212 if(!event->GetPrimaryVertex()) {
213 printf(" No vertex\n");
216 event->GetPrimaryVertex()->GetXYZ(fV1);
218 // call function which applies sigle-track selection and
219 // separetes positives and negatives
220 TObjArray trksP(trkEntries/2);
221 Int_t *trkEntryP = new Int_t[trkEntries];
222 TObjArray trksN(trkEntries/2);
223 Int_t *trkEntryN = new Int_t[trkEntries];
224 TTree *trkTree = new TTree();
225 SelectTracks(event,trksP,trkEntryP,nTrksP,
226 trksN,trkEntryN,nTrksN);
228 printf(" pos. tracks: %d neg .tracks: %d\n",nTrksP,nTrksN);
233 // LOOP ON POSITIVE TRACKS
234 for(iTrkP=0; iTrkP<nTrksP; iTrkP++) {
235 if(iTrkP%1==0) printf(" Processing positive track number %d of %d\n",iTrkP,nTrksP);
237 // get track from track array
238 postrack = (AliESDtrack*)trksP.UncheckedAt(iTrkP);
239 trkNum[0] = trkEntryP[iTrkP];
241 // LOOP ON NEGATIVE TRACKS
242 for(iTrkN=0; iTrkN<nTrksN; iTrkN++) {
244 // get track from tracks array
245 negtrack = (AliESDtrack*)trksN.UncheckedAt(iTrkN);
246 trkNum[1] = trkEntryN[iTrkN];
249 // ----------- DCA MINIMIZATION ------------------
251 // find the DCA and propagate the tracks to the DCA
252 Double_t b=event->GetMagneticField();
253 AliESDtrack nt(*negtrack), pt(*postrack);
254 dca = nt.PropagateToDCA(&pt,b);
256 // define the AliESDv0 object
257 AliESDv0 vertex2(nt,trkNum[0],pt,trkNum[1]);
259 // get position of the secondary vertex
260 vertex2.GetXYZ(v2[0],v2[1],v2[2]);
261 vertex2.GetPPxPyPz(mom[0],mom[1],mom[2]);
262 vertex2.GetNPxPyPz(mom[3],mom[4],mom[5]);
263 // impact parameters of the tracks w.r.t. the primary vertex
264 // d0[0] = 10000.*pt.GetD(fV1[0],fV1[1],b);
265 // d0[1] = -10000.*nt.GetD(fV1[0],fV1[1],b);
268 // no vertexing if DeltaMass > fMassCut
269 if(fVertexOnTheFly) {
271 if(SelectInvMass(mom)) {
272 // primary vertex from *other* tracks in the event
273 vertexer1->SetFieldkG(event->GetMagneticField());
274 skipped[0] = trkEntryP[iTrkP];
275 skipped[1] = trkEntryN[iTrkN];
276 vertexer1->SetSkipTracks(2,skipped);
277 AliESDVertex *vertex1onfly =
278 (AliESDVertex*)vertexer1->FindPrimaryVertex(event);
279 if(vertex1onfly->GetNContributors()>0) goodVtx1 = kTRUE;
280 vertex1onfly->GetXYZ(fV1);
281 //vertex1onfly->PrintStatus();
286 // impact parameters of the tracks w.r.t. the primary vertex
287 d0[0] = 10000.*pt.GetD(fV1[0],fV1[1],b);
288 d0[1] = -10000.*nt.GetD(fV1[0],fV1[1],b);
290 // create the object AliBtoJPSItoEle
291 AliBtoJPSItoEle theB(ev,trkNum,fV1,v2,dca,mom,d0);
293 if(goodVtx1 && theB.Select(fBCuts,okB)) {
294 // get PID info from ESD
295 AliESDtrack *t0 = (AliESDtrack*)event->GetTrack(trkNum[0]);
297 t0->GetESDpid(esdpid0);
298 if(t0->GetStatus()&AliESDtrack::kTOFpid) {
299 tofmass[0] = CalculateTOFmass(t0->GetP(),
300 t0->GetIntegratedLength(),
305 AliESDtrack *t1 = (AliESDtrack*)event->GetTrack(trkNum[1]);
307 t1->GetESDpid(esdpid1);
308 if(t1->GetStatus()&AliESDtrack::kTOFpid) {
309 tofmass[1] = CalculateTOFmass(t1->GetP(),
310 t1->GetIntegratedLength(),
316 theB.SetPIDresponse(esdpid0,esdpid1);
317 theB.SetTOFmasses(tofmass);
320 ioBtoJPSItoEle=&theB;
328 } // loop on negative tracks
330 } // loop on positive tracks
336 printf(" Number of B candidates: %d\n",nBrec1ev);
337 } // loop on events in file
340 printf("\n+++\n+++ Total number of events: %d\n+++\n",nTotEv);
341 printf("\n+++\n+++ Total number of B candidates: %d\n+++\n",nBrec);
347 // create a copy of this class to be written to output file
348 AliBtoJPSItoEleAnalysis *copy = (AliBtoJPSItoEleAnalysis*)this->Clone("BtoJPSItoEleAnalysis");
350 // add PDG codes to decay tracks in found candidates (in simulation mode)
351 // and store tree in the output file
353 TFile *outroot = new TFile(outName1.Data(),"recreate");
359 printf(" Now adding information from simulation (PDG codes) ...\n");
360 TTree *treeBsim = new TTree("TreeB","Tree with B candidates");
361 SimulationInfo(treeB,treeBsim);
363 TFile *outroot = new TFile(outName1.Data(),"recreate");
370 // write to a file the total number of events
371 FILE *outdat = fopen(outName2.Data(),"w");
372 fprintf(outdat,"%d\n",nTotEv);
377 //-----------------------------------------------------------------------------
378 void AliBtoJPSItoEleAnalysis::PrintStatus() const {
379 // Print parameters being used
381 printf("Preselections:\n");
382 printf(" fPtCut = %f GeV\n",fPtCut);
383 printf(" fd0Cut = %f micron\n",fd0Cut);
384 printf(" fMassCut = %f GeV\n",fMassCut);
385 printf(" fPidCut > %f \n",fPidCut);
386 if(fVertexOnTheFly) printf("Primary vertex on the fly\n");
388 printf("Simulation mode\n");
389 if(fOnlySignal && !(fOnlyPrimaryJpsi)) printf(" Only signal goes to file\n");
390 if(fOnlyPrimaryJpsi && !(fOnlySignal)) printf(" Only primary Jpsi go to file\n");
391 if(fOnlyPrimaryJpsi && fOnlySignal) printf(" Both signal and primary Jpsi go to file\n");
393 printf("Cuts on candidates:\n");
394 printf(" |M-MJPsi| [GeV] < %f\n",fBCuts[0]);
395 printf(" dca [micron] < %f\n",fBCuts[1]);
396 printf(" cosThetaStar < %f\n",fBCuts[2]);
397 printf(" pTP [GeV] > %f\n",fBCuts[3]);
398 printf(" pTN [GeV] > %f\n",fBCuts[4]);
399 printf(" |d0P| [micron] < %f\n",fBCuts[5]);
400 printf(" |d0N| [micron] < %f\n",fBCuts[6]);
401 printf(" d0d0 [micron^2] < %f\n",fBCuts[7]);
402 printf(" cosThetaPoint > %f\n",fBCuts[8]);
406 //-----------------------------------------------------------------------------
407 Bool_t AliBtoJPSItoEleAnalysis::SelectInvMass(const Double_t p[6]) const {
408 // Apply preselection in the invariant mass of the pair
410 Double_t mJPsi = 3.096916;
411 Double_t mel = 0.00510998902;
414 Double_t mom2[2],momTot2;
416 mom2[0] = p[0]*p[0] + p[1]*p[1] + p[2]*p[2];
417 mom2[1] = p[3]*p[3] + p[4]*p[4] + p[5]*p[5];
419 momTot2 = (p[0]+p[3])*(p[0]+p[3])+
420 (p[1]+p[4])*(p[1]+p[4])+
421 (p[2]+p[5])*(p[2]+p[5]);
424 energy[1] = TMath::Sqrt(mel*mel+mom2[1]);
425 energy[0] = TMath::Sqrt(mel*mel+mom2[0]);
427 Double_t minvJPsi = TMath::Sqrt((energy[0]+energy[1])*(energy[0]+energy[1])-momTot2);
429 if(TMath::Abs(minvJPsi-mJPsi) < fMassCut) return kTRUE;
432 //-----------------------------------------------------------------------------
433 void AliBtoJPSItoEleAnalysis::SelectTracks(AliESDEvent *event,
434 TObjArray &trksP,Int_t *trkEntryP,Int_t &nTrksP,
435 TObjArray &trksN,Int_t *trkEntryN,Int_t &nTrksN) const {
436 // Create two TObjArrays with positive and negative tracks and
437 // apply single-track preselection
441 Int_t entr = event->GetNumberOfTracks();
443 // transfer ITS tracks from ESD to arrays and to a tree
444 for(Int_t i=0; i<entr; i++) {
446 AliESDtrack *esdtrack = event->GetTrack(i);
447 UInt_t status = esdtrack->GetStatus();
449 if(!(status&AliESDtrack::kITSin)) continue;
451 // single track selection
452 if(!SingleTrkCuts(*esdtrack,event->GetMagneticField())) continue;
454 if(esdtrack->GetSign()<0) { // negative track
455 trksN.AddLast(esdtrack);
456 trkEntryN[nTrksN] = i;
458 } else { // positive track
459 trksP.AddLast(esdtrack);
460 trkEntryP[nTrksP] = i;
464 } // loop on ESD tracks
468 //-----------------------------------------------------------------------------
469 void AliBtoJPSItoEleAnalysis::SetBCuts(Double_t cut0,Double_t cut1,
470 Double_t cut2,Double_t cut3,Double_t cut4,
471 Double_t cut5,Double_t cut6,
472 Double_t cut7,Double_t cut8) {
473 // Set the cuts for B selection
486 //-----------------------------------------------------------------------------
487 void AliBtoJPSItoEleAnalysis::SetBCuts(const Double_t cuts[9]) {
488 // Set the cuts for B selection
490 for(Int_t i=0; i<9; i++) fBCuts[i] = cuts[i];
494 //-----------------------------------------------------------------------------
496 AliBtoJPSItoEleAnalysis::SingleTrkCuts(const AliESDtrack& trk, Double_t b) const {
497 // Check if track passes some kinematical cuts
498 // Magnetic field "b" (kG)
500 if(TMath::Abs(1./trk.GetParameter()[4]) < fPtCut)
502 if(TMath::Abs(10000.*trk.GetD(fV1[0],fV1[1],b)) < fd0Cut)
504 //select only tracks with the "combined PID"
505 UInt_t status = trk.GetStatus();
506 if ((status&AliESDtrack::kESDpid)==0) return kTRUE;
509 if(r[0] < fPidCut) return kFALSE;
513 //----------------------------------------------------------------------------
514 void AliBtoJPSItoEleAnalysis::MakeTracksRefFile(AliRun *gAlice,
515 Int_t evFirst,Int_t evLast) const {
516 // Create a file with simulation info for the reconstructed tracks
518 TFile *outFile = TFile::Open("BTracksRefFile.root","recreate");
519 TFile *esdFile = TFile::Open("AliESDs.root");
521 AliMC *mc = gAlice->GetMCApp();
529 AliESDEvent* event = new AliESDEvent();
530 TTree* tree = (TTree*) esdFile->Get("esdTree");
531 event->ReadFromTree(tree);
532 // loop on events in file
533 for(Int_t iEvent=evFirst; iEvent<tree->GetEntries(); iEvent++) {
534 if(iEvent>evLast) break;
535 tree->GetEvent(iEvent);
536 Int_t ev = (Int_t)event->GetEventNumberInFile();
538 gAlice->GetEvent(ev);
540 Int_t nentr=(Int_t)event->GetNumberOfTracks();
542 // Tree for true track parameters
544 sprintf(ttname,"Tree_Ref_%d",ev);
545 TTree *reftree = new TTree(ttname,"Tree with true track params");
546 reftree->Branch("rectracks",&reftrk,"lab/I:pdg:mumlab:mumpdg:gmumlab:gmumpdg:mumprongs:Vx/F:Vy:Vz:Px:Py:Pz");
547 // reftree->Branch("rectracks",&reftrk,"lab/I:pdg:mumlab:mumpdg:Vx/F:Vy:Vz:Px:Py:Pz");
549 for(Int_t i=0; i<nentr; i++) {
550 AliESDtrack *esdtrack = (AliESDtrack*)event->GetTrack(i);
551 label = TMath::Abs(esdtrack->GetLabel());
553 part = (TParticle*)mc->Particle(label);
555 reftrk.pdg = part->GetPdgCode();
556 reftrk.mumlab = part->GetFirstMother();
557 if(part->GetFirstMother()>=0) {
558 mumpart = (TParticle*)gAlice->GetMCApp()->Particle(part->GetFirstMother());
559 reftrk.mumpdg = mumpart->GetPdgCode();
560 reftrk.mumprongs = mumpart->GetNDaughters();
561 reftrk.gmumlab = mumpart->GetFirstMother();
562 if(mumpart->GetFirstMother()>=0) {
563 gmumpart = (TParticle*)gAlice->GetMCApp()->Particle(mumpart->GetFirstMother());
564 reftrk.gmumpdg = gmumpart->GetPdgCode();
570 reftrk.gmumlab=part->GetFirstMother(); // If it hasn't any mother, then it has neither Gmother!
571 // reftrk.gmumlab=-1; // If it hasn't any mother, then it has neither Gmother!
573 reftrk.Vx = part->Vx();
574 reftrk.Vy = part->Vy();
575 reftrk.Vz = part->Vz();
576 reftrk.Px = part->Px();
577 reftrk.Py = part->Py();
578 reftrk.Pz = part->Pz();
595 //-----------------------------------------------------------------------------
596 void AliBtoJPSItoEleAnalysis::SimulationInfo(TTree *treeBin,TTree *treeBout) const {
597 // add pdg codes to candidate decay tracks (for sim)
599 TString refFileName("BTracksRefFile.root");
600 if(fSim && gSystem->AccessPathName(refFileName.Data(),kFileExists)) {
601 printf("AliBtoJPSItoEleAnalysis::SimulationInfo: no reference file found!\n");
604 TFile *refFile = TFile::Open(refFileName.Data());
606 Char_t refTreeName[100];
608 Int_t pdg[2],mumpdg[2],mumlab[2],gmumpdg[2],gmumlab[2];
611 // read-in reference tree for event 0 (the only event for Pb-Pb)
612 sprintf(refTreeName,"Tree_Ref_%d",0);
613 TTree *refTree0 = (TTree*)refFile->Get(refTreeName);
614 refTree0->SetBranchAddress("rectracks",&reftrk);
616 AliBtoJPSItoEle *theB = 0;
617 treeBin->SetBranchAddress("BtoJPSItoEle",&theB);
618 treeBout->Branch("BtoJPSItoEle","AliBtoJPSItoEle",&theB,200000,0);
620 Int_t entries = (Int_t)treeBin->GetEntries();
622 for(Int_t i=0; i<entries; i++) {
623 if(i%100==0) printf(" done %d candidates of %d\n",i,entries);
625 treeBin->GetEvent(i);
626 event = theB->EventNo();
628 if(event==0) { // always true for Pb-Pb (avoid to read-in tree every time)
629 refTree0->GetEvent(theB->GetTrkNum(0));
631 mumpdg[0] = reftrk.mumpdg;
632 mumlab[0] = reftrk.mumlab;
633 gmumpdg[0] = reftrk.gmumpdg;
634 gmumlab[0] = reftrk.gmumlab;
635 refTree0->GetEvent(theB->GetTrkNum(1));
637 mumpdg[1] = reftrk.mumpdg;
638 mumlab[1] = reftrk.mumlab;
639 gmumpdg[1] = reftrk.gmumpdg;
640 gmumlab[1] = reftrk.gmumlab;
642 sprintf(refTreeName,"Tree_Ref_%d",event);
643 TTree *refTree = (TTree*)refFile->Get(refTreeName);
644 refTree->SetBranchAddress("rectracks",&reftrk);
645 refTree->GetEvent(theB->GetTrkNum(0));
647 mumpdg[0] = reftrk.mumpdg;
648 mumlab[0] = reftrk.mumlab;
649 gmumpdg[0] = reftrk.gmumpdg;
650 gmumlab[0] = reftrk.gmumlab;
651 refTree->GetEvent(theB->GetTrkNum(1));
653 mumpdg[1] = reftrk.mumpdg;
654 mumlab[1] = reftrk.mumlab;
655 gmumpdg[1] = reftrk.gmumpdg;
656 gmumlab[1] = reftrk.gmumlab;
660 theB->SetPdgCodes(pdg);
661 theB->SetMumPdgCodes(mumpdg);
662 theB->SetGMumPdgCodes(gmumpdg);
664 if(gmumpdg[0]==gmumpdg[1] && // Both GrandMothers are of the same sign
665 (TMath::Abs(gmumpdg[0])==521 || TMath::Abs(gmumpdg[0])==511 || // GrandMother Bplus/Bminus or B0/B0bar
666 TMath::Abs(gmumpdg[0])==523 || TMath::Abs(gmumpdg[0])==513 || // B0s/B0sbar
667 TMath::Abs(gmumpdg[0])==515 || TMath::Abs(gmumpdg[0])==525 || //
668 TMath::Abs(gmumpdg[0])==531 || TMath::Abs(gmumpdg[0])==533 || //
669 TMath::Abs(gmumpdg[0])==535 || TMath::Abs(gmumpdg[0])==541 || //
670 TMath::Abs(gmumpdg[0])==543 || TMath::Abs(gmumpdg[0])==545 || //
671 TMath::Abs(gmumpdg[0])==10521 || TMath::Abs(gmumpdg[0])==10511 || // all possible
672 TMath::Abs(gmumpdg[0])==10523 || TMath::Abs(gmumpdg[0])==10513 || // B mesons
673 TMath::Abs(gmumpdg[0])==20523 || TMath::Abs(gmumpdg[0])==20513 || //
674 TMath::Abs(gmumpdg[0])==10531 || TMath::Abs(gmumpdg[0])==10533 || //
675 TMath::Abs(gmumpdg[0])==20533 || TMath::Abs(gmumpdg[0])==10541 || //
676 TMath::Abs(gmumpdg[0])==20543 || TMath::Abs(gmumpdg[0])==10543 || //
677 TMath::Abs(gmumpdg[0])==4122 || TMath::Abs(gmumpdg[0])==4222 || // All possible B baryons
678 TMath::Abs(gmumpdg[0])==4212 || TMath::Abs(gmumpdg[0])==4112 || // All possible B baryons
679 TMath::Abs(gmumpdg[0])==4224 || TMath::Abs(gmumpdg[0])==4214 || // All possible B baryons
680 TMath::Abs(gmumpdg[0])==4114 || TMath::Abs(gmumpdg[0])==4232 || // All possible B baryons
681 TMath::Abs(gmumpdg[0])==4132 || TMath::Abs(gmumpdg[0])==4322 || // All possible B baryons
682 TMath::Abs(gmumpdg[0])==4312 || TMath::Abs(gmumpdg[0])==4324 || // All possible B baryons
683 TMath::Abs(gmumpdg[0])==4314 || TMath::Abs(gmumpdg[0])==4332 || // All possible B baryons
684 TMath::Abs(gmumpdg[0])==4334 || TMath::Abs(gmumpdg[0])==4412 || // All possible B baryons
685 TMath::Abs(gmumpdg[0])==4422 || TMath::Abs(gmumpdg[0])==4414 || // All possible B baryons
686 TMath::Abs(gmumpdg[0])==4424 || TMath::Abs(gmumpdg[0])==4432 || // All possible B baryons
687 TMath::Abs(gmumpdg[0])==4434 || TMath::Abs(gmumpdg[0])==4444 // All possible B baryons
689 mumpdg[0]==443 && mumpdg[1]== 443 &&
690 mumlab[0]==mumlab[1] &&
691 reftrk.mumprongs==2 &&
692 pdg[0]==-11 && pdg[1]==11
695 else if ( // here consider the case of primary J/psi
696 gmumlab[0]<0 && gmumlab[1]<0 &&
697 mumpdg[0]==443 && mumpdg[1]== 443 &&
698 mumlab[0]==mumlab[1] &&
699 reftrk.mumprongs==2 &&
700 pdg[0]==-11 && pdg[1]==11
701 ) theB->SetJpsiPrimary();
703 // if(!fOnlySignal || theB->IsSignal()) treeBout->Fill();
705 // write it out 1) always if you have not asked for onlySignal or OnlyPrimaryJpsi (or both)
706 // or 2) if you have asked for Signal and it is Signal
707 // or 3) if you have asked for Primary Jpsi and it is a Primary Jpsi
708 if ( (!fOnlySignal && !fOnlyPrimaryJpsi) || (fOnlySignal && theB->IsSignal())
709 || (fOnlyPrimaryJpsi && theB->IsJpsiPrimary()) ) treeBout->Fill();