X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=ZDC%2FAliZDCDigitizer.cxx;h=d8b43c933dc682b885712eb99b8ff7fe19433c60;hb=0d76e8c6989aae84f17591621a687c63b0b9bd4b;hp=58cfb8199e3c092db15df5a34ac75db9858e4bbe;hpb=e24fa7ff740b209990c521159f4a79dcf4eadb0f;p=u%2Fmrichter%2FAliRoot.git diff --git a/ZDC/AliZDCDigitizer.cxx b/ZDC/AliZDCDigitizer.cxx index 58cfb8199e3..d8b43c933dc 100644 --- a/ZDC/AliZDCDigitizer.cxx +++ b/ZDC/AliZDCDigitizer.cxx @@ -36,6 +36,7 @@ #include "AliGenHijingEventHeader.h" #include "AliRunDigitizer.h" #include "AliRunLoader.h" +#include "AliLoader.h" #include "AliGRPObject.h" #include "AliCDBManager.h" #include "AliCDBEntry.h" @@ -47,8 +48,6 @@ class AliCDBStorage; class AliZDCPedestals; -class AliZDCEnCalib; -class AliZDCTowerCalib; ClassImp(AliZDCDigitizer) @@ -59,12 +58,11 @@ AliZDCDigitizer::AliZDCDigitizer() : fIsSignalInADCGate(kFALSE), fFracLostSignal(0.), fPedData(0), - fEnCalibData(0), - fTowCalibData(0), - fSpectators2Track(kFALSE) + fSpectators2Track(kFALSE), + fBeamEnergy(0.) { -// Default constructor - + // Default constructor + for(Int_t i=0; i<2; i++) fADCRes[i]=0.; } //____________________________________________________________________________ @@ -74,13 +72,15 @@ AliZDCDigitizer::AliZDCDigitizer(AliRunDigitizer* manager): fIsSignalInADCGate(kFALSE), fFracLostSignal(0.), fPedData(GetPedData()), - fEnCalibData(GetEnCalibData()), - fTowCalibData(GetTowCalibData()), - fSpectators2Track(kFALSE) + fSpectators2Track(kFALSE), + fBeamEnergy(0.) { // Get calibration data if(fIsCalibration!=0) printf("\n\t AliZDCDigitizer -> Creating calibration data (pedestals)\n"); - + for(Int_t i=0; i<5; i++){ + for(Int_t j=0; j<5; j++) + fPMGain[i][j] = 0.; + } } //____________________________________________________________________________ @@ -98,13 +98,12 @@ AliZDCDigitizer::AliZDCDigitizer(const AliZDCDigitizer &digitizer): fIsSignalInADCGate(digitizer.fIsSignalInADCGate), fFracLostSignal(digitizer.fFracLostSignal), fPedData(digitizer.fPedData), - fEnCalibData(digitizer.fEnCalibData), - fTowCalibData(digitizer.fTowCalibData), - fSpectators2Track(digitizer.fSpectators2Track) + fSpectators2Track(digitizer.fSpectators2Track), + fBeamEnergy(digitizer.fBeamEnergy) { // Copy constructor - for(Int_t i=0; i<6; i++){ + for(Int_t i=0; i<5; i++){ for(Int_t j=0; j<5; j++){ fPMGain[i][j] = digitizer.fPMGain[i][j]; } @@ -120,6 +119,7 @@ Bool_t AliZDCDigitizer::Init() // Initialize the digitizer AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/GRP/Data"); + if(!entry) AliFatal("No calibration data loaded!"); AliGRPObject* grpData = 0x0; if(entry){ TMap* m = dynamic_cast(entry->GetObject()); // old GRP entry @@ -134,46 +134,73 @@ Bool_t AliZDCDigitizer::Init() entry->SetOwner(0); AliCDBManager::Instance()->UnloadFromCache("GRP/GRP/Data"); } - if(!grpData) AliError("No GRP entry found in OCDB!"); + if(!grpData){ + AliError("No GRP entry found in OCDB! \n "); + return kFALSE; + } TString beamType = grpData->GetBeamType(); if(beamType==AliGRPObject::GetInvalidString()){ AliError("\t UNKNOWN beam type from GRP obj -> PMT gains not set in ZDC digitizer!!!\n"); } - Float_t beamEnergy = grpData->GetBeamEnergy(); - if(beamEnergy==AliGRPObject::GetInvalidFloat()){ + fBeamEnergy = grpData->GetBeamEnergy(); + if(fBeamEnergy==AliGRPObject::GetInvalidFloat()){ AliWarning("GRP/GRP/Data entry: missing value for the beam energy ! Using 0."); AliError("\t UNKNOWN beam type from GRP obj -> PMT gains not set in ZDC digitizer!!!\n"); - beamEnergy = 0.; + fBeamEnergy = 0.; } - - if((beamType.CompareTo("P-P")) == 0){ + + if((beamType.CompareTo("P-P")) == 0 || (beamType.CompareTo("p-p")) == 0){ //PTM gains rescaled to beam energy for p-p - if(beamEnergy != 0){ + // New correction coefficients for PMT gains needed + // to reproduce experimental spectra (from Grazia Jul 2010) + if(fBeamEnergy != 0){ for(Int_t j = 0; j < 5; j++){ - fPMGain[0][j] = (661.444/beamEnergy+0.000740671)*10000000; - fPMGain[1][j] = (864.350/beamEnergy+0.002344)*10000000; - fPMGain[3][j] = (1.32312-0.000101515*beamEnergy)*10000000; - fPMGain[4][j] = fPMGain[0][j]; - fPMGain[5][j] = fPMGain[1][j] ; + fPMGain[0][j] = 1.515831*(661.444/fBeamEnergy+0.000740671)*10000000; + fPMGain[1][j] = 0.674234*(864.350/fBeamEnergy+0.00234375)*10000000; + fPMGain[3][j] = 1.350938*(661.444/fBeamEnergy+0.000740671)*10000000; + fPMGain[4][j] = 0.678597*(864.350/fBeamEnergy+0.00234375)*10000000; } + fPMGain[2][1] = 0.869654*(1.32312-0.000101515*fBeamEnergy)*10000000; + fPMGain[2][2] = 1.030883*(1.32312-0.000101515*fBeamEnergy)*10000000; + // + AliInfo(Form(" PMT gains for p-p @ %1.0f+%1.0f GeV: ZN(%1.0f), ZP(%1.0f), ZEM(%1.0f)\n", + fBeamEnergy, fBeamEnergy, fPMGain[0][0], fPMGain[1][0], fPMGain[2][1])); + } + else{ // for RELDIS simulation + Float_t scalGainFactor = 0.5; + for(Int_t j = 0; j < 5; j++){ + fPMGain[0][j] = 50000./scalGainFactor; // ZNC + fPMGain[1][j] = 100000./scalGainFactor; // ZPC + fPMGain[2][j] = 100000./scalGainFactor; // ZEM + fPMGain[3][j] = 50000./scalGainFactor; // ZNA + fPMGain[4][j] = 100000./scalGainFactor; // ZPA + } + // + AliInfo(Form(" PMT gains for RELDIS simulation: ZN(%1.0f), ZP(%1.0f), ZEM(%1.0f)\n", + fPMGain[0][0], fPMGain[1][0], fPMGain[2][1])); } } else if((beamType.CompareTo("A-A")) == 0){ - // PTM gains for Pb-Pb @ 2.7_2.7 A TeV *************** + // PTM gains for Pb-Pb @ 2.7+2.7 A TeV *************** + // rescaled for Pb-Pb @ 1.38+1.38 A TeV *************** + // Values corrected after 2010 Pb-Pb data taking (7/2/2011 - Ch.) + // Experimental data compared to EMD simulation for single nucleon peaks: + // ZN gains must be divided by 4, ZP gains by 10! + Float_t scalGainFactor = fBeamEnergy/2760.; for(Int_t j = 0; j < 5; j++){ - fPMGain[0][j] = 50000.; - fPMGain[1][j] = 100000.; - fPMGain[2][j] = 100000.; - fPMGain[3][j] = 50000.; - fPMGain[4][j] = 100000.; - fPMGain[5][j] = 100000.; + fPMGain[0][j] = 50000./(4*scalGainFactor); // ZNC + fPMGain[1][j] = 100000./(5*scalGainFactor); // ZPC + fPMGain[2][j] = 100000./scalGainFactor; // ZEM + fPMGain[3][j] = 50000./(4*scalGainFactor); // ZNA + fPMGain[4][j] = 100000./(5*scalGainFactor); // ZPA } + AliInfo(Form(" PMT gains for Pb-Pb @ %1.0f+%1.0f A GeV: ZN(%1.0f), ZP(%1.0f), ZEM(%1.0f)\n", + fBeamEnergy, fBeamEnergy, fPMGain[0][0], fPMGain[1][0], fPMGain[2][1])); } - - - // ADC Caen V965 + + // ADC Caen V965 fADCRes[0] = 0.0000008; // ADC Resolution high gain: 200 fC/adcCh fADCRes[1] = 0.0000064; // ADC Resolution low gain: 25 fC/adcCh @@ -212,7 +239,7 @@ void AliZDCDigitizer::Exec(Option_t* /*option*/) } // impact parameter and number of spectators - Float_t impPar = -1; + Float_t impPar = 0; Int_t specNTarg = 0, specPTarg = 0; Int_t specNProj = 0, specPProj = 0; Float_t signalTime0 = 0.; @@ -257,6 +284,7 @@ void AliZDCDigitizer::Exec(Option_t* /*option*/) fFracLostSignal = (sdigit.GetTrackTime()-30)*(sdigit.GetTrackTime()-30)/280.; } } + Float_t sdSignal = sdigit.GetLightPM(); if(fIsSignalInADCGate == kFALSE){ AliDebug(2,Form("\t Signal time %f -> fraction %f of ZDC signal for det.(%d, %d) out of ADC gate\n", @@ -265,9 +293,10 @@ void AliZDCDigitizer::Exec(Option_t* /*option*/) } pm[(sdigit.GetSector(0))-1][sdigit.GetSector(1)] += sdigit.GetLightPM(); - /*printf("\n\t Detector %d, Tower %d -> pm[%d][%d] = %.0f \n", + //Ch. debug + /*printf("\t Detector %d, Tower %d -> pm[%d][%d] = %.0f \n", sdigit.GetSector(0), sdigit.GetSector(1),sdigit.GetSector(0)-1, - sdigit.GetSector(1), pm[sdigit.GetSector(0)-1][sdigit.GetSector(1)]); // Chiara debugging! + sdigit.GetSector(1), pm[sdigit.GetSector(0)-1][sdigit.GetSector(1)]); */ } @@ -288,10 +317,10 @@ void AliZDCDigitizer::Exec(Option_t* /*option*/) specPProj = ((AliGenHijingEventHeader*) genHeader)->ProjSpectatorsp(); specNTarg = ((AliGenHijingEventHeader*) genHeader)->TargSpectatorsn(); specPTarg = ((AliGenHijingEventHeader*) genHeader)->TargSpectatorsp(); - printf("\n\t AliZDCDigitizer: b = %1.2f fm\n" + /*printf("\n\t AliZDCDigitizer: b = %1.2f fm\n" " \t PROJ.: #spectator n %d, #spectator p %d\n" - " \t TARG.: #spectator n %d, #spectator p %d\n\n", - impPar, specNProj, specPProj, specNTarg, specPTarg); + " \t TARG.: #spectator n %d, #spectator p %d\n", + impPar, specNProj, specPProj, specNTarg, specPTarg);*/ } } @@ -303,13 +332,13 @@ void AliZDCDigitizer::Exec(Option_t* /*option*/) Int_t freeSpecNTarg, freeSpecPTarg; Fragmentation(impPar, specNTarg, specPTarg, freeSpecNTarg, freeSpecPTarg); SpectatorSignal(1, freeSpecNProj, pm); - //printf(" AliZDCDigitizer -> Signal for %d PROJ free spectator n",freeSpecNProj); +// printf("\t AliZDCDigitizer -> Adding signal for %d PROJ free spectator n",freeSpecNProj); SpectatorSignal(2, freeSpecPProj, pm); - //printf(" and %d free spectator p added\n",freeSpecPProj); +// printf(" and %d free spectator p\n",freeSpecPProj); SpectatorSignal(3, freeSpecNTarg, pm); - //printf(" AliZDCDigitizer -> Signal for %d TARG free spectator n",freeSpecNTarg); +// printf("\t AliZDCDigitizer -> Adding signal for %d TARG free spectator n",freeSpecNTarg); SpectatorSignal(4, freeSpecPTarg, pm); - //printf("and %d free spectator p added\n",freeSpecPTarg); +// printf(" and %d free spectator p\n\n",freeSpecPTarg); } @@ -343,12 +372,13 @@ void AliZDCDigitizer::Exec(Option_t* /*option*/) digi[res] = Phe2ADCch(sector[0], sector[1], pm[sector[0]-1][sector[1]], res) + Pedestal(sector[0], sector[1], res); } - /*printf("\t DIGIT added -> det %d quad %d - digi[0,1] = [%d, %d]\n", - sector[0], sector[1], digi[0], digi[1]); // Chiara debugging! - */ - // new(pdigit) AliZDCDigit(sector, digi); treeD->Fill(); + + //Ch. debug + //printf("\t DIGIT added -> det %d quad %d - digi[0,1] = [%d, %d]\n", + // sector[0], sector[1], digi[0], digi[1]); // Chiara debugging! + } } // Loop over detector // Adding in-time digits for 2 reference PTM signals (after signal ch.) @@ -365,11 +395,12 @@ void AliZDCDigitizer::Exec(Option_t* /*option*/) for(Int_t res=0; res<2; res++){ sigRef[res] += Pedestal(sectorRef[0], sectorRef[1], res); } - /*printf("\t RefDigit added -> det = %d, quad = %d - digi[0,1] = [%d, %d]\n", - sectorRef[0], sectorRef[1], sigRef[0], sigRef[1]); // Chiara debugging! - */ new(pdigit) AliZDCDigit(sectorRef, sigRef); treeD->Fill(); + + //Ch. debug + //printf("\t RefDigit added -> det = %d, quad = %d - digi[0,1] = [%d, %d]\n", + // sectorRef[0], sectorRef[1], sigRef[0], sigRef[1]); // Chiara debugging! } // // --- Adding digits for out-of-time channels after signal digits @@ -379,12 +410,12 @@ void AliZDCDigitizer::Exec(Option_t* /*option*/) for(Int_t res=0; res<2; res++){ digioot[res] = Pedestal(sector[0], sector[1], res); // out-of-time ADCs } - /*printf("\t DIGIToot added -> det = %d, quad = %d - digi[0,1] = [%d, %d]\n", - sector[0], sector[1], digioot[0], digioot[1]); // Chiara debugging! - */ - // new(pdigit) AliZDCDigit(sector, digioot); treeD->Fill(); + + //Ch. debug + //printf("\t DIGIToot added -> det = %d, quad = %d - digi[0,1] = [%d, %d]\n", + // sector[0], sector[1], digioot[0], digioot[1]); // Chiara debugging! } } // Adding out-of-time digits for 2 reference PTM signals (after out-of-time ch.) @@ -394,11 +425,12 @@ void AliZDCDigitizer::Exec(Option_t* /*option*/) for(Int_t res=0; res<2; res++){ sigRefoot[res] = Pedestal(sectorRef[0], sectorRef[1], res); } - /*printf("\t RefDigitoot added -> det = %d, quad = %d - digi[0,1] = [%d, %d]\n", - sectorRef[0], sectorRef[1], sigRefoot[0], sigRefoot[1]); // Chiara debugging! - */ new(pdigit) AliZDCDigit(sectorRef, sigRefoot); treeD->Fill(); + //Ch. debug + //printf("\t RefDigitoot added -> det = %d, quad = %d - digi[0,1] = [%d, %d]\n", + // sectorRef[0], sectorRef[1], sigRefoot[0], sigRefoot[1]); // Chiara debugging! + } //printf("\t AliZDCDigitizer -> TreeD has %d entries\n",(Int_t) treeD->GetEntries()); @@ -421,14 +453,18 @@ void AliZDCDigitizer::Fragmentation(Float_t impPar, Int_t specN, Int_t specP, Int_t nAlpha = frag.GetNalpha(); // Attach neutrons + frag.AttachNeutrons(); Int_t ztot = frag.GetZtot(); Int_t ntot = frag.GetNtot(); - frag.AttachNeutrons(); + + // Removing fragments and alpha pcs freeSpecN = specN-ntot-2*nAlpha; freeSpecP = specP-ztot-2*nAlpha; + // Removing deuterons Int_t ndeu = (Int_t) (freeSpecN*frag.DeuteronNumber()); freeSpecN -= ndeu; + freeSpecP -= ndeu; // if(freeSpecN<0) freeSpecN=0; if(freeSpecP<0) freeSpecP=0; @@ -441,27 +477,62 @@ void AliZDCDigitizer::SpectatorSignal(Int_t SpecType, Int_t numEvents, { // add signal of the spectators - TString hfn; - if(SpecType == 1) { // --- Signal for projectile spectator neutrons - hfn = "$ALICE_ROOT/ZDC/ZNCSignal.root"; - } - else if(SpecType == 2) { // --- Signal for projectile spectator protons - hfn = "$ALICE_ROOT/ZDC/ZPCSignal.root"; + TFile *specSignalFile = TFile::Open("$ALICE_ROOT/ZDC/SpectatorSignal.root"); + if(!specSignalFile || !specSignalFile->IsOpen()) { + AliError((" Opening file $ALICE_ROOT/ZDC/SpectatorSignal.root failed\n")); + return; } - else if(SpecType == 3) { // --- Signal for target spectator neutrons - hfn = "$ALICE_ROOT/ZDC/ZNASignal.root"; + + TNtuple* zdcSignal=0x0; + + Float_t sqrtS = 2*fBeamEnergy; + // + if(TMath::Abs(sqrtS-5500) < 100.){ + specSignalFile->cd("energy5500"); + // + if(SpecType == 1) { // --- Signal for projectile spectator neutrons + specSignalFile->GetObject("energy5500/ZNCSignal;1",zdcSignal); + if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZNCSignal from SpectatorSignal.root file"); + } + else if(SpecType == 2) { // --- Signal for projectile spectator protons + specSignalFile->GetObject("energy5500/ZPCSignal;1",zdcSignal); + if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZPCSignal from SpectatorSignal.root file"); + } + else if(SpecType == 3) { // --- Signal for target spectator neutrons + specSignalFile->GetObject("energy5500/ZNASignal;1",zdcSignal); + if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZNASignal from SpectatorSignal.root file"); + } + else if(SpecType == 4) { // --- Signal for target spectator protons + specSignalFile->GetObject("energy5500/ZPASignal;1",zdcSignal); + if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZPASignal from SpectatorSignal.root file"); + } } - else if(SpecType == 4) { // --- Signal for target spectator protons - hfn = "$ALICE_ROOT/ZDC/ZPASignal.root"; + else if(TMath::Abs(sqrtS-2760) < 100.){ + specSignalFile->cd("energy2760"); + // + if(SpecType == 1) { // --- Signal for projectile spectator neutrons + specSignalFile->GetObject("energy2760/ZNCSignal;1",zdcSignal); + if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZNCSignal from SpectatorSignal.root file"); + } + else if(SpecType == 2) { // --- Signal for projectile spectator protons + specSignalFile->GetObject("energy2760/ZPCSignal;1",zdcSignal); + if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZPCSignal from SpectatorSignal.root file"); + } + else if(SpecType == 3) { // --- Signal for target spectator neutrons + specSignalFile->GetObject("energy2760/ZNASignal;1",zdcSignal); + if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZNASignal from SpectatorSignal.root file"); + } + else if(SpecType == 4) { // --- Signal for target spectator protons + specSignalFile->GetObject("energy2760/ZPASignal;1",zdcSignal); + if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZPASignal from SpectatorSignal.root file"); + } } - TFile* file = TFile::Open(hfn); - if(!file || !file->IsOpen()) { - AliError((Form(" Opening file %s failed\n",hfn.Data()))); + if(!zdcSignal){ + printf("\n No spectator signal available for ZDC digitization\n"); return; - } - - TNtuple* zdcSignal = (TNtuple*) file->Get("ZDCSignal"); + } + Int_t nentries = (Int_t) zdcSignal->GetEntries(); Float_t *entry; @@ -507,8 +578,8 @@ void AliZDCDigitizer::SpectatorSignal(Int_t SpecType, Int_t numEvents, } }while(ievClose(); - delete file; + specSignalFile->Close(); + delete specSignalFile; } @@ -518,7 +589,9 @@ Int_t AliZDCDigitizer::Phe2ADCch(Int_t Det, Int_t Quad, Float_t Light, { // Evaluation of the ADC channel corresponding to the light yield Light Int_t vADCch = (Int_t) (Light * fPMGain[Det-1][Quad] * fADCRes[Res]); - //printf("\t Phe2ADCch -> det %d quad %d - phe %.0f ADC %d\n", Det,Quad,Light,vADCch); + // Ch. debug + //printf("\t Phe2ADCch -> det %d quad %d - PMGain[%d][%d] %1.0f phe %1.2f ADC %d\n", + // Det,Quad,Det-1,Quad,fPMGain[Det-1][Quad],Light,vADCch); return vADCch; } @@ -598,32 +671,3 @@ AliZDCPedestals* AliZDCDigitizer::GetPedData() const return calibdata; } -//_____________________________________________________________________________ -AliZDCEnCalib* AliZDCDigitizer::GetEnCalibData() const -{ - - // Getting calibration object for ZDC set - - AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/EnergyCalib"); - if(!entry) AliFatal("No calibration data loaded!"); - - AliZDCEnCalib *calibdata = dynamic_cast (entry->GetObject()); - if(!calibdata) AliFatal("Wrong calibration object in calibration file!"); - - return calibdata; -} - -//_____________________________________________________________________________ -AliZDCTowerCalib* AliZDCDigitizer::GetTowCalibData() const -{ - - // Getting calibration object for ZDC set - - AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/TowerCalib"); - if(!entry) AliFatal("No calibration data loaded!"); - - AliZDCTowerCalib *calibdata = dynamic_cast (entry->GetObject()); - if(!calibdata) AliFatal("Wrong calibration object in calibration file!"); - - return calibdata; -}