X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PHOS%2FAliPHOSPIDv1.cxx;h=3df63d818cda0711ae9564de795f023d5ccae14d;hb=8ebd7df4127272d8dc8e711539dd45778cf247f6;hp=ba0ee55bf087a1f2bf7172083a4df6cd7599810c;hpb=999f9a8fac5455f15d5d2b5eef6016cb9779b9b2;p=u%2Fmrichter%2FAliRoot.git diff --git a/PHOS/AliPHOSPIDv1.cxx b/PHOS/AliPHOSPIDv1.cxx index ba0ee55bf08..3df63d818cd 100644 --- a/PHOS/AliPHOSPIDv1.cxx +++ b/PHOS/AliPHOSPIDv1.cxx @@ -18,6 +18,27 @@ /* History of cvs commits: * * $Log$ + * Revision 1.113 2007/08/07 14:12:03 kharlov + * Quality assurance added (Yves Schutz) + * + * Revision 1.112 2007/07/11 13:43:30 hristov + * New class AliESDEvent, backward compatibility with the old AliESD (Christian) + * + * Revision 1.111 2007/05/04 14:49:29 policheh + * AliPHOSRecPoint inheritance from AliCluster + * + * Revision 1.110 2007/04/24 10:08:03 kharlov + * Vertex extraction from GenHeader + * + * Revision 1.109 2007/04/18 09:34:05 kharlov + * Geometry bug fixes + * + * Revision 1.108 2007/04/16 09:03:37 kharlov + * Incedent angle correction fixed + * + * Revision 1.107 2007/04/02 15:00:16 cvetan + * No more calls to gAlice in the reconstruction + * * Revision 1.106 2007/04/01 15:40:15 kharlov * Correction for actual vertex position implemented * @@ -117,18 +138,20 @@ //#include "AliLog.h" #include "AliPHOS.h" #include "AliPHOSPIDv1.h" -#include "AliPHOSGetter.h" -#include "AliESD.h" +#include "AliESDEvent.h" #include "AliESDVertex.h" +#include "AliPHOSTrackSegment.h" +#include "AliPHOSEmcRecPoint.h" +#include "AliPHOSRecParticle.h" ClassImp( AliPHOSPIDv1) //____________________________________________________________________________ AliPHOSPIDv1::AliPHOSPIDv1() : + AliPHOSPID(), fBayesian(kFALSE), fDefaultInit(kFALSE), fWrite(kFALSE), - fNEvent(0), fFileNamePrincipalPhoton(), fFileNamePrincipalPi0(), fFileNameParameters(), @@ -137,8 +160,8 @@ AliPHOSPIDv1::AliPHOSPIDv1() : fX(0), fPPhoton(0), fPPi0(0), - fRecParticlesInRun(0), fParameters(0), + fVtx(0.), fTFphoton(0), fTFpiong(0), fTFkaong(0), @@ -164,7 +187,6 @@ AliPHOSPIDv1::AliPHOSPIDv1(const AliPHOSPIDv1 & pid ) : fBayesian(kFALSE), fDefaultInit(kFALSE), fWrite(kFALSE), - fNEvent(0), fFileNamePrincipalPhoton(), fFileNamePrincipalPi0(), fFileNameParameters(), @@ -173,8 +195,8 @@ AliPHOSPIDv1::AliPHOSPIDv1(const AliPHOSPIDv1 & pid ) : fX(0), fPPhoton(0), fPPi0(0), - fRecParticlesInRun(0), fParameters(0), + fVtx(0.), fTFphoton(0), fTFpiong(0), fTFkaong(0), @@ -191,17 +213,15 @@ AliPHOSPIDv1::AliPHOSPIDv1(const AliPHOSPIDv1 & pid ) : { // ctor InitParameters() ; - Init() ; } //____________________________________________________________________________ -AliPHOSPIDv1::AliPHOSPIDv1(const TString alirunFileName, const TString eventFolderName) : - AliPHOSPID(alirunFileName, eventFolderName), +AliPHOSPIDv1::AliPHOSPIDv1(AliPHOSGeometry *geom): + AliPHOSPID(geom), fBayesian(kFALSE), fDefaultInit(kFALSE), fWrite(kFALSE), - fNEvent(0), fFileNamePrincipalPhoton(), fFileNamePrincipalPi0(), fFileNameParameters(), @@ -210,8 +230,8 @@ AliPHOSPIDv1::AliPHOSPIDv1(const TString alirunFileName, const TString eventFold fX(0), fPPhoton(0), fPPi0(0), - fRecParticlesInRun(0), fParameters(0), + fVtx(0.), fTFphoton(0), fTFpiong(0), fTFkaong(0), @@ -229,7 +249,6 @@ AliPHOSPIDv1::AliPHOSPIDv1(const TString alirunFileName, const TString eventFold //ctor with the indication on where to look for the track segments InitParameters() ; - Init() ; fDefaultInit = kFALSE ; } @@ -253,38 +272,14 @@ AliPHOSPIDv1::~AliPHOSPIDv1() delete fTFhhadronl; delete fDFmuon; } -//____________________________________________________________________________ -const TString AliPHOSPIDv1::BranchName() const -{ - - return GetName() ; -} -//____________________________________________________________________________ -void AliPHOSPIDv1::Init() -{ - // Make all memory allocations that are not possible in default constructor - // Add the PID task to the list of PHOS tasks - - AliPHOSGetter * gime = AliPHOSGetter::Instance() ; - if(!gime) - gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName.Data()) ; - - if ( !gime->PID() ) - gime->PostPID(this) ; -} - //____________________________________________________________________________ void AliPHOSPIDv1::InitParameters() { // Initialize PID parameters fWrite = kTRUE ; - fRecParticlesInRun = 0 ; - fNEvent = 0 ; - fRecParticlesInRun = 0 ; fBayesian = kTRUE ; SetParameters() ; // fill the parameters matrix from parameters file - SetEventRange(0,-1) ; // initialisation of response function parameters // Tof @@ -485,12 +480,10 @@ void AliPHOSPIDv1::InitParameters() } //________________________________________________________________________ -void AliPHOSPIDv1::Exec(Option_t *option) +void AliPHOSPIDv1::TrackSegments2RecParticles(Option_t *option) { // Steering method to perform particle reconstruction and identification // for the event range from fFirstEvent to fLastEvent. - // This range is optionally set by SetEventRange(). - // if fLastEvent=-1 (by default), then process events until the end. if(strstr(option,"tim")) gBenchmark->Start("PHOSPID"); @@ -500,44 +493,26 @@ void AliPHOSPIDv1::Exec(Option_t *option) return ; } + if(fTrackSegments && //Skip events, where no track segments made + fTrackSegments->GetEntriesFast()) { - AliPHOSGetter * gime = AliPHOSGetter::Instance() ; - - if (fLastEvent == -1) - fLastEvent = gime->MaxEvent() - 1 ; - else - fLastEvent = TMath::Min(fLastEvent,gime->MaxEvent()); - Int_t nEvents = fLastEvent - fFirstEvent + 1; - - Int_t ievent ; - for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) { - gime->Event(ievent,"TR") ; - if(gime->TrackSegments() && //Skip events, where no track segments made - gime->TrackSegments()->GetEntriesFast()) { - - GetVertex() ; - MakeRecParticles() ; - - if(fBayesian) - MakePID() ; + GetVertex() ; + MakeRecParticles() ; + + if(fBayesian) + MakePID() ; - WriteRecParticles(); - if(strstr(option,"deb")) - PrintRecParticles(option) ; - //increment the total number of rec particles per run - fRecParticlesInRun += gime->RecParticles()->GetEntriesFast() ; - } + if(strstr(option,"deb")) + PrintRecParticles(option) ; } + if(strstr(option,"deb")) PrintRecParticles(option); if(strstr(option,"tim")){ gBenchmark->Stop("PHOSPID"); - AliInfo(Form("took %f seconds for PID %f seconds per event", - gBenchmark->GetCpuTime("PHOSPID"), - gBenchmark->GetCpuTime("PHOSPID")/nEvents)) ; + AliInfo(Form("took %f seconds for PID", + gBenchmark->GetCpuTime("PHOSPID"))); } - if(fWrite) - Unload(); } //________________________________________________________________________ @@ -770,7 +745,7 @@ Float_t AliPHOSPIDv1::GetParameterToCalculateEllipse(TString particle, TString //{ // // Calculates the distance between the EMC RecPoint and the PPSD RecPoint // -// const AliPHOSGeometry * geom = AliPHOSGetter::Instance()->PHOSGeometry() ; +// AliPHOSGeometry * geom = AliPHOSGeometry::GetInstance(); // TVector3 vecEmc ; // TVector3 vecCpv ; // if(cpv){ @@ -893,11 +868,6 @@ TVector3 AliPHOSPIDv1::GetMomentumDirection(AliPHOSEmcRecPoint * emc, AliPHOSCpv // However because of the poor position resolution of PPSD the direction is always taken as if we were // in case 1. - TVector3 dir(0,0,0) ; - TMatrixF dummy ; - - emc->GetGlobalPosition(dir, dummy) ; - TVector3 local ; emc->GetLocalPosition(local) ; @@ -917,27 +887,28 @@ TVector3 AliPHOSPIDv1::GetMomentumDirection(AliPHOSEmcRecPoint * emc, AliPHOSCpv Float_t depthzOld = 0.; Float_t energy = emc->GetEnergy() ; if (energy > 0 && vInc.Y()!=0.) { - depthxOld = ( para * TMath::Log(energy) + parb ) * vInc.X()/vInc.Y() ; - depthzOld = ( para * TMath::Log(energy) + parb ) * vInc.Z()/vInc.Y() ; + depthxOld = ( para * TMath::Log(energy) + parb ) * vInc.X()/TMath::Abs(vInc.Y()) ; + depthzOld = ( para * TMath::Log(energy) + parb ) * vInc.Z()/TMath::Abs(vInc.Y()) ; } else{ AliError("Cluster with zero energy \n"); - } + } //Apply Real vertex phosgeom->GetIncidentVector(fVtx,emc->GetPHOSMod(),x,z,vInc) ; Float_t depthx = 0.; Float_t depthz = 0.; if (energy > 0 && vInc.Y()!=0.) { - depthx = ( para * TMath::Log(energy) + parb ) * vInc.X()/vInc.Y() ; - depthz = ( para * TMath::Log(energy) + parb ) * vInc.Z()/vInc.Y() ; + depthx = ( para * TMath::Log(energy) + parb ) * vInc.X()/TMath::Abs(vInc.Y()) ; + depthz = ( para * TMath::Log(energy) + parb ) * vInc.Z()/TMath::Abs(vInc.Y()) ; } - dir.SetXYZ(dir.X()-(depthxOld-depthx)*TMath::Sin(dir.Phi()), - dir.Y()-(depthxOld-depthx)*TMath::Cos(dir.Phi()), - dir.Z()+depthzOld-depthz) ; + //Correct for the vertex position and shower depth + Double_t xd=x+(depthxOld-depthx) ; + Double_t zd=z+(depthzOld-depthz) ; + TVector3 dir(0,0,0) ; + phosgeom->Local2Global(emc->GetPHOSMod(),xd,zd,dir) ; - //Correct for the vertex position - dir = dir - fVtx ; + dir-=fVtx ; dir.SetMag(1.) ; return dir ; @@ -1023,17 +994,13 @@ void AliPHOSPIDv1::MakePID() const Int_t kSPECIES = AliPID::kSPECIESN ; - AliPHOSGetter * gime = AliPHOSGetter::Instance() ; + Int_t nparticles = fRecParticles->GetEntriesFast() ; - Int_t nparticles = gime->RecParticles()->GetEntriesFast() ; - - TObjArray * emcRecPoints = gime->EmcRecPoints() ; - TObjArray * cpvRecPoints = gime->CpvRecPoints() ; - TClonesArray * trackSegments = gime->TrackSegments() ; - if ( !emcRecPoints || !cpvRecPoints || !trackSegments ) { + if ( !fEMCRecPoints || !fCPVRecPoints || !fTrackSegments ) { AliFatal("RecPoints or TrackSegments not found !") ; } - TIter next(trackSegments) ; + + TIter next(fTrackSegments) ; AliPHOSTrackSegment * ts ; Int_t index = 0 ; @@ -1057,7 +1024,7 @@ void AliPHOSPIDv1::MakePID() AliPHOSEmcRecPoint * emc = 0 ; if(ts->GetEmcIndex()>=0) - emc = (AliPHOSEmcRecPoint *) emcRecPoints->At(ts->GetEmcIndex()) ; + emc = (AliPHOSEmcRecPoint *) fEMCRecPoints->At(ts->GetEmcIndex()) ; // AliPHOSCpvRecPoint * cpv = 0 ; // if(ts->GetCpvIndex()>=0) @@ -1305,7 +1272,7 @@ void AliPHOSPIDv1::MakePID() for(index = 0 ; index < nparticles ; index ++) { - AliPHOSRecParticle * recpar = gime->RecParticle(index) ; + AliPHOSRecParticle * recpar = static_cast(fRecParticles->At(index)); //Conversion electron? @@ -1364,34 +1331,29 @@ void AliPHOSPIDv1::MakeRecParticles() { // Makes a RecParticle out of a TrackSegment - AliPHOSGetter * gime = AliPHOSGetter::Instance() ; - TObjArray * emcRecPoints = gime->EmcRecPoints() ; - TObjArray * cpvRecPoints = gime->CpvRecPoints() ; - TClonesArray * trackSegments = gime->TrackSegments() ; - if ( !emcRecPoints || !cpvRecPoints || !trackSegments ) { + if ( !fEMCRecPoints || !fCPVRecPoints || !fTrackSegments ) { AliFatal("RecPoints or TrackSegments not found !") ; } - TClonesArray * recParticles = gime->RecParticles() ; - recParticles->Clear(); + fRecParticles->Clear(); - TIter next(trackSegments) ; + TIter next(fTrackSegments) ; AliPHOSTrackSegment * ts ; Int_t index = 0 ; AliPHOSRecParticle * rp ; while ( (ts = (AliPHOSTrackSegment *)next()) ) { // cout<<">>>>>>>>>>>>>>>PCA Index "<At(index) ; + new( (*fRecParticles)[index] ) AliPHOSRecParticle() ; + rp = (AliPHOSRecParticle *)fRecParticles->At(index) ; rp->SetTrackSegment(index) ; rp->SetIndexInList(index) ; AliPHOSEmcRecPoint * emc = 0 ; if(ts->GetEmcIndex()>=0) - emc = (AliPHOSEmcRecPoint *) emcRecPoints->At(ts->GetEmcIndex()) ; + emc = (AliPHOSEmcRecPoint *) fEMCRecPoints->At(ts->GetEmcIndex()) ; AliPHOSCpvRecPoint * cpv = 0 ; if(ts->GetCpvIndex()>=0) - cpv = (AliPHOSCpvRecPoint *) cpvRecPoints->At(ts->GetCpvIndex()) ; + cpv = (AliPHOSCpvRecPoint *) fCPVRecPoints->At(ts->GetCpvIndex()) ; Int_t track = 0 ; track = ts->GetTrackIndex() ; @@ -1493,11 +1455,10 @@ void AliPHOSPIDv1::MakeRecParticles() rp->SetLastDaughter(-1); rp->SetPolarisation(0,0,0); //Set the position in global coordinate system from the RecPoint - AliPHOSGeometry * geom = gime->PHOSGeometry() ; - AliPHOSTrackSegment * ts = gime->TrackSegment(rp->GetPHOSTSIndex()) ; - AliPHOSEmcRecPoint * erp = gime->EmcRecPoint(ts->GetEmcIndex()) ; + AliPHOSTrackSegment * ts = static_cast(fTrackSegments->At(rp->GetPHOSTSIndex())); + AliPHOSEmcRecPoint * erp = static_cast(fEMCRecPoints->At(ts->GetEmcIndex())); TVector3 pos ; - geom->GetGlobal(erp, pos) ; + fGeom->GetGlobalPHOS(erp, pos) ; rp->SetPos(pos); index++ ; } @@ -1528,23 +1489,17 @@ void AliPHOSPIDv1::PrintRecParticles(Option_t * option) { // Print table of reconstructed particles - AliPHOSGetter *gime = AliPHOSGetter::Instance() ; - - TClonesArray * recParticles = gime->RecParticles() ; - TString message ; - message = "\nevent " ; - message += gime->EventNumber(); - message += " found " ; - message += recParticles->GetEntriesFast(); + message = " found " ; + message += fRecParticles->GetEntriesFast(); message += " RecParticles\n" ; if(strstr(option,"all")) { // printing found TS message += "\n PARTICLE Index \n" ; Int_t index ; - for (index = 0 ; index < recParticles->GetEntries() ; index++) { - AliPHOSRecParticle * rp = (AliPHOSRecParticle * ) recParticles->At(index) ; + for (index = 0 ; index < fRecParticles->GetEntries() ; index++) { + AliPHOSRecParticle * rp = (AliPHOSRecParticle * ) fRecParticles->At(index) ; message += "\n" ; message += rp->Name().Data() ; message += " " ; @@ -1704,39 +1659,6 @@ void AliPHOSPIDv1::SetParameterToCalculateEllipse(TString particle, TString par (*fParameters)(p,i) = par ; } -//____________________________________________________________________________ -void AliPHOSPIDv1::Unload() -{ - //Unloads RecPoints, Tracks and RecParticles - AliPHOSGetter * gime = AliPHOSGetter::Instance() ; - gime->PhosLoader()->UnloadRecPoints() ; - gime->PhosLoader()->UnloadTracks() ; - gime->PhosLoader()->UnloadRecParticles() ; -} - -//____________________________________________________________________________ -void AliPHOSPIDv1::WriteRecParticles() -{ - //It writes reconstructed particles and pid to file - - AliPHOSGetter *gime = AliPHOSGetter::Instance() ; - - TClonesArray * recParticles = gime->RecParticles() ; - recParticles->Expand(recParticles->GetEntriesFast() ) ; - if(fWrite){ - TTree * treeP = gime->TreeP(); - - //First rp - Int_t bufferSize = 32000 ; - TBranch * rpBranch = treeP->Branch("PHOSRP",&recParticles,bufferSize); - rpBranch->SetTitle(BranchName()); - - rpBranch->Fill() ; - - gime->WriteRecParticles("OVERWRITE"); - gime->WritePID("OVERWRITE"); - } -} //____________________________________________________________________________ void AliPHOSPIDv1::GetVertex(void) { //extract vertex either using ESD or generator @@ -1744,12 +1666,14 @@ void AliPHOSPIDv1::GetVertex(void) //Try to extract vertex from data if(fESD){ const AliESDVertex *esdVtx = fESD->GetVertex() ; - if(esdVtx){ + if(esdVtx && esdVtx->GetChi2()!=0.){ fVtx.SetXYZ(esdVtx->GetXv(),esdVtx->GetYv(),esdVtx->GetZv()) ; return ; } } - + + // Use vertex diamond from CDB GRP folder if the one from ESD is missing + // PLEASE FIX IT AliWarning("Can not read vertex from data, use fixed \n") ; fVtx.SetXYZ(0.,0.,0.) ;