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 **************************************************************************/
19 /***********************************************************************
20 * This code defines the reconstructed cascades visualized with EVE
22 * Ludovic Gaudichet (gaudichet@to.infn.it)
23 ************************************************************************/
27 #include "MCHelixLine.hi"
29 #include <TPolyLine3D.h>
30 #include <TPolyMarker3D.h>
34 #include <Reve/RGTopFrame.h>
45 /***********************************************************************
49 ************************************************************************/
51 const Float_t Cascade::fgkMassPion2 = 0.13956995*0.13956995;
52 const Float_t Cascade::fgkMassKaon2 = 0.493677*0.493677;
53 const Float_t Cascade::fgkMassProton2 = 0.93827231*0.93827231;
54 const Float_t Cascade::fgkMassLambda2 = 1.115683*1.115683;
56 ClassImp(Reve::Cascade)
84 fCasCosPointingAngle(999),
87 fPolyLinePos.SetLineColor(2); // red
88 fPolyLineNeg.SetLineColor(7); // light blue
89 fPolyLineBach.SetLineColor(4); // blue
98 Cascade::Cascade(TrackRnrStyle* rs) :
123 fCasCosPointingAngle(999),
126 fMarkerColor = fRnrStyle->GetColor();
127 fPolyLineV0.SetLineColor(fMarkerColor);
128 fPolyLinePos.SetLineColor(2); // red
129 fPolyLineNeg.SetLineColor(7); // light blue
130 fPolyLineBach.SetLineColor(4); // blue
132 fMainColorPtr = &fMarkerColor;
141 for (vpPathMark_i i=fPathMarksNeg.begin(); i!=fPathMarksNeg.end(); ++i)
143 for (vpPathMark_i i=fPathMarksPos.begin(); i!=fPathMarksPos.end(); ++i)
145 for (vpPathMark_i i=fPathMarksBach.begin(); i!=fPathMarksBach.end(); ++i)
148 void Cascade::Reset(TPolyLine3D* polyLine) {
149 //polyLine->SetPolyLine(n_points);
150 polyLine->SetPolyLine(0);
154 //______________________________________________________________________
155 void Cascade::SetDecayLength(Float_t primx, Float_t primy, Float_t primz) {
158 Float_t dx = fV_decay.x-primx;
159 Float_t dy = fV_decay.y-primy;
160 Float_t dz = fV_decay.z-primz;
162 fCasDecayLength = sqrt(dx*dx+dy*dy+dz*dz);
163 // This is probably wrong but I can only do this for now
164 Float_t distNorm = fCasDecayLength/GetMomentum();
165 fV_birth.x = fV_decay.x - distNorm*GetPx();
166 fV_birth.y = fV_decay.y - distNorm*GetPy();
167 fV_birth.z = fV_decay.z - distNorm*GetPz();
169 fV_bach.x = fV_decay.x;
170 fV_bach.y = fV_decay.y;
171 fV_bach.z = fV_decay.z;
175 //______________________________________________________________________
176 void Cascade::MakeTrack(vpPathMark_t& pathMark, Reve::Vector& vtx, Reve::Vector& p,
177 Int_t charge, Float_t beta, TPolyLine3D& polyLine) {
179 TrackRnrStyle& RS((fRnrStyle != 0) ? *fRnrStyle : TrackRnrStyle::fgDefStyle);
181 Float_t px = p.x, py = p.y, pz = p.z;
189 std::vector<MCVertex> track_points;
190 Bool_t decay = kFALSE;
192 if ((TMath::Abs(vtx.z) > RS.fMaxZ) || (vtx.x*vtx.x + vtx.y*vtx.y > RS.fMaxR*RS.fMaxR))
195 if (TMath::Abs(RS.fMagField) > 1e-5) {
197 // Charged particle in magnetic field
199 Float_t a = RS.fgkB2C * RS.fMagField * charge;
201 MCHelix helix(fRnrStyle, &mc_v0, TMath::C()*beta, &track_points, a); //m->cm
202 helix.Init(TMath::Sqrt(px*px+py*py), pz);
204 if(!pathMark.empty()){
205 for(std::vector<Reve::PathMark*>::iterator i=pathMark.begin();
206 i!=pathMark.end(); ++i) {
208 Reve::PathMark* pm = *i;
210 if(RS.fFitDaughters && pm->type == Reve::PathMark::Daughter){
211 if(TMath::Abs(pm->V.z) > RS.fMaxZ
212 || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR )
215 //printf("%s fit daughter \n", fName.Data());
216 helix.LoopToVertex(p.x, p.y, p.z, pm->V.x, pm->V.y, pm->V.z);
221 if(RS.fFitDecay && pm->type == Reve::PathMark::Decay){
223 if(TMath::Abs(pm->V.z) > RS.fMaxZ
224 || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR )
226 helix.LoopToVertex(p.x, p.y, p.z, pm->V.x, pm->V.y, pm->V.z);
234 if(!decay || RS.fFitDecay == kFALSE){
235 helix.LoopToBounds(px,py,pz);
236 // printf("%s loop to bounds \n",fName.Data() );
241 // Neutral particle or no field
243 MCLine line(fRnrStyle, &mc_v0, TMath::C()*beta, &track_points);
245 if(!pathMark.empty()) {
246 for(std::vector<Reve::PathMark*>::iterator i=pathMark.begin();
247 i!=pathMark.end(); ++i) {
248 Reve::PathMark* pm = *i;
250 if(RS.fFitDaughters && pm->type == Reve::PathMark::Daughter){
251 if(TMath::Abs(pm->V.z) > RS.fMaxZ
252 || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR )
254 line.GotoVertex(pm->V.x, pm->V.y, pm->V.z);
260 if(RS.fFitDecay && pm->type == Reve::PathMark::Decay){
261 if(TMath::Abs(pm->V.z) > RS.fMaxZ
262 || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR )
264 line.GotoVertex(pm->V.x, pm->V.y, pm->V.z);
272 if(!decay || RS.fFitDecay == kFALSE)
273 line.GotoBounds(px,py,pz);
278 for(std::vector<MCVertex>::iterator i=track_points.begin();
279 i!=track_points.end(); ++i) {
280 polyLine.SetNextPoint(i->x,i->y, i->z);
285 //______________________________________________________________________
286 void Cascade::MakeV0path() {
289 mc_v0.x = (fV_neg.x+fV_pos.x)/2;
290 mc_v0.y = (fV_neg.y+fV_pos.y)/2;
291 mc_v0.z = (fV_neg.z+fV_pos.z)/2;
294 std::vector<MCVertex> track_points;
295 MCLine line(fRnrStyle, &mc_v0, TMath::C()*0.99, &track_points);
297 line.GotoVertex(fV_decay.x,fV_decay.y,fV_decay.z);
300 for(std::vector<MCVertex>::iterator i=track_points.begin();
301 i!=track_points.end(); ++i) {
302 fPolyLineV0.SetNextPoint(i->x,i->y, i->z);
307 //______________________________________________________________________
308 void Cascade::MakeCasPath() {
311 mc_v0.x = fV_birth.x;
312 mc_v0.y = fV_birth.y;
313 mc_v0.z = fV_birth.z;
316 std::vector<MCVertex> track_points;
317 MCLine line(fRnrStyle, &mc_v0, TMath::C()*0.99, &track_points);
319 line.GotoVertex(fV_decay.x,fV_decay.y,fV_decay.z);
321 Reset(&fPolyLineCas);
322 for(std::vector<MCVertex>::iterator i=track_points.begin();
323 i!=track_points.end(); ++i) {
324 fPolyLineCas.SetNextPoint(i->x,i->y, i->z);
329 //______________________________________________________________________
330 void Cascade::MakeCascade()
332 SetNextPoint(fV_neg.x, fV_neg.y, fV_neg.z);
333 SetNextPoint(fV_decay.x, fV_decay.y, fV_decay.z);
335 MakeTrack(fPathMarksNeg, fV_neg, fP_neg, -1, fBeta_neg, fPolyLineNeg);
336 MakeTrack(fPathMarksPos, fV_pos, fP_pos, 1, fBeta_pos, fPolyLinePos);
338 MakeTrack(fPathMarksBach, fV_bach, fP_bach, 1, fBeta_bach, fPolyLineBach);
340 MakeTrack(fPathMarksBach, fV_bach, fP_bach, -1, -fBeta_bach, fPolyLineBach);
348 //______________________________________________________________________
349 Float_t Cascade::GetCasAlphaArmenteros() const
351 Float_t px = GetPx(), py = GetPy(), pz = GetPz();
352 Float_t posXcas, negXcas;
355 posXcas = fP_bach.x*px + fP_bach.y*py + fP_bach.z*pz;
356 negXcas = (fP_neg.x+fP_pos.x)*px + (fP_neg.y+fP_pos.y)*py + (fP_neg.z+fP_pos.z)*pz;
358 posXcas = (fP_neg.x+fP_pos.x)*px + (fP_neg.y+fP_pos.y)*py + (fP_neg.z+fP_pos.z)*pz;
359 negXcas = fP_bach.x*px + fP_bach.y*py + fP_bach.z*pz;
362 if (posXcas + negXcas > 1.e-39)
363 return (posXcas - negXcas)/(posXcas + negXcas);
368 //______________________________________________________________________
369 Float_t Cascade::GetCasPtArmenteros() const
371 Float_t px = GetPx(), py = GetPy(), pz = GetPz();
372 Float_t p2 = px*px + py*py + pz*pz;
373 if (p2 < 1.e-39) return -999;
375 Float_t posXcas, posP2;
378 posXcas = fP_bach.x*px + fP_bach.y*py + fP_bach.z*pz;
381 posXcas = (fP_neg.x+fP_pos.x)*px + (fP_neg.y+fP_pos.y)*py + (fP_neg.z+fP_pos.z)*pz;
384 return sqrt( posP2 - posXcas*posXcas/p2 );
389 /***********************************************************************
393 ************************************************************************/
395 ClassImp(Reve::CascadeList)
398 //______________________________________________________________________
399 CascadeList::CascadeList(TrackRnrStyle* rs) :
404 fRnrV0Daughters(kTRUE),
417 //______________________________________________________________________
418 CascadeList::CascadeList(const Text_t* name, TrackRnrStyle* rs) :
423 fRnrV0Daughters(kTRUE),
437 //______________________________________________________________________
438 void CascadeList::Init()
441 if (fRnrStyle== 0) fRnrStyle = new TrackRnrStyle;
442 SetMainColorPtr(&fRnrStyle->fColor);
445 fMin[0] = 0; fMax[0] = 5; // Xi mass
446 fMin[1] = 0; fMax[1] = 5; // Omega mass
447 fMin[2] = 0; fMax[2] = 1e5; // Index
448 fMin[3] = 0.8; fMax[3] = 1; // cosPointingAngle
449 fMin[4] = 0; fMax[4] = 5; // bachV0DCA
450 fMin[5] = 0; fMax[5] = 100; // radius
451 fMin[6] = 0; fMax[6] = 10; // Pt
452 fMin[7] = -2; fMax[7] = 2; // PseudoRapidity
453 fMin[8] = 0; fMax[8] = 10; // negPt
454 fMin[9] = -2; fMax[9] = 2; // negEta
455 fMin[10] = 0; fMax[10] = 10; // posPt
456 fMin[11] = -2; fMax[11] = 2; // posEta
457 fMin[12] = 0; fMax[12] = 10; // bachPt
458 fMin[13] = -2; fMax[13] = 2; // backEta
461 fHist[0] = new TH1F(ch,ch, 100, fMin[0], fMax[0]);
463 fHist[1] = new TH1F(ch,ch, 100, fMin[1], fMax[1]);
465 fHist[2] = new TH1F(ch,ch, 100, fMin[2], fMax[2]);
467 ch = "cosPointingAngle";
468 fHist[3] = new TH1F(ch,ch, 100, fMin[3], fMax[3]);
470 fHist[4] = new TH1F(ch,ch, 100, fMin[4], fMax[4]);
472 fHist[5] = new TH1F(ch,ch, 100, fMin[5], fMax[5]);
474 fHist[6] = new TH1F(ch,ch, 100, fMin[6], fMax[6]);
475 ch = "PseudoRapidity";
476 fHist[7] = new TH1F(ch,ch, 100, fMin[7], fMax[7]);
479 fHist[8] = new TH1F(ch,ch, 100, fMin[8], fMax[8]);
481 fHist[9] = new TH1F(ch,ch, 100, fMin[9], fMax[9]);
483 fHist[10] = new TH1F(ch,ch, 100, fMin[10], fMax[10]);
485 fHist[11] = new TH1F(ch,ch, 100, fMin[11], fMax[11]);
487 fHist[12] = new TH1F(ch,ch, 100, fMin[12], fMax[12]);
489 fHist[13] = new TH1F(ch,ch, 100, fMin[13], fMax[13]);
495 ch = "ArmenterosPodolansky";
496 fHist2D[0] = new TH2F(ch,ch, 70, fMinX[0], fMaxX[0], 70,
499 for (Int_t i=0; i<fgkNcutVar; i++) {
500 fHist[i]->GetXaxis()->SetLabelSize(0.07);
501 fHist[i]->GetYaxis()->SetLabelSize(0.07);
502 fHist[i]->SetStats(0);
504 for (Int_t i=0; i<fgkNcutVar2D; i++) {
505 fHist2D[i]->GetXaxis()->SetLabelSize(0.07);
506 fHist2D[i]->GetYaxis()->SetLabelSize(0.07);
507 fHist2D[i]->SetStats(0);
512 //______________________________________________________________________
513 void CascadeList::Paint(Option_t* option) {
517 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
518 if((*i)->GetRnrSelf()) {
519 ((Cascade*)(*i))->PaintBachelor(option);
524 if(fRnrV0Daughters) {
525 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
526 if((*i)->GetRnrSelf()) {
527 ((Cascade*)(*i))->PaintV0Daughters(option);
533 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
534 if((*i)->GetRnrSelf()) {
535 ((Cascade*)(*i))->PaintV0Path(option);
541 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
542 if((*i)->GetRnrSelf()) {
543 ((Cascade*)(*i))->Paint(option);
549 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
550 if((*i)->GetRnrSelf()) {
551 ((Cascade*)(*i))->PaintCasPath(option);
556 } // end if(fRnrSelf)
560 //______________________________________________________________________
562 void CascadeList::AddElement(RenderElement* el)
564 static const Exc_t eH("CascadeList::AddElement ");
565 if (dynamic_cast<Cascade*>(el) == 0)
566 throw(eH + "new element not a Cascade.");
567 RenderElement::AddElement(el);
572 void CascadeList::SetRnrV0vtx(Bool_t rnr)
578 void CascadeList::SetRnrV0path(Bool_t rnr)
584 void CascadeList::SetRnrV0Daughters(Bool_t rnr)
586 fRnrV0Daughters = rnr;
591 void CascadeList::SetRnrCasPath(Bool_t rnr)
597 void CascadeList::SetRnrCasVtx(Bool_t rnr)
603 void CascadeList::SetRnrBachelor(Bool_t rnr)
610 //______________________________________________________________________
612 void CascadeList::MakeCascades()
614 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
615 ((Cascade*)(*i))->MakeCascade();
620 //_________________________________________________________________________
621 void CascadeList::AdjustHist(Int_t iHist) {
623 if ((iHist<0)||(iHist>=fgkNcutVar)) return;
624 if (! fHist[iHist]) return;
626 TString name = fHist[iHist]->GetName();
627 Int_t nBin = fHist[iHist]->GetXaxis()->GetNbins();
629 fHist[iHist] = new TH1F(name.Data(), name.Data(), nBin, GetMin(iHist),
631 fHist[iHist]->GetXaxis()->SetLabelSize(0.07);
632 fHist[iHist]->GetYaxis()->SetLabelSize(0.07);
633 fHist[iHist]->SetStats(0);
637 //______________________________________________________________________
638 void CascadeList::UnFill(Cascade* cas) {
641 Int_t bin = fHist[0]->GetXaxis()->FindBin(cas->GetXiMass());
642 fHist[0]->SetBinContent( bin, fHist[0]->GetBinContent(bin)-1 );
644 bin = fHist[1]->GetXaxis()->FindBin( cas->GetOmegaMass() );
645 fHist[1]->SetBinContent( bin, fHist[1]->GetBinContent(bin)-1 );
648 bin = fHist[2]->GetXaxis()->FindBin( cas->GetESDIndex() );
649 fHist[2]->SetBinContent( bin, fHist[2]->GetBinContent(bin)-1 );
651 bin = fHist[3]->GetXaxis()->FindBin( cas->GetCasCosPointingAngle() );
652 fHist[3]->SetBinContent( bin, fHist[3]->GetBinContent(bin)-1 );
654 bin = fHist[4]->GetXaxis()->FindBin( cas->GetDCA_v0_Bach() );
655 fHist[4]->SetBinContent( bin, fHist[4]->GetBinContent(bin)-1 );
657 bin = fHist[5]->GetXaxis()->FindBin( cas->GetRadius() );
658 fHist[5]->SetBinContent( bin, fHist[5]->GetBinContent(bin)-1 );
660 bin = fHist[6]->GetXaxis()->FindBin( cas->GetPt() );
661 fHist[6]->SetBinContent( bin, fHist[6]->GetBinContent(bin)-1 );
663 bin = fHist[7]->GetXaxis()->FindBin( cas->GetPseudoRapidity() );
664 fHist[7]->SetBinContent( bin, fHist[7]->GetBinContent(bin)-1 );
667 bin = fHist[8]->GetXaxis()->FindBin( cas->GetNegPt() );
668 fHist[8]->SetBinContent( bin, fHist[8]->GetBinContent(bin)-1 );
670 bin = fHist[9]->GetXaxis()->FindBin( cas->GetNegPseudoRapidity() );
671 fHist[9]->SetBinContent( bin, fHist[9]->GetBinContent(bin)-1 );
673 bin = fHist[10]->GetXaxis()->FindBin( cas->GetPosPt() );
674 fHist[10]->SetBinContent( bin, fHist[10]->GetBinContent(bin)-1 );
676 bin = fHist[11]->GetXaxis()->FindBin( cas->GetPosPseudoRapidity() );
677 fHist[11]->SetBinContent( bin, fHist[11]->GetBinContent(bin)-1 );
679 bin = fHist[12]->GetXaxis()->FindBin( cas->GetBachPt() );
680 fHist[12]->SetBinContent( bin, fHist[12]->GetBinContent(bin)-1 );
682 bin = fHist[13]->GetXaxis()->FindBin( cas->GetBachPseudoRapidity() );
683 fHist[13]->SetBinContent( bin, fHist[13]->GetBinContent(bin)-1 );
686 bin = fHist2D[0]->GetXaxis()->FindBin( cas->GetCasAlphaArmenteros() );
687 Int_t binY = fHist2D[0]->GetYaxis()->FindBin( cas->GetCasPtArmenteros() );
688 fHist2D[0]->SetBinContent( bin, binY, fHist2D[0]->GetBinContent(bin,binY)-1 );
692 //______________________________________________________________________
693 void CascadeList::Filter(Cascade* cas) {
695 Float_t xiMass = cas->GetXiMass();
696 if ((xiMass<fMin[0])||(xiMass>fMax[0])) return;
698 Float_t omegaMass = cas->GetOmegaMass();
699 if ( (omegaMass<fMin[1])||(omegaMass>fMax[1]) ) return;
702 Float_t index = cas->GetESDIndex();
703 if ( (index<fMin[2])||(index>fMax[2]) ) return;
705 Float_t cosPointingAngle = cas->GetCasCosPointingAngle();
706 if ( (cosPointingAngle<fMin[3])||(cosPointingAngle>fMax[3]) ) return;
708 Float_t bachV0DCA = cas->GetDCA_v0_Bach();
709 if ( (bachV0DCA<fMin[4])||(bachV0DCA>fMax[4]) ) return;
711 Float_t radius = cas->GetRadius();
712 if ( (radius<fMin[5])||(radius>fMax[5]) ) return;
714 Float_t pt = cas->GetPt();
715 if ( (pt<fMin[6])||(pt>fMax[6]) ) return;
717 Float_t pseudoRapidity = cas->GetPseudoRapidity();
718 if ( (pseudoRapidity<fMin[7])||(pseudoRapidity>fMax[7]) ) return;
720 Float_t negPt = cas->GetNegPt();
721 if ( (negPt<fMin[8])||(negPt>fMax[8]) ) return;
723 Float_t negEta = cas->GetNegPseudoRapidity();
724 if ( (negEta<fMin[9])||(negEta>fMax[9]) ) return;
726 Float_t posPt = cas->GetPosPt();
727 if ( (posPt<fMin[10])||(posPt>fMax[10]) ) return;
729 Float_t posEta = cas->GetPosPseudoRapidity();
730 if ( (posEta<fMin[11])||(posEta>fMax[11]) ) return;
732 Float_t bachPt = cas->GetBachPt();
733 if ( (bachPt<fMin[12])||(bachPt>fMax[12]) ) return;
735 Float_t bachEta = cas->GetBachPseudoRapidity();
736 if ( (bachEta<fMin[13])||(bachEta>fMax[13]) ) return;
738 cas->SetRnrSelf(kTRUE);
739 fHist[0]->Fill(xiMass);
740 fHist[1]->Fill(omegaMass);
741 fHist[2]->Fill(index);
742 fHist[3]->Fill(cosPointingAngle);
743 fHist[4]->Fill(bachV0DCA);
744 fHist[5]->Fill(radius);
746 fHist[7]->Fill(pseudoRapidity);
747 fHist[8]->Fill(negPt);
748 fHist[9]->Fill(negEta);
749 fHist[10]->Fill(posPt);
750 fHist[11]->Fill(posEta);
751 fHist[12]->Fill(bachPt);
752 fHist[13]->Fill(bachEta);
754 fHist2D[0]->Fill(cas->GetCasAlphaArmenteros(), cas->GetCasPtArmenteros() );
757 //______________________________________________________________________
758 void CascadeList::FilterAll() {
760 for (Int_t i=0; i<fgkNcutVar; i++)
763 for (Int_t i=0; i<fgkNcutVar2D; i++)
767 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
769 myCas = (Cascade*)(*i);
775 //______________________________________________________________________
776 void CascadeList::GetCasIndexRange(Int_t &imin, Int_t &imax) {
780 List_i i = fChildren.begin();
781 myCas = (Cascade*)(*i);
782 index = myCas->GetESDIndex();
786 for(; i!=fChildren.end(); ++i) {
788 myCas = (Cascade*)(*i);
789 index = myCas->GetESDIndex();
790 if (index<imin) imin = index;
791 if (index>imax) imax = index;
795 //______________________________________________________________________
796 void CascadeList::XiMassFilter(Float_t min, Float_t max) {
806 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
808 myCas = (Cascade*)(*i);
809 val = myCas->GetXiMass();
810 wasSelected = myCas->GetRnrSelf();
811 isSelected = ( (val>=min) && (val<=max) );
816 myCas->SetRnrSelf(isSelected);
819 if (isSelected) Filter(myCas);
824 //______________________________________________________________________
825 void CascadeList::OmegaMassFilter(Float_t min, Float_t max) {
835 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
837 myCas = (Cascade*)(*i);
838 val = myCas->GetOmegaMass();
839 wasSelected = myCas->GetRnrSelf();
840 isSelected = ( (val>=min) && (val<=max) );
845 myCas->SetRnrSelf(isSelected);
848 if (isSelected) Filter(myCas);
853 //______________________________________________________________________
854 void CascadeList::IndexFilter(Float_t min, Float_t max) {
864 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
866 myCas = (Cascade*)(*i);
867 val = myCas->GetESDIndex();
868 wasSelected = myCas->GetRnrSelf();
869 isSelected = ( (val>=min) && (val<=max) );
874 myCas->SetRnrSelf(isSelected);
877 if (isSelected) Filter(myCas);
882 //______________________________________________________________________
883 void CascadeList::CosPointingFilter(Float_t min, Float_t max) {
893 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
895 myCas = (Cascade*)(*i);
896 val = myCas->GetCasCosPointingAngle();
897 wasSelected = myCas->GetRnrSelf();
898 isSelected = ( (val>=min) && (val<=max) );
903 myCas->SetRnrSelf(isSelected);
906 if (isSelected) Filter(myCas);
912 //______________________________________________________________________
913 void CascadeList::BachV0DCAFilter(Float_t min, Float_t max) {
923 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
925 myCas = (Cascade*)(*i);
926 val = myCas->GetDCA_v0_Bach();
927 wasSelected = myCas->GetRnrSelf();
928 isSelected = ( (val>=min) && (val<=max) );
933 myCas->SetRnrSelf(isSelected);
936 if (isSelected) Filter(myCas);
942 //______________________________________________________________________
943 void CascadeList::RadiusFilter(Float_t min, Float_t max) {
953 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
955 myCas = (Cascade*)(*i);
956 val = myCas->GetRadius();
957 wasSelected = myCas->GetRnrSelf();
958 isSelected = ( (val>=min) && (val<=max) );
963 myCas->SetRnrSelf(isSelected);
966 if (isSelected) Filter(myCas);
972 //______________________________________________________________________
973 void CascadeList::PtFilter(Float_t min, Float_t max) {
983 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
985 myCas = (Cascade*)(*i);
986 val = myCas->GetPt();
987 wasSelected = myCas->GetRnrSelf();
988 isSelected = ( (val>=min) && (val<=max) );
993 myCas->SetRnrSelf(isSelected);
996 if (isSelected) Filter(myCas);
1002 //______________________________________________________________________
1003 void CascadeList::PseudoRapFilter(Float_t min, Float_t max) {
1013 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
1015 myCas = (Cascade*)(*i);
1016 val = myCas->GetPseudoRapidity();
1017 wasSelected = myCas->GetRnrSelf();
1018 isSelected = ( (val>=min) && (val<=max) );
1023 myCas->SetRnrSelf(isSelected);
1026 if (isSelected) Filter(myCas);
1032 //______________________________________________________________________
1033 void CascadeList::NegPtFilter(Float_t min, Float_t max) {
1043 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
1045 myCas = (Cascade*)(*i);
1046 val = myCas->GetNegPt();
1047 wasSelected = myCas->GetRnrSelf();
1048 isSelected = ( (val>=min) && (val<=max) );
1053 myCas->SetRnrSelf(isSelected);
1056 if (isSelected) Filter(myCas);
1062 //______________________________________________________________________
1063 void CascadeList::NegEtaFilter(Float_t min, Float_t max) {
1073 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
1075 myCas = (Cascade*)(*i);
1076 val = myCas->GetNegPseudoRapidity();
1077 wasSelected = myCas->GetRnrSelf();
1078 isSelected = ( (val>=min) && (val<=max) );
1083 myCas->SetRnrSelf(isSelected);
1086 if (isSelected) Filter(myCas);
1092 //______________________________________________________________________
1093 void CascadeList::PosPtFilter(Float_t min, Float_t max) {
1103 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
1105 myCas = (Cascade*)(*i);
1106 val = myCas->GetPosPt();
1107 wasSelected = myCas->GetRnrSelf();
1108 isSelected = ( (val>=min) && (val<=max) );
1113 myCas->SetRnrSelf(isSelected);
1116 if (isSelected) Filter(myCas);
1121 //______________________________________________________________________
1122 void CascadeList::PosEtaFilter(Float_t min, Float_t max) {
1132 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
1134 myCas = (Cascade*)(*i);
1135 val = myCas->GetPosPseudoRapidity();
1136 wasSelected = myCas->GetRnrSelf();
1137 isSelected = ( (val>=min) && (val<=max) );
1142 myCas->SetRnrSelf(isSelected);
1145 if (isSelected) Filter(myCas);
1151 //______________________________________________________________________
1152 void CascadeList::BachPtFilter(Float_t min, Float_t max) {
1162 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
1164 myCas = (Cascade*)(*i);
1165 val = myCas->GetBachPt();
1166 wasSelected = myCas->GetRnrSelf();
1167 isSelected = ( (val>=min) && (val<=max) );
1172 myCas->SetRnrSelf(isSelected);
1175 if (isSelected) Filter(myCas);
1180 //______________________________________________________________________
1181 void CascadeList::BachEtaFilter(Float_t min, Float_t max) {
1191 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
1193 myCas = (Cascade*)(*i);
1194 val = myCas->GetBachPseudoRapidity();
1195 wasSelected = myCas->GetRnrSelf();
1196 isSelected = ( (val>=min) && (val<=max) );
1201 myCas->SetRnrSelf(isSelected);
1204 if (isSelected) Filter(myCas);