#include "AliTPCCorrection.h"
#include "AliGRPObject.h"
#include "AliTPCTransform.h"
-
-
+#include "AliAnalysisManager.h"
-void AliTPCcalibV0::ProcessESD(AliESDEvent *esd, AliStack *stack){
+void AliTPCcalibV0::ProcessESD(AliESDEvent *esd){
//
//
//
if (TMath::Abs(AliTracker::GetBz())<1) return;
DumpToTree(esd);
DumpToTreeHPT(esd);
- return;
- //
- MakeV0s();
- if (stack) {
- fStack = stack;
- MakeMC();
- }else{
- fStack =0;
- }
}
void AliTPCcalibV0::DumpToTreeHPT(AliESDEvent *esd){
//
if (TMath::Abs(AliTracker::GetBz())<1) return;
const Int_t kMinCluster=110;
- const Float_t kMinPt =3.;
+ const Float_t kMinPt =4.;
AliESDfriend *esdFriend=static_cast<AliESDfriend*>(esd->FindListObject("AliESDfriend"));
- if (!esdFriend) {
- Printf("ERROR: esdFriend not available");
- return;
- }
+// if (!esdFriend) {
+// Printf("ERROR: esdFriend not available");
+// return;
+// }
//
Int_t ntracks=esd->GetNumberOfTracks();
for (Int_t i=0;i<ntracks;++i) {
Bool_t isOK=kFALSE;
AliESDtrack *track = esd->GetTrack(i);
if (track->GetTPCncls()<kMinCluster) continue;
- AliESDfriendTrack *friendTrack = esdFriend->GetTrack(i);
- if (!friendTrack) continue;
if (TMath::Abs(AliTracker::GetBz())>1){ // cut on momenta if measured
if (track->Pt()>kMinPt) isOK=kTRUE;
}
if (TMath::Abs(dvertex[0]/TMath::Sqrt(cvertex[0]+0.01))>10) isAccepted=kFALSE;
if (TMath::Abs(dvertex[1]/TMath::Sqrt(TMath::Abs(cvertex[2]+0.01)))>10) isAccepted=kFALSE;
if (!isAccepted) isOK=kFALSE;
+ }
+ if ( track->GetTPCsignal()>100 && track->GetInnerParam()->Pt()>1 ){
+ if (track->IsOn(AliESDtrack::kITSin)||track->IsOn(AliESDtrack::kTRDout)||track->IsOn(AliESDtrack::kTOFin))
+ isOK=kTRUE;
+ if (isOK){
+ TString filename(AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile()->GetName());
+ Int_t eventNumber = esd->GetEventNumberInFile();
+ Bool_t hasFriend=(esdFriend) ? (esdFriend->GetTrack(i)!=0):0;
+ Bool_t hasITS=(track->GetNcls(0)>2);
+ printf("DUMPIONTrack:%s|%f|%d|%d|%d\n",filename.Data(),track->GetInnerParam()->Pt()*track->GetTPCsignal()/50., eventNumber,hasFriend,hasITS);
+ }
}
-
+ if (!isOK) continue;
+ TString filename(AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile()->GetName());
+ Int_t eventNumber = esd->GetEventNumberInFile();
+ Bool_t hasFriend=(esdFriend) ? (esdFriend->GetTrack(i)!=0):0;
+ Bool_t hasITS=(track->GetNcls(0)>2);
+ printf("DUMPHPTTrack:%s|%f|%d|%d|%d\n",filename.Data(),track->Pt(), eventNumber,hasFriend,hasITS);
+ //
+ if (!esdFriend) continue;
+ AliESDfriendTrack *friendTrack = esdFriend->GetTrack(i);
+ if (!friendTrack) continue;
+
if (!isOK) continue;
//
Int_t nV0s = fESD->GetNumberOfV0s();
const Int_t kMinCluster=110;
const Double_t kDownscale=0.01;
- const Float_t kMinR =0;
const Float_t kMinPt =1.0;
const Float_t kMinMinPt =0.7;
AliESDfriend *esdFriend=static_cast<AliESDfriend*>(esd->FindListObject("AliESDfriend"));
- if (!esdFriend) {
- Printf("ERROR: esdFriend not available");
- return;
- }
//
for (Int_t ivertex=0; ivertex<nV0s; ivertex++){
if (TMath::Max(track0->Pt(),track1->Pt())>kMinPt) isOK=kTRUE;
if (gRandom->Rndm()<kDownscale) isOK=kTRUE;
if (!isOK) continue;
+ //
+ TString filename(AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile()->GetName());
+ Int_t eventNumber = esd->GetEventNumberInFile();
+ Bool_t hasITS=(track0->GetNcls(0)+ track1->GetNcls(0)>4);
+ printf("DUMPHPTV0:%s|%f|%d|%d|%d\n",filename.Data(), (TMath::Min(track0->Pt(),track1->Pt())), eventNumber,(esdFriend!=0), hasITS);
+ //
+ if (!esdFriend) continue;
+ //
+
//
AliESDfriendTrack *ftrack0 = esdFriend->GetTrack(v0->GetIndex(0));
if (!ftrack0) continue;
// Add the content of tree:
// Notice automatic copy of tree in ROOT does not work for such complicated tree
//
+ return ;
AliESDv0 * v0 = new AliESDv0;
Double_t kMinPt=0.8;
AliESDtrack * track0 = 0; // negative track
}
//
treeInput->GetEntry(i);
- ftrack0->GetCalibContainer()->SetOwner(kTRUE);
- ftrack1->GetCalibContainer()->SetOwner(kTRUE);
+ //ftrack0->GetCalibContainer()->SetOwner(kTRUE);
+ //ftrack1->GetCalibContainer()->SetOwner(kTRUE);
Bool_t isOK=kTRUE;
if (v0->GetOnFlyStatus()==kFALSE) isOK=kFALSE;
if (track0->GetTPCncls()<100) isOK=kFALSE;
// Add the content of tree:
// Notice automatic copy of tree in ROOT does not work for such complicated tree
//
+ return ;
AliESDtrack *track = 0;
AliESDfriendTrack *friendTrack = 0;
AliTPCseed *seed = 0;
}
-void AliTPCcalibV0::MakeMC(){
- //
- // MC comparison
- // 1. Select interesting particles
- // 2. Assign the recosntructed particles
- //
- //1. Select interesting particles
- const Float_t kMinP = 0.2;
- const Float_t kMinPt = 0.1;
- const Float_t kMaxR = 0.5;
- const Float_t kMaxTan = 1.2;
- const Float_t kMaxRad = 150;
- //
- if (!fParticles) fParticles = new TObjArray;
- TParticle *part=0;
- //
- Int_t entries = fStack->GetNtrack();
- for (Int_t ipart=0; ipart<entries; ipart++){
- part = fStack->Particle(ipart);
- if (!part) continue;
- if (part->P()<kMinP) continue;
- if (part->R()>kMaxR) continue;
- if (TMath::Abs(TMath::Tan(part->Theta()-TMath::Pi()*0.5))>kMaxTan) continue;
- Bool_t isInteresting = kFALSE;
- if (part->GetPdgCode()==22) isInteresting =kTRUE;
- if (part->GetPdgCode()==310) isInteresting =kTRUE;
- if (part->GetPdgCode()==111) isInteresting =kTRUE;
- if (TMath::Abs(part->GetPdgCode()==3122)) isInteresting =kTRUE;
-
- //
- if (!isInteresting) continue;
- fParticles->AddLast(new TParticle(*part));
- }
- if (fParticles->GetEntries()<1) {
- return;
- }
- //
- //
- //
- Int_t sentries=fParticles->GetEntries();;
- for (Int_t ipart=0; ipart<sentries; ipart++){
- part = (TParticle*)fParticles->At(ipart);
- TParticle *p0 = 0;
- TParticle *p1 = 0;
-
- Int_t nold =0;
- Int_t nnew =0;
- Int_t id0 = part->GetDaughter(0);
- Int_t id1 = part->GetDaughter(1);
- if (id0>=fStack->GetNtrack() ) id0*=-1;
- if (id1>=fStack->GetNtrack() ) id1*=-1;
- Bool_t findable = kTRUE;
- if (id0<0 || id1<0) findable = kFALSE;
- Int_t charge =0;
- if (findable){
- p0 = fStack->Particle(id0);
- if (p0->R()>kMaxRad) findable = kFALSE;
- if (p0->Pt()<kMinPt) findable = kFALSE;
- if (p0->Vz()>250) findable= kFALSE;
- if (TMath::Abs(TMath::Tan(p0->Theta()-TMath::Pi()*0.5))>2) findable=kFALSE;
- if (fPdg->GetParticle(p0->GetPdgCode())==0) findable =kFALSE;
- else
- if (fPdg->GetParticle(p0->GetPdgCode())->Charge()==0) charge++;
-
- p1 = fStack->Particle(id1);
- if (p1->R()>kMaxRad) findable = kFALSE;
- if (p1->Pt()<kMinPt) findable = kFALSE;
- if (TMath::Abs(p1->Vz())>250) findable= kFALSE;
- if (TMath::Abs(TMath::Tan(p1->Theta()-TMath::Pi()*0.5))>2) findable=kFALSE;
- if (fPdg->GetParticle(p1->GetPdgCode())==0) findable = kFALSE;
- else
- if (fPdg->GetParticle(p1->GetPdgCode())->Charge()==0) charge++;
-
- }
- // (*fDebugStream)<<"MC0"<<
- // "P.="<<part<<
- // "findable="<<findable<<
- // "id0="<<id0<<
- // "id1="<<id1<<
- // "\n";
- if (!findable) continue;
- Float_t minpt = TMath::Min(p0->Pt(), p1->Pt());
- Int_t type=-1;
-
- //
- //
- AliKFVertex primVtx(*(fESD->GetPrimaryVertex()));
- for (Int_t ivertex=0; ivertex<fESD->GetNumberOfV0s(); ivertex++){
- AliESDv0 * v0 = fESD->GetV0(ivertex);
- // select coresponding track
- AliESDtrack * trackN = fESD->GetTrack(v0->GetIndex(0));
- if (TMath::Abs(trackN->GetLabel())!=id0 && TMath::Abs(trackN->GetLabel())!=id1) continue;
- AliESDtrack * trackP = fESD->GetTrack(v0->GetIndex(1));
- if (TMath::Abs(trackP->GetLabel())!=id0 && TMath::Abs(trackP->GetLabel())!=id1) continue;
- TParticle *pn = fStack->Particle(TMath::Abs(trackN->GetLabel()));
- TParticle *pp = fStack->Particle(TMath::Abs(trackP->GetLabel()));
- //
- //
- if ( v0->GetOnFlyStatus()) nnew++;
- if (!v0->GetOnFlyStatus()) nold++;
- if (part->GetPdgCode()==22 && TMath::Abs(pn->GetPdgCode())==11 && TMath::Abs(pp->GetPdgCode())==11)
- type =1;
- if (part->GetPdgCode()==310 && TMath::Abs(pn->GetPdgCode())==211 && TMath::Abs(pp->GetPdgCode())==211)
- type =0;
- if (part->GetPdgCode()==3122){
- if (TMath::Abs(pn->GetPdgCode())==210 ) type=2;
- else type=3;
- }
- AliKFParticle *v0kf = Fit(primVtx,v0,pn->GetPdgCode(),pp->GetPdgCode());
- v0kf->SetProductionVertex( primVtx );
- AliKFParticle *v0kfc = Fit(primVtx,v0,pn->GetPdgCode(),pp->GetPdgCode());
- v0kfc->SetProductionVertex( primVtx );
- v0kfc->SetMassConstraint(fPdg->GetParticle(part->GetPdgCode())->Mass());
- Float_t chi2 = v0kf->GetChi2();
- Float_t chi2C = v0kf->GetChi2();
- //
- //
- TTreeSRedirector *cstream = GetDebugStreamer();
- if (cstream){
- (*cstream)<<"MCRC"<<
- "P.="<<part<<
- "type="<<type<<
- "chi2="<<chi2<<
- "chi2C="<<chi2C<<
- "minpt="<<minpt<<
- "id0="<<id0<<
- "id1="<<id1<<
- "Pn.="<<pn<<
- "Pp.="<<pp<<
- "tn.="<<trackN<<
- "tp.="<<trackP<<
- "nold.="<<nold<<
- "nnew.="<<nnew<<
- "v0.="<<v0<<
- "v0kf.="<<v0kf<<
- "v0kfc.="<<v0kfc<<
- "\n";
- delete v0kf;
- delete v0kfc;
- //
- }
-
- if (cstream){
- (*cstream)<<"MC"<<
- "P.="<<part<<
- "charge="<<charge<<
- "type="<<type<<
- "minpt="<<minpt<<
- "id0="<<id0<<
- "id1="<<id1<<
- "P0.="<<p0<<
- "P1.="<<p1<<
- "nold="<<nold<<
- "nnew="<<nnew<<
- "\n";
- }
- }
- fParticles->Delete();
- }
-}
-
-void AliTPCcalibV0::MakeFitTreeTrack(const TObjArray * corrArray, Double_t ptCut, Int_t run){
+void AliTPCcalibV0::MakeFitTreeTrack(const TObjArray * corrArray, Double_t ptCut, Int_t /*run*/){
//
// Make a fit tree
//
Int_t ncorr=0;
if (corrArray) ncorr = corrArray->GetEntries();
AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ;
- AliTPCParam *param = AliTPCcalibDB::Instance()->GetParameters();
- AliGRPObject* grp = AliTPCcalibDB::Instance()->GetGRP(run);
- Double_t time=0.5*(grp->GetTimeStart() +grp->GetTimeEnd());
+ // AliTPCParam *param = AliTPCcalibDB::Instance()->GetParameters();
+// AliGRPObject* grp = AliTPCcalibDB::Instance()->GetGRP(run);
+// Double_t time=0.5*(grp->GetTimeStart() +grp->GetTimeEnd());
//
//
//
}
}
//
- Double_t xref=134;
AliExternalTrackParam* paramInner=0;
AliExternalTrackParam* paramOuter=0;
AliExternalTrackParam* paramIO=0;
//
//Connect input
- const Int_t kMinNcl=120;
TFile f("TPCV0Objects.root");
AliTPCcalibV0 *v0TPC = (AliTPCcalibV0*) f.Get("v0TPC");
TTree * treeInput = v0TPC->GetV0Tree();
Int_t ncorr=0;
if (corrArray) ncorr = corrArray->GetEntries();
AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ;
- AliTPCParam *param = AliTPCcalibDB::Instance()->GetParameters();
- AliGRPObject* grp = AliTPCcalibDB::Instance()->GetGRP(run);
- Double_t time=0.5*(grp->GetTimeStart() +grp->GetTimeEnd());
Double_t massK0= pdg.GetParticle("K0")->Mass();
Double_t massLambda= pdg.GetParticle("Lambda0")->Mass();
Double_t massPion=pdg.GetParticle("pi+")->Mass();
Double_t massProton=pdg.GetParticle("proton")->Mass();
- Double_t pdgPion=pdg.GetParticle("pi+")->PdgCode();
- Double_t pdgProton=pdg.GetParticle("proton")->PdgCode();
- Double_t mass0=0;
- Double_t mass1=0;
+ Int_t pdgPion=pdg.GetParticle("pi+")->PdgCode();
+ Int_t pdgProton=pdg.GetParticle("proton")->PdgCode();
+ Double_t rmass0=0;
+ Double_t rmass1=0;
Double_t massV0=0;
Int_t pdg0=0;
Int_t pdg1=0;
if (TMath::Abs(v0->GetEffMass(4,2)-massLambda)<0.01) {isLambda=1; v0Type=2;} //select Lambda
if (TMath::Abs(v0->GetEffMass(2,4)-massLambda)<0.01) {isAntiLambda=1;v0Type=3;} //select Anti Lambda
if (isK0+isLambda+isAntiLambda!=1) continue;
- mass0=massPion;
- mass1=massPion;
+ rmass0=massPion;
+ rmass1=massPion;
pdg0=pdgPion;
pdg1=pdgPion;
- if (isLambda) {mass0=massProton; pdg0=pdgProton;}
- if (isAntiLambda) {mass1=massProton; pdg1=pdgProton;}
+ if (isLambda) {rmass0=massProton; pdg0=pdgProton;}
+ if (isAntiLambda) {rmass1=massProton; pdg1=pdgProton;}
massV0=massK0;
if (isK0==0) massV0=massLambda;
//
AliTPCCorrection *corr =0;
if (icorr>=0) corr = (AliTPCCorrection*)corrArray->At(icorr);
//
- AliExternalTrackParam * trackInner0 = RefitTrack(seed0, corr,160,85,mass0);
- AliExternalTrackParam * trackIO0 = RefitTrack(seed0, corr,245,85,mass0);
- AliExternalTrackParam * trackInner1 = RefitTrack(seed1, corr,160,85,mass1);
- AliExternalTrackParam * trackIO1 = RefitTrack(seed1, corr,245,85,mass1);
+ AliExternalTrackParam * trackInner0 = RefitTrack(seed0, corr,160,85,rmass0);
+ AliExternalTrackParam * trackIO0 = RefitTrack(seed0, corr,245,85,rmass0);
+ AliExternalTrackParam * trackInner1 = RefitTrack(seed1, corr,160,85,rmass1);
+ AliExternalTrackParam * trackIO1 = RefitTrack(seed1, corr,245,85,rmass1);
if (!trackInner0) isOK=kFALSE;
if (!trackInner1) isOK=kFALSE;
if (!trackIO0) isOK=kFALSE;
if (!trackIO0->Rotate(alpha)) isOK=kFALSE;
if (!trackIO1->Rotate(alpha)) isOK=kFALSE;
//
- if (!AliTracker::PropagateTrackToBxByBz(trackInner0, radius, mass0, 1, kFALSE)) isOK=kFALSE;
- if (!AliTracker::PropagateTrackToBxByBz(trackInner1, radius, mass1, 1, kFALSE)) isOK=kFALSE;
- if (!AliTracker::PropagateTrackToBxByBz(trackIO0, radius, mass0, 1, kFALSE)) isOK=kFALSE;
- if (!AliTracker::PropagateTrackToBxByBz(trackIO1, radius, mass1, 1, kFALSE)) isOK=kFALSE;
+ if (!AliTracker::PropagateTrackToBxByBz(trackInner0, radius, rmass0, 1, kFALSE)) isOK=kFALSE;
+ if (!AliTracker::PropagateTrackToBxByBz(trackInner1, radius, rmass1, 1, kFALSE)) isOK=kFALSE;
+ if (!AliTracker::PropagateTrackToBxByBz(trackIO0, radius, rmass0, 1, kFALSE)) isOK=kFALSE;
+ if (!AliTracker::PropagateTrackToBxByBz(trackIO1, radius, rmass1, 1, kFALSE)) isOK=kFALSE;
if (!isOK) continue;
arrayT0.AddAt(trackIO0->Clone(),icorr+1);
arrayT1.AddAt(trackIO1->Clone(),icorr+1);
- Int_t charge=(trackIO0->GetSign());
+ Int_t charge=TMath::Nint(trackIO0->GetSign());
AliKFParticle pin0( *trackInner0, pdg0*charge);
AliKFParticle pin1( *trackInner1, -pdg1*charge);
AliKFParticle pio0( *trackIO0, pdg0*charge);
Int_t dtype=30; // id for V0
Int_t ptype=5; // id for invariant mass
// Int_t id=TMath::Nint(100.*(param0->Pt()-param1->Pt())/(param0->Pt()+param1->Pt())); // K0s V0 asymetry
- Int_t id=1000.*(param0->Pt()-param1->Pt()); // K0s V0 asymetry
+ Int_t id=Int_t(1000.*(param0->Pt()-param1->Pt())); // K0s V0 asymetry
Double_t gx,gy,gz, px,py,pz;
Double_t pt = v0->Pt();
v0->GetXYZ(gx,gy,gz);
for (Int_t icorr=0; icorr<ncorr; icorr++){
AliTPCCorrection *corr =0;
if (icorr>=0) corr = (AliTPCCorrection*)corrArray->At(icorr);
- AliKFParticle* pin= (AliKFParticle*)arrayV0in.At(icorr+1);
+ // AliKFParticle* pin= (AliKFParticle*)arrayV0in.At(icorr+1);
AliKFParticle* pio= (AliKFParticle*)arrayV0io.At(icorr+1);
AliExternalTrackParam *par0=(AliExternalTrackParam *)arrayT0.At(icorr+1);
AliExternalTrackParam *par1=(AliExternalTrackParam *)arrayT1.At(icorr+1);
//
-void AliTPCcalibV0::MakeV0s(){
- //
- //
- //
- const Int_t kMinCluster=40;
- const Float_t kMinR =0;
- if (! fV0s) fV0s = new TObjArray(10);
- fV0s->Clear();
- //
- // Old V0 finder
- //
- for (Int_t ivertex=0; ivertex<fESD->GetNumberOfV0s(); ivertex++){
- AliESDv0 * v0 = fESD->GetV0(ivertex);
- if (v0->GetOnFlyStatus()) continue;
- fV0s->AddLast(v0);
- }
- ProcessV0(0);
- fV0s->Clear(0);
- //
- // MI V0 finder
- //
- for (Int_t ivertex=0; ivertex<fESD->GetNumberOfV0s(); ivertex++){
- AliESDv0 * v0 = fESD->GetV0(ivertex);
- if (!v0->GetOnFlyStatus()) continue;
- fV0s->AddLast(v0);
- }
- ProcessV0(1);
- fV0s->Clear();
- return;
- //
- // combinatorial
- //
- Int_t ntracks = fESD->GetNumberOfTracks();
- for (Int_t itrack0=0; itrack0<ntracks; itrack0++){
- AliESDtrack * track0 = fESD->GetTrack(itrack0);
- if (track0->GetSign()>0) continue;
- if ( track0->GetTPCNcls()<kMinCluster) continue;
- if (track0->GetKinkIndex(0)>0) continue;
- //
- for (Int_t itrack1=0; itrack1<ntracks; itrack1++){
- AliESDtrack * track1 = fESD->GetTrack(itrack1);
- if (track1->GetSign()<0) continue;
- if ( track1->GetTPCNcls()<kMinCluster) continue;
- if (track1->GetKinkIndex(0)>0) continue;
- //
- // AliExternalTrackParam param0(*track0);
- // AliExternalTrackParam param1(*track1);
- AliV0 vertex;
- vertex.SetParamN(*track0);
- vertex.SetParamP(*track1);
- Float_t xyz[3];
- xyz[0] = fESD->GetPrimaryVertex()->GetXv();
- xyz[1] = fESD->GetPrimaryVertex()->GetYv();
- xyz[2] = fESD->GetPrimaryVertex()->GetZv();
- vertex.Update(xyz);
- if (vertex.GetRr()<kMinR) continue;
- if (vertex.GetDcaV0Daughters()>1.) continue;
- if (vertex.GetDcaV0Daughters()>0.3*vertex.GetRr()) continue;
- // if (vertex.GetPointAngle()<0.9) continue;
- vertex.SetIndex(0,itrack0);
- vertex.SetIndex(1,itrack1);
- fV0s->AddLast(new AliV0(vertex));
- }
- }
- ProcessV0(2);
- for (Int_t i=0;i<fV0s->GetEntries(); i++) delete fV0s->At(i);
- fV0s->Clear();
-}
-
-
-
-
-
-
-
-
-void AliTPCcalibV0::ProcessV0(Int_t ftype){
- //
- // Obsolete
- //
- if (! fGammas) fGammas = new TObjArray(10);
- fGammas->Clear();
- Int_t nV0s = fV0s->GetEntries();
- if (nV0s==0) return;
- AliKFVertex primVtx(*(fESD->GetPrimaryVertex()));
- //
- for (Int_t ivertex=0; ivertex<nV0s; ivertex++){
- AliESDv0 * v0 = (AliESDv0*)fV0s->At(ivertex);
- AliESDtrack * trackN = fESD->GetTrack(v0->GetIndex(0)); // negative track
- AliESDtrack * trackP = fESD->GetTrack(v0->GetIndex(1)); // positive track
-
- const AliExternalTrackParam * paramInNeg = trackN->GetInnerParam();
- const AliExternalTrackParam * paramInPos = trackP->GetInnerParam();
-
- if (!paramInPos) continue; // in case the inner paramters do not exist
- if (!paramInNeg) continue;
- //
- //
- //
- AliKFParticle *v0K0 = Fit(primVtx,v0,-211,211);
- AliKFParticle *v0Gamma = Fit(primVtx,v0,11,-11);
- AliKFParticle *v0Lambda42 = Fit(primVtx,v0,-2212,211);
- AliKFParticle *v0Lambda24 = Fit(primVtx,v0,-211,2212);
- //Set production vertex
- v0K0->SetProductionVertex( primVtx );
- v0Gamma->SetProductionVertex( primVtx );
- v0Lambda42->SetProductionVertex( primVtx );
- v0Lambda24->SetProductionVertex( primVtx );
- Double_t massK0, massGamma, massLambda42,massLambda24, massSigma;
- v0K0->GetMass( massK0,massSigma);
- v0Gamma->GetMass( massGamma,massSigma);
- v0Lambda42->GetMass( massLambda42,massSigma);
- v0Lambda24->GetMass( massLambda24,massSigma);
- Float_t chi2K0 = v0K0->GetChi2()/v0K0->GetNDF();
- Float_t chi2Gamma = v0Gamma->GetChi2()/v0Gamma->GetNDF();
- Float_t chi2Lambda42 = v0Lambda42->GetChi2()/v0Lambda42->GetNDF();
- Float_t chi2Lambda24 = v0Lambda24->GetChi2()/v0Lambda24->GetNDF();
- //
- // Mass Contrained params
- //
- AliKFParticle *v0K0C = Fit(primVtx,v0,-211,211);
- AliKFParticle *v0GammaC = Fit(primVtx,v0,11,-11);
- AliKFParticle *v0Lambda42C = Fit(primVtx,v0,-2212,211); //lambdaBar
- AliKFParticle *v0Lambda24C = Fit(primVtx,v0,-211,2212); //lambda
- //
- v0K0C->SetProductionVertex( primVtx );
- v0GammaC->SetProductionVertex( primVtx );
- v0Lambda42C->SetProductionVertex( primVtx );
- v0Lambda24C->SetProductionVertex( primVtx );
-
- v0K0C->SetMassConstraint(fPdg->GetParticle(310)->Mass());
- v0GammaC->SetMassConstraint(0);
- v0Lambda42C->SetMassConstraint(fPdg->GetParticle(-3122)->Mass());
- v0Lambda24C->SetMassConstraint(fPdg->GetParticle(3122)->Mass());
- //
- Double_t timeK0, sigmaTimeK0;
- Double_t timeLambda42, sigmaTimeLambda42;
- Double_t timeLambda24, sigmaTimeLambda24;
- v0K0C->GetLifeTime(timeK0, sigmaTimeK0);
- //v0K0Gamma->GetLifeTime(timeK0, sigmaTimeK0);
- v0Lambda42C->GetLifeTime(timeLambda42, sigmaTimeLambda42);
- v0Lambda24C->GetLifeTime(timeLambda24, sigmaTimeLambda24);
-
-
- //
- Float_t chi2K0C = v0K0C->GetChi2()/v0K0C->GetNDF();
- if (chi2K0C<0) chi2K0C=100;
- Float_t chi2GammaC = v0GammaC->GetChi2()/v0GammaC->GetNDF();
- if (chi2GammaC<0) chi2GammaC=100;
- Float_t chi2Lambda42C = v0Lambda42C->GetChi2()/v0Lambda42C->GetNDF();
- if (chi2Lambda42C<0) chi2Lambda42C=100;
- Float_t chi2Lambda24C = v0Lambda24C->GetChi2()/v0Lambda24C->GetNDF();
- if (chi2Lambda24C<0) chi2Lambda24C=100;
- //
- Float_t minChi2C=99;
- Int_t type =-1;
- if (chi2K0C<minChi2C) { minChi2C= chi2K0C; type=0;}
- if (chi2GammaC<minChi2C) { minChi2C= chi2GammaC; type=1;}
- if (chi2Lambda42C<minChi2C) { minChi2C= chi2Lambda42C; type=2;}
- if (chi2Lambda24C<minChi2C) { minChi2C= chi2Lambda24C; type=3;}
- Float_t minChi2=99;
- Int_t type0 =-1;
- if (chi2K0<minChi2) { minChi2= chi2K0; type0=0;}
- if (chi2Gamma<minChi2) { minChi2= chi2Gamma; type0=1;}
- if (chi2Lambda42<minChi2) { minChi2= chi2Lambda42; type0=2;}
- if (chi2Lambda24<minChi2) { minChi2= chi2Lambda24; type0=3;}
-
- // 0 is negative particle; 1 is positive particle
- Float_t betaGamma0 = 0;
- Float_t betaGamma1 = 0;
-
- switch (type) {
- case 0:
- betaGamma0 = paramInNeg->GetP()/fPdg->GetParticle(-211)->Mass();
- betaGamma1 = paramInPos->GetP()/fPdg->GetParticle(211)->Mass();
- break;
- case 1:
- betaGamma0 = paramInNeg->GetP()/fPdg->GetParticle(11)->Mass();
- betaGamma1 = paramInPos->GetP()/fPdg->GetParticle(-11)->Mass();
- break;
- case 2:
- betaGamma0 = paramInNeg->GetP()/fPdg->GetParticle(-2212)->Mass();
- betaGamma1 = paramInPos->GetP()/fPdg->GetParticle(211)->Mass();
- break;
- case 3:
- betaGamma0 = paramInNeg->GetP()/fPdg->GetParticle(-211)->Mass();
- betaGamma1 = paramInPos->GetP()/fPdg->GetParticle(2212)->Mass();
- break;
- }
-
- // cuts and histogram filling
- Int_t numCand = 0; // number of particle types which have a chi2 < 10*minChi2
-
- if (minChi2C < 2 && ftype == 1) {
- //
- if (chi2K0C < 10*minChi2C) numCand++;
- if (chi2GammaC < 10*minChi2C) numCand++;
- if (chi2Lambda42C < 10*minChi2C) numCand++;
- if (chi2Lambda24C < 10*minChi2C) numCand++;
- //
- }
-
- //
- //
- // write output tree
- if (minChi2>50) continue;
- TTreeSRedirector *cstream = GetDebugStreamer();
- if (cstream){
- (*cstream)<<"V0"<<
- "ftype="<<ftype<<
- "v0.="<<v0<<
- "trackN.="<<trackN<<
- "trackP.="<<trackP<<
- //
- "betaGamma0="<<betaGamma0<<
- "betaGamma1="<<betaGamma1<<
- //
- "type="<<type<<
- "chi2C="<<minChi2C<<
- "v0K0.="<<v0K0<<
- "v0Gamma.="<<v0Gamma<<
- "v0Lambda42.="<<v0Lambda42<<
- "v0Lambda24.="<<v0Lambda24<<
- //
- "chi20K0.="<<chi2K0<<
- "chi2Gamma.="<<chi2Gamma<<
- "chi2Lambda42.="<<chi2Lambda42<<
- "chi2Lambda24.="<<chi2Lambda24<<
- //
- "chi20K0c.="<<chi2K0C<<
- "chi2Gammac.="<<chi2GammaC<<
- "chi2Lambda42c.="<<chi2Lambda42C<<
- "chi2Lambda24c.="<<chi2Lambda24C<<
- //
- "v0K0C.="<<v0K0C<<
- "v0GammaC.="<<v0GammaC<<
- "v0Lambda42C.="<<v0Lambda42C<<
- "v0Lambda24C.="<<v0Lambda24C<<
- //
- "massK0="<<massK0<<
- "massGamma="<<massGamma<<
- "massLambda42="<<massLambda42<<
- "massLambda24="<<massLambda24<<
- //
- "timeK0="<<timeK0<<
- "timeLambda42="<<timeLambda42<<
- "timeLambda24="<<timeLambda24<<
- "\n";
- }
- if (type==1) fGammas->AddLast(v0);
- //
- //
- //
- delete v0K0;
- delete v0Gamma;
- delete v0Lambda42;
- delete v0Lambda24;
- delete v0K0C;
- delete v0GammaC;
- delete v0Lambda42C;
- delete v0Lambda24C;
- }
- ProcessPI0();
-}
-
-
-
-void AliTPCcalibV0::ProcessPI0(){
- //
- //
- //
- Int_t nentries = fGammas->GetEntries();
- if (nentries<2) return;
- //
- Double_t m0[3], m1[3];
- AliKFVertex primVtx(*(fESD->GetPrimaryVertex()));
- for (Int_t i0=0; i0<nentries; i0++){
- AliESDv0 *v00 = (AliESDv0*)fGammas->At(i0);
- v00->GetPxPyPz (m0[0], m0[1], m0[2]);
- AliKFParticle *p00 = Fit(primVtx, v00, 11,-11);
- p00->SetProductionVertex( primVtx );
- p00->SetMassConstraint(0);
- //
- for (Int_t i1=i0; i1<nentries; i1++){
- AliESDv0 *v01 = (AliESDv0*)fGammas->At(i1);
- v01->GetPxPyPz (m1[0], m1[1], m1[2]);
- AliKFParticle *p01 = Fit(primVtx, v01, 11,-11);
- p01->SetProductionVertex( primVtx );
- p01->SetMassConstraint(0);
- if (v00->GetIndex(0) != v01->GetIndex(0) &&
- v00->GetIndex(1) != v01->GetIndex(1)){
- AliKFParticle pi0( *p00,*p01);
- pi0.SetProductionVertex(primVtx);
- Double_t n1 = TMath::Sqrt (m0[0]*m0[0] + m0[1]*m0[1] + m0[2]*m0[2]);
- Double_t n2 = TMath::Sqrt (m1[0]*m1[0] + m1[1]*m1[1] + m1[2]*m1[2]);
- Double_t mass = TMath::Sqrt(2.*(n1*n2 - (m0[0]*m1[0] + m0[1]*m1[1] + m0[2]*m1[2])));
- TTreeSRedirector *cstream = GetDebugStreamer();
- if (cstream){
- (*cstream)<<"PI0"<<
- "v00.="<<v00<<
- "v01.="<<v01<<
- "mass="<<mass<<
- "p00.="<<p00<<
- "p01.="<<p01<<
- "pi0.="<<&pi0<<
- "\n";
- }
- }
- delete p01;
- }
- delete p00;
- }
-}
-
-
-
AliKFParticle * AliTPCcalibV0::Fit(AliKFVertex & /*primVtx*/, AliESDv0 *v0, Int_t PDG1, Int_t PDG2){
new_bins[i] = factor * new_bins[i-1];
}
axis->Set(bins, new_bins);
- delete [] new_bins;
-
+ delete [] new_bins;
}
-
+void AliTPCcalibV0::FilterV0s(AliESDEvent* event){
+ //
+ //
+ TDatabasePDG pdg;
+ const Double_t kChi2Cut=20;
+ const Double_t kMinR=2;
+ const Double_t ptCut=0.2;
+ const Int_t kMinNcl=110;
+ //
+ Int_t nv0 = event->GetNumberOfV0s();
+ AliESDVertex *vertex= (AliESDVertex *)event->GetPrimaryVertex();
+ AliKFVertex kfvertex=*vertex;
+ //
+ for (Int_t iv0=0;iv0<nv0;iv0++){
+ AliESDv0 *v0 = event->GetV0(iv0);
+ if (!v0) continue;
+ if (v0->GetPindex()<0) continue;
+ if (v0->GetNindex()<0) continue;
+ if (TMath::Max(v0->GetPindex(), v0->GetNindex())>event->GetNumberOfTracks()) continue;
+ //
+ //
+ AliExternalTrackParam pp=(v0->GetParamP()->GetSign()>0) ? (*(v0->GetParamP())):(*(v0->GetParamN()));
+ AliExternalTrackParam pn=(v0->GetParamP()->GetSign()>0) ? (*(v0->GetParamN())):(*(v0->GetParamP()));
+ AliKFParticle kfp1( pp, 211 );
+ AliKFParticle kfp2( pn, -211 );
+ AliKFParticle *v0KFK0 = new AliKFParticle(kfp1,kfp2);
+ AliKFParticle *v0KFK0CV = new AliKFParticle(*v0KFK0);
+ v0KFK0CV->SetProductionVertex(kfvertex);
+ v0KFK0CV->TransportToProductionVertex();
+ AliKFParticle *v0KFK0CVM = new AliKFParticle(*v0KFK0CV);
+ v0KFK0CVM->SetMassConstraint(pdg.GetParticle("K_S0")->Mass());
+ Double_t chi2K0 = v0KFK0CV->GetChi2();
+ Double_t chi2K0M= v0KFK0CVM->GetChi2();
+ Double_t massK0=v0KFK0CV->GetMass();
+ if (chi2K0>kChi2Cut) continue;
+ if (v0->GetRr()>kMinR) continue;
+ //
+ Double_t maxPt = TMath::Max(v0->GetParamP()->Pt(), v0->GetParamN()->Pt());
+ Double_t effMass22=v0->GetEffMass(2,2);
+ Double_t effMass42=v0->GetEffMass(4,2);
+ Double_t effMass24=v0->GetEffMass(2,4);
+ Bool_t isV0= kFALSE;
+ isV0|=TMath::Abs(effMass22-pdg.GetParticle("K_S0")->Mass())<0.1;
+ isV0|=TMath::Abs(effMass42-pdg.GetParticle("Lambda0")->Mass())<0.1;
+ isV0|=TMath::Abs(effMass24-pdg.GetParticle("Lambda0")->Mass())<0.1;
+
+ Double_t sign= v0->GetParamP()->GetSign()* v0->GetParamN()->GetSign();
+ if (sign<0&&v0->GetOnFlyStatus()>0.5&&maxPt>ptCut&&isV0){
+ AliESDtrack * trackP = event->GetTrack(v0->GetPindex());
+ AliESDtrack * trackN = event->GetTrack(v0->GetNindex());
+ if (!trackN) continue;
+ if (!trackP) continue;
+ Int_t nclP= (Int_t)trackP->GetTPCClusterInfo(2,1);
+ Int_t nclN= (Int_t)trackN->GetTPCClusterInfo(2,1);
+ if (TMath::Min(nclP,nclN)<kMinNcl) continue;
+ Double_t eta = TMath::Max(TMath::Abs(trackP->Eta()), TMath::Abs(trackN->Eta()));
+ Double_t ncls = TMath::Min(TMath::Abs(trackP->GetNcls(0)), TMath::Abs(trackN->GetNcls(0)));
+ if (eta<0.8&&ncls>2){
+ // printf("%d\t%f\t%f\t%d\t%d\t%f\t%f\n",i, v0->Pt(), maxPt, v0->GetNindex(),v0->GetPindex(),v0->GetRr(),effMass22);
+ (*fDebugStreamer)<<"v0tree"<<
+ "v0.="<<v0<<
+ "tp.="<<trackP<<
+ "tm.="<<trackN<<
+ //
+ "v.="<<vertex<<
+ "ncls="<<ncls<<
+ "maxPt="<<maxPt<<
+ "\n";
+ }
+ }
+ }
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff