X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PHOS%2FAliPHOSClusterizerv1.cxx;h=840ffd05454774e42831d4c4172e7c48c6b2a4ed;hb=32abf0b9c3fa3a93907a70b3b950f8a87a0bf667;hp=1ed11b2adc8f73f2405a604ee2d7cf8212518437;hpb=78b28884c904b2adf03dfba0c309f35f2da6f4e8;p=u%2Fmrichter%2FAliRoot.git diff --git a/PHOS/AliPHOSClusterizerv1.cxx b/PHOS/AliPHOSClusterizerv1.cxx index 1ed11b2adc8..840ffd05454 100644 --- a/PHOS/AliPHOSClusterizerv1.cxx +++ b/PHOS/AliPHOSClusterizerv1.cxx @@ -180,7 +180,6 @@ #include "AliPHOSCpvRecPoint.h" #include "AliPHOSDigit.h" #include "AliPHOSDigitizer.h" -#include "AliPHOSCalibrationDB.h" #include "AliCDBManager.h" #include "AliCDBStorage.h" #include "AliCDBEntry.h" @@ -198,11 +197,17 @@ AliPHOSClusterizerv1::AliPHOSClusterizerv1() : fNumberOfEmcClusters(0), fNumberOfCpvClusters(0), fEmcClusteringThreshold(0), fCpvClusteringThreshold(0), fEmcLocMaxCut(0), fW0(0), fCpvLocMaxCut(0), - fW0CPV(0), fEmcTimeGate(0), fEcoreRadius(0) + fW0CPV(0), + fTimeGateLowAmp(0.), fTimeGateLow(0.), fTimeGateHigh(0.), + fEcoreRadius(0) { // default ctor (to be used mainly by Streamer) - fDefaultInit = kTRUE ; + fDefaultInit = kTRUE ; + + for(Int_t i=0; i<53760; i++){ + fDigitsUsed[i]=0 ; + } } //____________________________________________________________________________ @@ -213,10 +218,16 @@ AliPHOSClusterizerv1::AliPHOSClusterizerv1(AliPHOSGeometry *geom) : fNumberOfEmcClusters(0), fNumberOfCpvClusters(0), fEmcClusteringThreshold(0), fCpvClusteringThreshold(0), fEmcLocMaxCut(0), fW0(0), fCpvLocMaxCut(0), - fW0CPV(0), fEmcTimeGate(0), fEcoreRadius(0) + fW0CPV(0), + fTimeGateLowAmp(0.), fTimeGateLow(0.), fTimeGateHigh(0.), + fEcoreRadius(0) { // ctor with the indication of the file where header Tree and digits Tree are stored + for(Int_t i=0; i<53760; i++){ + fDigitsUsed[i]=0 ; + } + Init() ; fDefaultInit = kFALSE ; } @@ -275,6 +286,8 @@ Bool_t AliPHOSClusterizerv1::FindFit(AliPHOSEmcRecPoint * emcRP, AliPHOSDigit ** // Cluster will be fitted as a superposition of nPar/3 electromagnetic showers + if(!gMinuit) //it was deleted by someone else + gMinuit = new TMinuit(100) ; gMinuit->mncler(); // Reset Minuit's list of paramters gMinuit->SetPrintLevel(-1) ; // No Printout gMinuit->SetFCN(AliPHOSClusterizerv1::UnfoldingChiSquare) ; @@ -371,6 +384,8 @@ void AliPHOSClusterizerv1::Init() fgCalibData = new AliPHOSCalibData(-1); //use AliCDBManager's run number if (fgCalibData->GetCalibDataEmc() == 0) AliFatal("Calibration parameters for PHOS EMC not found. Stop reconstruction.\n"); + if (fgCalibData->GetCalibDataCpv() == 0) + AliFatal("Calibration parameters for PHOS CPV not found. Stop reconstruction.\n"); } @@ -384,8 +399,10 @@ void AliPHOSClusterizerv1::InitParameters() const AliPHOSRecoParam* recoParam = AliPHOSReconstructor::GetRecoParam(); if(!recoParam) AliFatal("Reconstruction parameters are not set!"); - fCpvClusteringThreshold = recoParam->GetEMCClusteringThreshold(); - fEmcClusteringThreshold = recoParam->GetCPVClusteringThreshold(); + recoParam->Print(); + + fEmcClusteringThreshold = recoParam->GetEMCClusteringThreshold(); + fCpvClusteringThreshold = recoParam->GetCPVClusteringThreshold(); fEmcLocMaxCut = recoParam->GetEMCLocalMaxCut(); fCpvLocMaxCut = recoParam->GetCPVLocalMaxCut(); @@ -393,7 +410,10 @@ void AliPHOSClusterizerv1::InitParameters() fW0 = recoParam->GetEMCLogWeight(); fW0CPV = recoParam->GetCPVLogWeight(); - fEmcTimeGate = 1.e-6 ; + fTimeGateLowAmp = recoParam->GetTimeGateAmpThresh() ; + fTimeGateLow = recoParam->GetTimeGateLow() ; + fTimeGateHigh = recoParam->GetTimeGateHigh() ; + fEcoreRadius = recoParam->GetEMCEcoreRadius(); fToUnfold = recoParam->EMCToUnfold() ; @@ -424,7 +444,7 @@ Int_t AliPHOSClusterizerv1::AreNeighbours(AliPHOSDigit * d1, AliPHOSDigit * d2)c if (( coldiff <= 1 ) && ( rowdiff <= 1 )){ //At least common vertex // if (( relid1[2]==relid2[2] && coldiff <= 1 ) || ( relid1[3]==relid2[3] && rowdiff <= 1 )){ //common side - if((relid1[1] != 0) || (TMath::Abs(d1->GetTime() - d2->GetTime() ) < fEmcTimeGate)) + if((relid1[1] != 0) || CheckTimeGate(d1->GetTime(),d1->GetEnergy(),d2->GetTime(),d2->GetEnergy())) return 1 ; } else { @@ -447,6 +467,18 @@ Int_t AliPHOSClusterizerv1::AreNeighbours(AliPHOSDigit * d1, AliPHOSDigit * d2)c return 0 ; } //____________________________________________________________________________ +Bool_t AliPHOSClusterizerv1::CheckTimeGate(Float_t t1, Float_t amp1, Float_t t2, Float_t amp2)const{ + //Check if two cells have reasonable time difference + //Note that at low amplitude time is defined up to 1 tick == 100 ns. + if(amp1( fEMCRecPoints->At(index) ); - rp->Purify(emcMinE) ; + static_cast( fEMCRecPoints->At(index) ); + rp->Purify(emcMinE,fDigitsArr) ; if(rp->GetMultiplicity()==0){ fEMCRecPoints->RemoveAt(index) ; delete rp ; @@ -506,7 +538,7 @@ void AliPHOSClusterizerv1::WriteRecPoints() fEMCRecPoints->Sort() ; // fEMCRecPoints->Expand(fEMCRecPoints->GetEntriesFast()) ; for(index = 0; index < fEMCRecPoints->GetEntries(); index++){ - dynamic_cast( fEMCRecPoints->At(index) )->SetIndexInList(index) ; + static_cast( fEMCRecPoints->At(index) )->SetIndexInList(index) ; } //For each rec.point set the distance to the nearest bad crystal (BVP) @@ -514,7 +546,7 @@ void AliPHOSClusterizerv1::WriteRecPoints() //Now the same for CPV for(index = 0; index < fCPVRecPoints->GetEntries(); index++){ - AliPHOSCpvRecPoint * rp = dynamic_cast( fCPVRecPoints->At(index) ); + AliPHOSCpvRecPoint * rp = static_cast( fCPVRecPoints->At(index) ); rp->EvalAll(fDigitsArr) ; rp->EvalAll(fW0CPV,fakeVtx,fDigitsArr) ; rp->EvalLocal2TrackingCSTransform(); @@ -522,7 +554,7 @@ void AliPHOSClusterizerv1::WriteRecPoints() fCPVRecPoints->Sort() ; for(index = 0; index < fCPVRecPoints->GetEntries(); index++) - dynamic_cast( fCPVRecPoints->At(index) )->SetIndexInList(index) ; + static_cast( fCPVRecPoints->At(index) )->SetIndexInList(index) ; fCPVRecPoints->Expand(fCPVRecPoints->GetEntriesFast()) ; @@ -541,12 +573,8 @@ void AliPHOSClusterizerv1::MakeClusters() fNumberOfEmcClusters = 0 ; //Mark all digits as unused yet - const Int_t maxNDigits = 1500; + const Int_t maxNDigits = 3584; // There is no clusters larger than PHOS module ;) Int_t nDigits=fDigitsArr->GetEntriesFast() ; - if (nDigits > maxNDigits) { - AliWarning(Form("Skip event with too high digit occupancy: nDigits=%d",nDigits)); - return; - } for(Int_t i=0; iGetEnergy() > fEmcClusteringThreshold ) || - ( IsInCpv(digit) && digit->GetEnergy() > fCpvClusteringThreshold ) ) { + if (( IsInEmc(digit) && Calibrate(digit->GetEnergy(),digit->GetId()) > fEmcClusteringThreshold ) || + ( IsInCpv(digit) && Calibrate(digit->GetEnergy(),digit->GetId()) > fCpvClusteringThreshold ) ) { Int_t iDigitInCluster = 0 ; - if ( IsInEmc(digit) ) { // start a new EMC RecPoint if(fNumberOfEmcClusters >= fEMCRecPoints->GetSize()) @@ -579,7 +606,7 @@ void AliPHOSClusterizerv1::MakeClusters() fEMCRecPoints->AddAt(new AliPHOSEmcRecPoint(""), fNumberOfEmcClusters) ; clu = static_cast( fEMCRecPoints->At(fNumberOfEmcClusters) ) ; fNumberOfEmcClusters++ ; - clu->AddDigit(*digit, digit->GetEnergy()) ; + clu->AddDigit(*digit, Calibrate(digit->GetEnergy(),digit->GetId()),CalibrateT(digit->GetTime(),digit->GetId(),digit->IsLG())) ; clusterdigitslist[iDigitInCluster] = digit->GetIndexInList() ; iDigitInCluster++ ; fDigitsUsed[i]=kTRUE ; @@ -591,12 +618,12 @@ void AliPHOSClusterizerv1::MakeClusters() fCPVRecPoints->AddAt(new AliPHOSCpvRecPoint(""), fNumberOfCpvClusters) ; clu = static_cast( fCPVRecPoints->At(fNumberOfCpvClusters) ) ; fNumberOfCpvClusters++ ; - clu->AddDigit(*digit, digit->GetEnergy()) ; + clu->AddDigit(*digit, Calibrate(digit->GetEnergy(),digit->GetId()),0.) ; // no timing information in CPV clusterdigitslist[iDigitInCluster] = digit->GetIndexInList() ; iDigitInCluster++ ; fDigitsUsed[i]=kTRUE ; } // else - + //Now scan remaining digits in list to find neigbours of our seed AliPHOSDigit * digitN ; @@ -605,6 +632,11 @@ void AliPHOSClusterizerv1::MakeClusters() digit = static_cast( fDigitsArr->At(clusterdigitslist[index]) ) ; index++ ; for(Int_t j=iFirst; j= maxNDigits) { + AliError(Form("The number of digits in cluster is more than %d, skip the rest of event", + maxNDigits)); + return; + } if(fDigitsUsed[j]) continue ; //look through remaining digits digitN = static_cast( fDigitsArr->At(j) ) ; @@ -616,22 +648,22 @@ void AliPHOSClusterizerv1::MakeClusters() case 0 : // not a neighbour break ; case 1 : // are neighbours - clu->AddDigit(*digitN, digitN->GetEnergy()); + clu->AddDigit(*digitN, Calibrate(digitN->GetEnergy(),digitN->GetId()),CalibrateT(digitN->GetTime(),digitN->GetId(),digit->IsLG())) ; clusterdigitslist[iDigitInCluster] = j ; iDigitInCluster++ ; fDigitsUsed[j]=kTRUE ; break ; case 2 : // too far from each other - goto endofloop; + goto endOfLoop; } // switch } - endofloop: ; //scanned all possible neighbours for this digit + endOfLoop: ; //scanned all possible neighbours for this digit } // loop over cluster } // energy theshold - } + } } @@ -650,7 +682,7 @@ void AliPHOSClusterizerv1::MakeUnfolding() Int_t index ; for(index = 0 ; index < numberofNotUnfolded ; index++){ - AliPHOSEmcRecPoint * emcRecPoint = dynamic_cast( fEMCRecPoints->At(index) ) ; + AliPHOSEmcRecPoint * emcRecPoint = static_cast( fEMCRecPoints->At(index) ) ; if(emcRecPoint->GetPHOSMod()> nModulesToUnfold) break ; @@ -688,12 +720,12 @@ void AliPHOSClusterizerv1::MakeUnfolding() Int_t index ; for(index = 0 ; index < numberofCpvNotUnfolded ; index++){ - AliPHOSRecPoint * recPoint = dynamic_cast( fCPVRecPoints->At(index) ) ; + AliPHOSRecPoint * recPoint = static_cast( fCPVRecPoints->At(index) ) ; if(recPoint->GetPHOSMod()> nModulesToUnfold) break ; - AliPHOSEmcRecPoint * emcRecPoint = dynamic_cast(recPoint) ; + AliPHOSEmcRecPoint * emcRecPoint = static_cast(recPoint) ; Int_t nMultipl = emcRecPoint->GetMultiplicity() ; AliPHOSDigit ** maxAt = new AliPHOSDigit*[nMultipl] ; @@ -769,7 +801,7 @@ void AliPHOSClusterizerv1::UnfoldCluster(AliPHOSEmcRecPoint * iniEmc, Int_t iparam ; Int_t iDigit ; for(iDigit = 0 ; iDigit < nDigits ; iDigit ++){ - digit = dynamic_cast( fDigitsArr->At(emcDigits[iDigit] ) ) ; + digit = static_cast( fDigitsArr->At(emcDigits[iDigit] ) ) ; fGeom->AbsToRelNumbering(digit->GetId(), relid) ; fGeom->RelPosInModule(relid, xDigit, zDigit) ; efit[iDigit] = 0; @@ -809,7 +841,7 @@ void AliPHOSClusterizerv1::UnfoldCluster(AliPHOSEmcRecPoint * iniEmc, fEMCRecPoints->Expand(2*fNumberOfEmcClusters) ; (*fEMCRecPoints)[fNumberOfEmcClusters] = new AliPHOSEmcRecPoint("") ; - emcRP = dynamic_cast( fEMCRecPoints->At(fNumberOfEmcClusters) ) ; + emcRP = static_cast( fEMCRecPoints->At(fNumberOfEmcClusters) ) ; fNumberOfEmcClusters++ ; emcRP->SetNExMax((Int_t)nPar/3) ; } @@ -818,19 +850,19 @@ void AliPHOSClusterizerv1::UnfoldCluster(AliPHOSEmcRecPoint * iniEmc, fCPVRecPoints->Expand(2*fNumberOfCpvClusters) ; (*fCPVRecPoints)[fNumberOfCpvClusters] = new AliPHOSCpvRecPoint("") ; - emcRP = dynamic_cast( fCPVRecPoints->At(fNumberOfCpvClusters) ) ; + emcRP = static_cast( fCPVRecPoints->At(fNumberOfCpvClusters) ) ; fNumberOfCpvClusters++ ; } Float_t eDigit ; for(iDigit = 0 ; iDigit < nDigits ; iDigit ++){ - digit = dynamic_cast( fDigitsArr->At( emcDigits[iDigit] ) ) ; + digit = static_cast( fDigitsArr->At( emcDigits[iDigit] ) ) ; fGeom->AbsToRelNumbering(digit->GetId(), relid) ; fGeom->RelPosInModule(relid, xDigit, zDigit) ; // ratio = epar * ShowerShape(xDigit - xpar,zDigit - zpar,vIncid) / efit[iDigit] ; ratio = epar * ShowerShape(xDigit - xpar,zDigit - zpar) / efit[iDigit] ; eDigit = emcEnergies[iDigit] * ratio ; - emcRP->AddDigit( *digit, eDigit ) ; + emcRP->AddDigit( *digit, eDigit,CalibrateT(digit->GetTime(),digit->GetId(),digit->IsLG()) ) ; } } @@ -845,17 +877,17 @@ void AliPHOSClusterizerv1::UnfoldingChiSquare(Int_t & nPar, Double_t * Grad, Dou // Calculates the Chi square for the cluster unfolding minimization // Number of parameters, Gradient, Chi squared, parameters, what to do - TList * toMinuit = dynamic_cast( gMinuit->GetObjectFit() ) ; + TList * toMinuit = static_cast( gMinuit->GetObjectFit() ) ; - AliPHOSEmcRecPoint * emcRP = dynamic_cast( toMinuit->At(0) ) ; - TClonesArray * digits = dynamic_cast( toMinuit->At(1) ) ; + AliPHOSEmcRecPoint * emcRP = static_cast( toMinuit->At(0) ) ; + TClonesArray * digits = static_cast( toMinuit->At(1) ) ; // A bit buggy way to get an access to the geometry // To be revised! - AliPHOSGeometry *geom = dynamic_cast(toMinuit->At(2)); + AliPHOSGeometry *geom = static_cast(toMinuit->At(2)); -// TVector3 * vtx = dynamic_cast(toMinuit->At(3)) ; //Vertex position +// TVector3 * vtx = static_cast(toMinuit->At(3)) ; //Vertex position - // AliPHOSEmcRecPoint * emcRP = dynamic_cast( gMinuit->GetObjectFit() ) ; // EmcRecPoint to fit + // AliPHOSEmcRecPoint * emcRP = static_cast( gMinuit->GetObjectFit() ) ; // EmcRecPoint to fit Int_t * emcDigits = emcRP->GetDigitsList() ; @@ -878,7 +910,7 @@ void AliPHOSClusterizerv1::UnfoldingChiSquare(Int_t & nPar, Double_t * Grad, Dou for( iDigit = 0 ; iDigit < nOdigits ; iDigit++) { - digit = dynamic_cast( digits->At( emcDigits[iDigit] ) ); + digit = static_cast( digits->At( emcDigits[iDigit] ) ); Int_t relid[4] ; Float_t xDigit ; @@ -974,7 +1006,7 @@ void AliPHOSClusterizerv1::Print(const Option_t *)const else message = " AliPHOSClusterizerv1 not initialized " ; - AliInfo(Form("%s, %s %s %s %s %s %s %s %s %s %s", message.Data(), + AliInfo(Form("%s, %s %s %s %s %f %f %f %f %f %f", message.Data(), taskName.Data(), GetTitle(), taskName.Data(), @@ -1042,9 +1074,9 @@ void AliPHOSClusterizerv1::SetDistancesToBadChannels() //Author: Boris Polichtchouk if(!fgCalibData->GetNumOfEmcBadChannels()) return; - AliInfo(Form("%d bad channel(s) found.\n",fgCalibData->GetNumOfEmcBadChannels())); Int_t badIds[8000]; + memset(badIds,0,8000*sizeof(Int_t)); fgCalibData->EmcBadChannelIds(badIds); AliPHOSEmcRecPoint* rp; @@ -1085,3 +1117,48 @@ void AliPHOSClusterizerv1::SetDistancesToBadChannels() } } +//================================================================================== +Float_t AliPHOSClusterizerv1::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 AliPHOSClusterizerv1::CalibrateT(Float_t time, Int_t absId,Bool_t isLG)const{ + // Calibrate time in EMC digit + + 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 + if(isLG) + time += fgCalibData->GetLGTimeShiftEmc(module,column,row); + else + time += fgCalibData->GetTimeShiftEmc(module,column,row); + return time ; + } +} +