X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PHOS%2FAliPHOSReconstructor.cxx;h=32120efb00c57acd58c9c0dd46b479b2821e83c0;hb=e079adac340cd38db0b469fe7b10ca2452a11102;hp=8792cebcb528ba5ecfb8444400c522c7839a34ea;hpb=48c5db5bcfa71b3891827afb90d9855268021b8c;p=u%2Fmrichter%2FAliRoot.git diff --git a/PHOS/AliPHOSReconstructor.cxx b/PHOS/AliPHOSReconstructor.cxx index 8792cebcb52..32120efb00c 100644 --- a/PHOS/AliPHOSReconstructor.cxx +++ b/PHOS/AliPHOSReconstructor.cxx @@ -22,6 +22,8 @@ // derived from STEER/AliReconstructor. // // --- ROOT system --- +#include "TGeoManager.h" +#include "TGeoMatrix.h" // --- Standard library --- @@ -37,51 +39,47 @@ #include "AliPHOSPIDv1.h" #include "AliPHOSTracker.h" #include "AliRawReader.h" +#include "AliPHOSCalibData.h" +#include "AliCDBEntry.h" +#include "AliCDBManager.h" #include "AliPHOSTrigger.h" #include "AliPHOSGeometry.h" -#include "AliPHOSRecoParam.h" -#include "AliPHOSRecoParamEmc.h" -#include "AliPHOSRecoParamCpv.h" #include "AliPHOSDigit.h" #include "AliPHOSTrackSegment.h" #include "AliPHOSEmcRecPoint.h" #include "AliPHOSRecParticle.h" -#include "AliPHOSRawDecoder.h" -#include "AliPHOSRawDecoderv1.h" -#include "AliPHOSRawDecoderv2.h" +#include "AliPHOSRawFitterv0.h" +#include "AliPHOSRawFitterv1.h" +#include "AliPHOSRawFitterv2.h" +#include "AliPHOSRawFitterv3.h" #include "AliPHOSRawDigiProducer.h" #include "AliPHOSPulseGenerator.h" ClassImp(AliPHOSReconstructor) Bool_t AliPHOSReconstructor::fgDebug = kFALSE ; -AliPHOSRecoParam* AliPHOSReconstructor::fgkRecoParamEmc =0; // EMC rec. parameters -AliPHOSRecoParam* AliPHOSReconstructor::fgkRecoParamCpv =0; // CPV rec. parameters TClonesArray* AliPHOSReconstructor::fgDigitsArray = 0; // Array of PHOS digits TObjArray* AliPHOSReconstructor::fgEMCRecPoints = 0; // Array of EMC rec.points +AliPHOSCalibData * AliPHOSReconstructor::fgCalibData = 0 ; + //____________________________________________________________________________ AliPHOSReconstructor::AliPHOSReconstructor() : fGeom(NULL),fClusterizer(NULL),fTSM(NULL),fPID(NULL) { // ctor - - if (!fgkRecoParamEmc) { - AliWarning("The Reconstruction parameters for EMC nonitialized - Used default one"); - fgkRecoParamEmc = AliPHOSRecoParamEmc::GetEmcDefaultParameters(); - } - - if (!fgkRecoParamCpv) { - AliWarning("The Reconstruction parameters for CPV nonitialized - Used default one"); - fgkRecoParamCpv = AliPHOSRecoParamCpv::GetCpvDefaultParameters(); - } - fGeom = AliPHOSGeometry::GetInstance("IHEP",""); fClusterizer = new AliPHOSClusterizerv1 (fGeom); fTSM = new AliPHOSTrackSegmentMakerv1(fGeom); fPID = new AliPHOSPIDv1 (fGeom); fgDigitsArray = new TClonesArray("AliPHOSDigit",100); fgEMCRecPoints= new TObjArray(100) ; + if (!fgCalibData) + fgCalibData = new AliPHOSCalibData(-1); //use AliCDBManager's run number + + AliInfo(Form("PHOS bad channel map contains %d bad channel(s).\n", + fgCalibData->GetNumOfEmcBadChannels())); + } //____________________________________________________________________________ @@ -104,6 +102,7 @@ void AliPHOSReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) c // segment maker needs access to the AliESDEvent object to retrieve the tracks reconstructed by // the global tracking. + fClusterizer->InitParameters(); fClusterizer->SetInput(digitsTree); fClusterizer->SetOutput(clustersTree); if ( Debug() ) @@ -129,6 +128,8 @@ void AliPHOSReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, else fTSM->Clusters2TrackSegments("") ; + fPID->SetEnergyCorrectionOn(GetRecoParam()->GetEMCEnergyCorrectionOn()); + fPID->SetInput(clustersTree, fTSM->GetTrackSegments()) ; fPID->SetESD(esd) ; if ( Debug() ) @@ -139,9 +140,6 @@ void AliPHOSReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, TClonesArray *recParticles = fPID->GetRecParticles(); Int_t nOfRecParticles = recParticles->GetEntriesFast(); - esd->SetNumberOfPHOSClusters(nOfRecParticles) ; - esd->SetFirstPHOSCluster(esd->GetNumberOfCaloClusters()) ; - AliDebug(2,Form("%d rec. particles, option %s",nOfRecParticles,GetOption())); // Read digits array @@ -236,10 +234,10 @@ void AliPHOSReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, // Add to CaloCells only EMC digits with non-zero energy for (Int_t idig = 0 ; idig < nDigits ; idig++) { const AliPHOSDigit * dig = (const AliPHOSDigit*)fgDigitsArray->At(idig); - if(dig->GetId() <= knEMC && dig->GetEnergy() > 0 ){ - //printf("i %d; id %d; amp %f; time %e\n", - //idignew,dig->GetId(),dig->GetEnergy(), dig->GetTime()); - phsCells.SetCell(idignew,dig->GetId(), dig->GetEnergy(), dig->GetTime()); + if(dig->GetId() <= knEMC && + Calibrate(dig->GetEnergy(),dig->GetId()) > GetRecoParam()->GetEMCMinE() ){ + phsCells.SetCell(idignew,dig->GetId(), Calibrate(dig->GetEnergy(),dig->GetId()), + CalibrateT(dig->GetTime(),dig->GetId())); idignew++; } } @@ -279,7 +277,7 @@ void AliPHOSReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, AliPHOSDigit *digit = static_cast(fgDigitsArray->At(digitsList[iCell])); absIdList[iCell] = (UShort_t)(digit->GetId()); if (digit->GetEnergy() > 0) - fracList[iCell] = rpElist[iCell]/digit->GetEnergy(); + fracList[iCell] = rpElist[iCell]/(Calibrate(digit->GetEnergy(),digit->GetId())); else fracList[iCell] = 0; } @@ -289,7 +287,7 @@ void AliPHOSReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, Int_t *primList = emcRP->GetPrimaries(primMult); Float_t energy; - if (fgkRecoParamEmc->Ecore2ESD()) + if (GetRecoParam()->EMCEcore2ESD()) energy = emcRP->GetCoreEnergy(); else energy = rp->Energy(); @@ -304,8 +302,9 @@ void AliPHOSReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, ec->SetM20(emcRP->GetM2z()) ; //second moment M2z ec->SetNExMax(emcRP->GetNExMax()); //number of local maxima ec->SetEmcCpvDistance(ts->GetCpvDistance("r")); //Only radius, what about separate x,z???? + ec->SetTrackDistance(ts->GetCpvDistance("x"),ts->GetCpvDistance("z")); ec->SetClusterChi2(-1); //not yet implemented - ec->SetTOF(emcRP->GetTime()); //Time of flight + ec->SetTOF(emcRP->GetTime()); //Time of flight - already calibrated in EMCRecPoint //Cells contributing to clusters ec->SetNCells(cellMult); @@ -318,15 +317,27 @@ void AliPHOSReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, //Array of MC indeces TArrayI arrayPrim(primMult,primList); ec->AddLabels(arrayPrim); - - //Array of tracks uncomment when available in future - //TArrayS arrayTrackMatched(1);// Only one track, temporal solution. - //arrayTrackMatched[0]= (Short_t)(matchedTrack[iClust]); - //ec->AddTracksMatched(arrayTrackMatched); - // add the track to the esd object - - esd->AddCaloCluster(ec); + //Matched ESD track + TArrayI arrayTrackMatched(1); + arrayTrackMatched[0]= ts->GetTrackIndex(); + ec->AddTracksMatched(arrayTrackMatched); + + Int_t index = esd->AddCaloCluster(ec); + + //Set pointer to this cluster in ESD track + Int_t nt=esd->GetNumberOfTracks(); + for (Int_t itr=0; itrGetTrack(itr); + if(!esdTrack->IsPHOS()) + continue ; + if(esdTrack->GetPHOScluster()==-recpart){ //we store negative cluster number + esdTrack->SetPHOScluster(index) ; +//no garatie that only one track matched this cluster +// break ; + } + } + delete ec; delete [] fracList; delete [] absIdList; @@ -334,6 +345,10 @@ void AliPHOSReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree, fgDigitsArray ->Delete(); fgEMCRecPoints->Delete(); recParticles ->Delete(); + + //Store PHOS misalignment matrixes + FillMisalMatrixes(esd) ; + } //____________________________________________________________________________ @@ -349,63 +364,142 @@ void AliPHOSReconstructor::ConvertDigits(AliRawReader* rawReader, TTree* digits // Converts raw data to // PHOS digits // Works on a single-event basis - rawReader->Reset() ; - AliPHOSRawDecoder * dc ; + AliPHOSRawFitterv0 * fitter ; - const TObjArray* maps = AliPHOSRecoParamEmc::GetMappings(); + const TObjArray* maps = AliPHOSRecoParam::GetMappings(); if(!maps) AliFatal("Cannot retrieve ALTRO mappings!!"); - AliAltroMapping *mapping[4]; - for(Int_t i = 0; i < 4; i++) { + AliAltroMapping *mapping[20]; + for(Int_t i = 0; i < 20; i++) { mapping[i] = (AliAltroMapping*)maps->At(i); } - if(strcmp(fgkRecoParamEmc->DecoderVersion(),"v1")==0) - dc=new AliPHOSRawDecoderv1(rawReader,mapping); - else - if(strcmp(fgkRecoParamEmc->DecoderVersion(),"v2")==0) - dc=new AliPHOSRawDecoderv2(rawReader,mapping); - else - dc=new AliPHOSRawDecoder(rawReader,mapping); + if (strcmp(GetRecoParam()->EMCFitterVersion(),"v1")==0) + fitter=new AliPHOSRawFitterv1(); + else if (strcmp(GetRecoParam()->EMCFitterVersion(),"v2")==0) + fitter=new AliPHOSRawFitterv2(); + else if (strcmp(GetRecoParam()->EMCFitterVersion(),"v3")==0) + fitter=new AliPHOSRawFitterv3(); + else + fitter=new AliPHOSRawFitterv0(); + + fitter->SubtractPedestals(GetRecoParam()->EMCSubtractPedestals()); + fitter->SetAmpOffset (GetRecoParam()->GetGlobalAltroOffset()); + fitter->SetAmpThreshold (GetRecoParam()->GetGlobalAltroThreshold()); - dc->SubtractPedestals(fgkRecoParamEmc->SubtractPedestals()); - TClonesArray *digits = new TClonesArray("AliPHOSDigit",1); digits->SetName("DIGITS"); Int_t bufsize = 32000; digitsTree->Branch("PHOS", &digits, bufsize); - AliPHOSRawDigiProducer pr(fgkRecoParamEmc,fgkRecoParamCpv); - pr.MakeDigits(digits,dc); + AliPHOSRawDigiProducer rdp(rawReader,mapping); - delete dc ; + rdp.SetEmcMinAmp(GetRecoParam()->GetEMCRawDigitThreshold()); // in ADC + rdp.SetCpvMinAmp(GetRecoParam()->GetCPVMinE()); + rdp.SetSampleQualityCut(GetRecoParam()->GetEMCSampleQualityCut()); + rdp.MakeDigits(digits,fitter); - //!!!!for debug!!! -/* - Int_t modMax=-111; - Int_t colMax=-111; - Int_t rowMax=-111; - Float_t eMax=-333; - //!!!for debug!!! + delete fitter ; - Int_t relId[4]; - for(Int_t iDigit=0; iDigitGetEntries(); iDigit++) { - AliPHOSDigit* digit = (AliPHOSDigit*)digits->At(iDigit); - if(digit->GetEnergy()>eMax) { - fGeom->AbsToRelNumbering(digit->GetId(),relId); - eMax=digit->GetEnergy(); - modMax=relId[0]; - rowMax=relId[2]; - colMax=relId[3]; + if (AliLog::GetGlobalDebugLevel() == 1) { + Int_t modMax=-111; + Int_t colMax=-111; + Int_t rowMax=-111; + Float_t eMax=-333; + //!!!for debug!!! + + Int_t relId[4]; + for(Int_t iDigit=0; iDigitGetEntries(); iDigit++) { + AliPHOSDigit* digit = (AliPHOSDigit*)digits->At(iDigit); + if(digit->GetEnergy()>eMax) { + fGeom->AbsToRelNumbering(digit->GetId(),relId); + eMax=digit->GetEnergy(); + modMax=relId[0]; + rowMax=relId[2]; + colMax=relId[3]; + } } + + AliDebug(1,Form("Digit with max. energy: modMax %d colMax %d rowMax %d eMax %f\n\n", + modMax,colMax,rowMax,eMax)); } - AliDebug(1,Form("Digit with max. energy: modMax %d colMax %d rowMax %d eMax %f\n\n", - modMax,colMax,rowMax,eMax)); -*/ digitsTree->Fill(); digits->Delete(); delete digits; } +//================================================================================== +Float_t AliPHOSReconstructor::Calibrate(Float_t amp, Int_t absId)const{ + // Calibrate EMC digit, i.e. multiply its Amp by a factor read from CDB + + const AliPHOSGeometry *geom = AliPHOSGeometry::GetInstance() ; + + //Determine rel.position of the cell absolute ID + Int_t relId[4]; + geom->AbsToRelNumbering(absId,relId); + Int_t module=relId[0]; + Int_t row =relId[2]; + Int_t column=relId[3]; + if(relId[1]){ //CPV + Float_t calibration = fgCalibData->GetADCchannelCpv(module,column,row); + return amp*calibration ; + } + else{ //EMC + Float_t calibration = fgCalibData->GetADCchannelEmc(module,column,row); + return amp*calibration ; + } +} +//================================================================================== +Float_t AliPHOSReconstructor::CalibrateT(Float_t time, Int_t absId)const{ + // Calibrate EMC digit, i.e. multiply its Amp by a factor read from CDB + + const AliPHOSGeometry *geom = AliPHOSGeometry::GetInstance() ; + + //Determine rel.position of the cell absolute ID + Int_t relId[4]; + geom->AbsToRelNumbering(absId,relId); + Int_t module=relId[0]; + Int_t row =relId[2]; + Int_t column=relId[3]; + if(relId[1]){ //CPV + return 0. ; + } + else{ //EMC + time += fgCalibData->GetTimeShiftEmc(module,column,row); + return time ; + } +} +//================================================================================== +void AliPHOSReconstructor::FillMisalMatrixes(AliESDEvent* esd)const{ + //Store PHOS matrixes in ESD Header + + //Check, if matrixes was already stored + for(Int_t mod=0 ;mod<5; mod++){ + if(esd->GetPHOSMatrix(mod)!=0) + return ; + } + + //Create and store matrixes + if(!gGeoManager){ + AliError("Can not store misal. matrixes: no gGeoManager! \n") ; + return ; + } + //Note, that owner of copied marixes will be header + char path[255] ; + TGeoHMatrix * m ; + for(Int_t mod=0; mod<5; mod++){ + sprintf(path,"/ALIC_1/PHOS_%d",mod+1) ; //In Geometry modules numbered 1,2,.,5 + if (gGeoManager->cd(path)){ + m = gGeoManager->GetCurrentMatrix() ; + esd->SetPHOSMatrix(new TGeoHMatrix(*m),mod) ; + } + else{ + esd->SetPHOSMatrix(NULL,mod) ; + } + } + +} + +