X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=ITS%2FAliITSsimulationSPD.cxx;h=ac65a13e6345d7eded82390dbd4553d5ab3e9b6e;hb=9618bdc71f5513c7ebd235a625c31b4ed069d1b5;hp=e0539444ef0ba5cfc788c6252531cd125b4380e3;hpb=cf6105a456c09f6fe863bd72c84b177a1c1f0be4;p=u%2Fmrichter%2FAliRoot.git diff --git a/ITS/AliITSsimulationSPD.cxx b/ITS/AliITSsimulationSPD.cxx index e0539444ef0..ac65a13e634 100644 --- a/ITS/AliITSsimulationSPD.cxx +++ b/ITS/AliITSsimulationSPD.cxx @@ -1,716 +1,1030 @@ /************************************************************************** - * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * - * * - * Author: The ALICE Off-line Project. * - * Contributors are mentioned in the code where appropriate. * - * * - * Permission to use, copy, modify and distribute this software and its * - * documentation strictly for non-commercial purposes is hereby granted * - * without fee, provided that the above copyright notice appears in all * - * copies and that both the copyright notice and this permission notice * - * appear in the supporting documentation. The authors make no claims * - * about the suitability of this software for any purpose. It is * - * provided "as is" without express or implied warranty. * - **************************************************************************/ +* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * +* * +* Author: The ALICE Off-line Project. * +* Contributors are mentioned in the code where appropriate. * +* * +* Permission to use, copy, modify and distribute this software and its * +* documentation strictly for non-commercial purposes is hereby granted * +* without fee, provided that the above copyright notice appears in all * +* copies and that both the copyright notice and this permission notice * +* appear in the supporting documentation. The authors make no claims * +* about the suitability of this software for any purpose. It is * +* provided "as is" without express or implied warranty. * +**************************************************************************/ /* -$Log$ -Revision 1.16 2002/06/19 16:02:22 hristov -Division by zero corrected - -Revision 1.15 2002/03/15 17:32:14 nilsen -Reintroduced SDigitization, and Digitization from SDigits, along with -functions InitSimulationModule, and FinishSDigitizModule. - -Revision 1.14 2001/11/23 13:04:07 barbera -Some protection added in case of high multiplicity - -Revision 1.13 2001/11/13 11:13:24 barbera -A protection against tracks with the same entrance and exit has been made more strict - -Revision 1.12 2001/10/04 22:44:31 nilsen -Major changes in supppor of PreDigits (SDigits). Changes made with will make -it easier to suppor expected changes in AliITSHit class. Added use of new -class AliITSpList. Both SPD and SDD have added effects of Dead Channels. Both -of these will require addtional work as data bases of detectors and the like -are developed. - +$Id$ */ -#include -#include + +#include +#include #include -#include #include -#include - -#include "AliRun.h" #include "AliITS.h" +#include "AliITSdigitSPD.h" #include "AliITShit.h" -#include "AliITSdigit.h" #include "AliITSmodule.h" -#include "AliITSMapA2.h" #include "AliITSpList.h" -#include "AliITSsimulationSPD.h" -#include "AliITSsegmentation.h" -#include "AliITSresponse.h" +#include "AliITSCalibrationSPD.h" #include "AliITSsegmentationSPD.h" -#include "AliITSresponseSPD.h" +#include "AliITSsimulationSPD.h" +#include "AliLog.h" +#include "AliRun.h" +#include "AliMagF.h" +#include "AliMathBase.h" +//#define DEBUG +using std::endl; +using std::cout; ClassImp(AliITSsimulationSPD) //////////////////////////////////////////////////////////////////////// -// Version: 0 -// Written by Rocco Caliandro -// from a model developed with T. Virgili and R.A. Fini -// June 15 2000 +// Version: 1 +// Modified by D. Elia, G.E. Bruno, H. Tydesjo +// Fast diffusion code by Bjorn S. Nilsen +// March-April 2006 +// October 2007: GetCalibrationObjects() removed // -// AliITSsimulationSPD is the simulation of SPDs +// Version: 0 +// Written by Boris Batyunya +// December 20 1999 // -//______________________________________________________________________ -AliITSsimulationSPD::AliITSsimulationSPD(){ - // Default constructor - - fResponse = 0; - fSegmentation = 0; - fHis = 0; - fMapA2 = 0; +// +// AliITSsimulationSPD is to do the simulation of SPDs. +// +//////////////////////////////////////////////////////////////////////// -/* - fThresh = 0.; - fSigma = 0.; - fCouplCol = 0.; - fCouplRow = 0.; */ +//______________________________________________________________________ +AliITSsimulationSPD::AliITSsimulationSPD(): +AliITSsimulation(), +fHis(0), +fSPDname(), +fCoupling(), +fLorentz(kFALSE), +fTanLorAng(0), +fStrobe(kTRUE), +fStrobeLenght(4), +fStrobePhase(-12.5e-9){ + // Default constructor. + // Inputs: + // none. + // Outputs: + // none. + // Return: + // A default constructed AliITSsimulationSPD class. + + AliDebug(1,Form("Calling default constructor")); +// Init(); } //______________________________________________________________________ -AliITSsimulationSPD::AliITSsimulationSPD(AliITSsegmentation *seg, - AliITSresponse *resp) { - // Standard constructor - - fResponse = 0; - fSegmentation = 0; - fHis = 0; - fMapA2 = 0; - -/* - fThresh = 0.; - fSigma = 0.; - fCouplCol = 0.; - fCouplRow = 0.*/ - Init((AliITSsegmentationSPD*)seg,(AliITSresponseSPD*)resp); +AliITSsimulationSPD::AliITSsimulationSPD(AliITSDetTypeSim *dettyp): +AliITSsimulation(dettyp), +fHis(0), +fSPDname(), +fCoupling(), +fLorentz(kFALSE), +fTanLorAng(0), +fStrobe(kTRUE), +fStrobeLenght(4), +fStrobePhase(-12.5e-9){ + // standard constructor + // Inputs: + // AliITSsegmentation *seg A pointer to the segmentation class + // to be used for this simulation + // AliITSCalibration *resp A pointer to the responce class to + // be used for this simulation + // Outputs: + // none. + // Return: + // A default constructed AliITSsimulationSPD class. + + AliDebug(1,Form("Calling standard constructor ")); + Init(); } //______________________________________________________________________ -void AliITSsimulationSPD::Init(AliITSsegmentationSPD *seg, - AliITSresponseSPD *resp) { - // Initilizes the variables of AliITSsimulation SPD. - - fHis = 0; - fResponse = resp; - fSegmentation = seg; - fMapA2 = new AliITSMapA2(fSegmentation); - fpList = new AliITSpList(GetNPixelsZ()+1,GetNPixelsX()+1); -/* - fResponse->Thresholds(fThresh,fSigma); - fResponse->GetNoiseParam(fCouplCol,fCouplRow); - fNPixelsZ = fSegmentation->Npz(); - fNPixelsX = fSegmentation->Npx(); -*/ +void AliITSsimulationSPD::Init(){ + // Initilization + // Inputs: + // none. + // Outputs: + // none. + // Return: + // none. + const Double_t kmictocm = 1.0e-4; // convert microns to cm. + + SetModuleNumber(0); + SetEventNumber(0); + SetMap(new AliITSpList(GetNPixelsZ(),GetNPixelsX())); + AliITSSimuParam* simpar = fDetType->GetSimuParam(); + AliITSsegmentationSPD* seg = (AliITSsegmentationSPD*)GetSegmentationModel(0); + Double_t bias = simpar->GetSPDBiasVoltage(); +// cout << "Bias Voltage --> " << bias << endl; // dom + simpar->SetDistanceOverVoltage(kmictocm*seg->Dy(),bias); +// set kind of coupling ("old" or "new") + char opt[20]; + simpar->GetSPDCouplingOption(opt); + char *old = strstr(opt,"old"); + if (old) { + fCoupling=2; + } else { + fCoupling=1; + } // end if + SetLorentzDrift(simpar->GetSPDLorentzDrift()); + if (fLorentz) SetTanLorAngle(simpar->GetSPDLorentzHoleWeight()); + //SetStrobeGeneration(kFALSE); + if (fStrobe) GenerateStrobePhase(); } //______________________________________________________________________ -AliITSsimulationSPD::~AliITSsimulationSPD() { - // destructor - - delete fMapA2; -// delete fpList; - - if (fHis) { - fHis->Delete(); - delete fHis; - } // end if +Bool_t AliITSsimulationSPD::SetTanLorAngle(Double_t WeightHole) { + // This function set the Tangent of the Lorentz angle. + // A weighted average is used for electrons and holes + // Input: Double_t WeightHole: wheight for hole: it should be in the range [0,1] + // output: Bool_t : kTRUE in case of success + // + if(!fDetType) { + AliError("AliITSsimulationSPD::SetTanLorAngle: AliITSDetTypeSim* fDetType not set "); + return kFALSE;} + if(WeightHole<0) { + WeightHole=0.; + AliWarning("AliITSsimulationSPD::SetTanLorAngle: You have asked for negative Hole weight"); + AliWarning("AliITSsimulationSPD::SetTanLorAngle: I'm going to use only electrons"); + } + if(WeightHole>1) { + WeightHole=1.; + AliWarning("AliITSsimulationSPD::SetTanLorAngle: You have asked for weight > 1"); + AliWarning("AliITSsimulationSPD::SetTanLorAngle: I'm going to use only holes"); + } + Double_t WeightEle=1.-WeightHole; + AliITSSimuParam* simpar = fDetType->GetSimuParam(); + AliMagF* fld = (AliMagF*)TGeoGlobalMagField::Instance()->GetField(); + if (!fld) AliFatal("The field is not initialized"); + Double_t bz = fld->SolenoidField(); + fTanLorAng = TMath::Tan(WeightHole*simpar->LorentzAngleHole(bz) + + WeightEle*simpar->LorentzAngleElectron(bz)); + return kTRUE; } //______________________________________________________________________ -AliITSsimulationSPD::AliITSsimulationSPD(const AliITSsimulationSPD &source){ - // Copy Constructor - - if(&source == this) return; - - this->fMapA2 = source.fMapA2; - this->fHis = source.fHis; -/* - this->fThresh = source.fThresh; - this->fSigma = source.fSigma; - this->fCouplCol = source.fCouplCol; - this->fCouplRow = source.fCouplRow; - this->fNPixelsX = source.fNPixelsX; - this->fNPixelsZ = source.fNPixelsZ; -*/ - return; +AliITSsimulationSPD::~AliITSsimulationSPD(){ + // destructor + // Inputs: + // none. + // Outputs: + // none. + // Return: + // none. + + if (fHis) { + fHis->Delete(); + delete fHis; + } // end if fHis } //______________________________________________________________________ -AliITSsimulationSPD& AliITSsimulationSPD::operator=(const AliITSsimulationSPD - &source) { - // Assignment operator +AliITSsimulationSPD::AliITSsimulationSPD(const + AliITSsimulationSPD + &s) : AliITSsimulation(s), +fHis(s.fHis), +fSPDname(s.fSPDname), +fCoupling(s.fCoupling), +fLorentz(s.fLorentz), +fTanLorAng(s.fTanLorAng), +fStrobe(s.fStrobe), +fStrobeLenght(s.fStrobeLenght), +fStrobePhase(s.fStrobePhase){ + // Copy Constructor + // Inputs: + // AliITSsimulationSPD &s The original class for which + // this class is a copy of + // Outputs: + // none. + // Return: - if(&source == this) return *this; - - this->fMapA2 = source.fMapA2; - this->fHis = source.fHis; -/* - this->fThresh = source.fThresh; - this->fSigma = source.fSigma; - this->fCouplCol = source.fCouplCol; - this->fCouplRow = source.fCouplRow; - this->fNPixelsX = source.fNPixelsX; - this->fNPixelsZ = source.fNPixelsZ; -*/ - return *this; -} +} //______________________________________________________________________ -void AliITSsimulationSPD::InitSimulationModule(Int_t module,Int_t event){ - // Creates maps to build the list of tracks for each sumable digit - // Inputs: - // Int_t module // Module number to be simulated - // Int_t event // Event number to be simulated - // Outputs: - // none. - // Return - // none. - - fModule = module; - fEvent = event; - fMapA2->ClearMap(); - fpList->ClearMap(); +AliITSsimulationSPD& AliITSsimulationSPD::operator=(const + AliITSsimulationSPD &s){ + // Assignment operator + // Inputs: + // AliITSsimulationSPD &s The original class for which + // this class is a copy of + // Outputs: + // none. + // Return: + + if(&s == this) return *this; + this->fHis = s.fHis; + fCoupling = s.fCoupling; + fSPDname = s.fSPDname; + fLorentz = s.fLorentz; + fTanLorAng = s.fTanLorAng; + fStrobe = s.fStrobe; + fStrobeLenght = s.fStrobeLenght; + fStrobePhase = s.fStrobePhase; + return *this; } +/* //______________________________________________________________________ -void AliITSsimulationSPD::FinishSDigitiseModule(){ - // Does the Sdigits to Digits work - // Inputs: - // none. - // Outputs: - // none. - // Return: - // none. - - SDigitsToDigits(fModule,fpList); +AliITSsimulation& AliITSsimulationSPD::operator=(const + AliITSsimulation &s){ + // Assignment operator + // Inputs: + // AliITSsimulationSPD &s The original class for which + // this class is a copy of + // Outputs: + // none. + // Return: + + if(&s == this) return *this; + Error("AliITSsimulationSPD","Not allowed to make a = with " + "AliITSsimulationSPD","Using default creater instead"); + + return *this; } +*/ //______________________________________________________________________ -void AliITSsimulationSPD::SDigitiseModule(AliITSmodule *mod, Int_t dummy0, - Int_t dummy1) { - // Sum digitize module - if (!(mod->GetNhits())) return; // if module has no hits then no Sdigits. - Int_t number = 10000; - Int_t *frowpixel = new Int_t[number]; - Int_t *fcolpixel = new Int_t[number]; - Double_t *fenepixel = new Double_t[number]; - - fModule = mod->GetIndex(); - - // Array of pointers to store the track index of the digits - // leave +1, otherwise pList crashes when col=256, row=192 - - HitsToAnalogDigits(mod,frowpixel,fcolpixel,fenepixel,fpList); - - WriteSDigits(fpList); - - // clean memory - delete[] frowpixel; - delete[] fcolpixel; - delete[] fenepixel; - fMapA2->ClearMap(); - fpList->ClearMap(); +void AliITSsimulationSPD::InitSimulationModule(Int_t module, Int_t event){ + // This function creates maps to build the list of tracks for each + // summable digit. Inputs defined by base class. + // Inputs: + // Int_t module // Module number to be simulated + // Int_t event // Event number to be simulated + // Outputs: + // none + // Returns: + // none + + AliDebug(1,Form("(module=%d,event=%d)",module,event)); + const Double_t kmictocm = 1.0e-4; // convert microns to cm. + AliITSSimuParam* simpar = fDetType->GetSimuParam(); + AliITSsegmentationSPD* seg = (AliITSsegmentationSPD*)GetSegmentationModel(0); + SetModuleNumber(module); + SetEventNumber(event); + simpar->SetDistanceOverVoltage(kmictocm*seg->Dy(),simpar->GetSPDBiasVoltage(module)); + ClearMap(); +} +//_____________________________________________________________________ +void AliITSsimulationSPD::SDigitiseModule(AliITSmodule *mod,Int_t, + Int_t event){ + // This function begins the work of creating S-Digits. Inputs defined + // by base class. + // Inputs: + // AliITSmodule *mod // module + // Int_t // not used + // Int_t event // Event number + // Outputs: + // none + // Return: + // test // test returns kTRUE if the module contained hits + // // test returns kFALSE if it did not contain hits + + AliDebug(1,Form("(mod=%p, ,event=%d)",mod,event)); + if(!(mod->GetNhits())){ + AliDebug(1,Form("In event %d module %d there are %d hits returning.", + event, mod->GetIndex(),mod->GetNhits())); + return;// if module has no hits don't create Sdigits + } // end if + SetModuleNumber(mod->GetIndex()); + if (fStrobe) if(event != GetEventNumber()) GenerateStrobePhase(); + SetEventNumber(event); + InitSimulationModule( GetModuleNumber() , event ); + // HitToSDigit(mod); + HitToSDigitFast(mod); + if (fDetType->GetSimuParam()->GetSPDAddNoisyFlag()) AddNoisyPixels(); + if (fDetType->GetSimuParam()->GetSPDRemoveDeadFlag()) RemoveDeadPixels(); + +// cout << "After Remove in SDigitiseModule !!!!!" << endl; // dom +// cout << "Module " << mod->GetIndex() << " Event " << event << endl; // dom + WriteSDigits(); + ClearMap(); } //______________________________________________________________________ -void AliITSsimulationSPD::DigitiseModule(AliITSmodule *mod, Int_t dummy0, - Int_t dummy1) { - // digitize module. Also need to digitize modules with only noise. - - Int_t number = 10000; - Int_t *frowpixel = new Int_t[number]; - Int_t *fcolpixel = new Int_t[number]; - Double_t *fenepixel = new Double_t[number]; - - // Array of pointers to store the track index of the digits - // leave +1, otherwise pList crashes when col=256, row=192 - fModule = mod->GetIndex(); - // noise setting - SetFluctuations(fpList,fModule); - - HitsToAnalogDigits(mod,frowpixel,fcolpixel,fenepixel,fpList); - - // apply mask to SPD module - SetMask(); - - CreateDigit(fModule,fpList); - - // clean memory - delete[] frowpixel; - delete[] fcolpixel; - delete[] fenepixel; - fMapA2->ClearMap(); - fpList->ClearMap(); +void AliITSsimulationSPD::WriteSDigits(){ + // This function adds each S-Digit to pList + // Inputs: + // none. + // Outputs: + // none. + // Return: + // none + Int_t ix, nix, iz, niz; + static AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS"); + + AliDebug(1,Form("Writing SDigits for module %d",GetModuleNumber())); +// cout << "WriteSDigits for module " << GetModuleNumber() << endl; // dom + GetMap()->GetMaxMapIndex(niz, nix); + for(iz=0; izGetSignalOnly(iz,ix)>0.0){ +// cout << " Signal gt 0 iz ix " << iz << ix << " Module " << GetModuleNumber() << endl; // dom + aliITS->AddSumDigit(*(GetMap()->GetpListItem(iz,ix))); + if(AliDebugLevel()>0) { + AliDebug(1,Form("%d, %d",iz,ix)); + cout << *(GetMap()->GetpListItem(iz,ix)) << endl; + } // end if GetDebug + } // end if GetMap()->GetSignalOnly(iz,ix)>0.0 + } // end for iz,ix + return; } //______________________________________________________________________ -void AliITSsimulationSPD::SDigitsToDigits(Int_t module,AliITSpList *pList) { - // sum digits to Digits. -// cout << "Entering AliITSsimulatinSPD::SDigitsToDigits for module="; -// cout << module << endl; - fModule = module; - - // noise setting - SetFluctuations(pList,module); - - fMapA2->ClearMap(); // since noise is in pList aready. Zero Map so that - // noise is not doubled when calling FillMapFrompList. - - FillMapFrompList(pList); - - // apply mask to SPD module - SetMask(); - - CreateDigit(module,pList); - - fMapA2->ClearMap(); - pList->ClearMap(); +void AliITSsimulationSPD::FinishSDigitiseModule(){ + // This function calls SDigitsToDigits which creates Digits from SDigits + // Inputs: + // none + // Outputs: + // none + // Return + // none + + AliDebug(1,"()"); +// cout << "FinishSDigitiseModule for module " << GetModuleNumber() << endl; // dom + FrompListToDigits(); // Charge To Signal both adds noise and + ClearMap(); + return; } //______________________________________________________________________ -void AliITSsimulationSPD::UpdateMapSignal(Int_t row,Int_t col,Int_t trk, - Int_t hit,Int_t mod,Double_t ene, - AliITSpList *pList) { - // updates the Map of signal, adding the energy (ene) released by - // the current track - - fMapA2->AddSignal(row,col,ene); - pList->AddSignal(row,col,trk,hit,mod,ene); +void AliITSsimulationSPD::DigitiseModule(AliITSmodule *mod,Int_t, + Int_t event){ + // This function creates Digits straight from the hits and then adds + // electronic noise to the digits before adding them to pList + // Each of the input variables is passed along to HitToSDigit + // Inputs: + // AliITSmodule *mod module + // Int_t Dummy. + // Int_t Dummy + // Outputs: + // none. + // Return: + // none. + + if (fStrobe) if(event != GetEventNumber()) GenerateStrobePhase(); + AliDebug(1,Form("(mod=%p,,0)",mod)); + // HitToSDigit(mod); + InitSimulationModule( mod->GetIndex(), event ); + HitToSDigitFast(mod); + + if (fDetType->GetSimuParam()->GetSPDAddNoisyFlag()) AddNoisyPixels(); + if (fDetType->GetSimuParam()->GetSPDRemoveDeadFlag()) RemoveDeadPixels(); +// cout << "After Remove in DigitiseModule in module " << mod->GetIndex() << endl; // dom + FrompListToDigits(); + ClearMap(); } //______________________________________________________________________ -void AliITSsimulationSPD::UpdateMapNoise(Int_t row,Int_t col,Int_t mod, - Double_t ene,AliITSpList *pList) { - // updates the Map of noise, adding the energy (ene) give my noise - - fMapA2->AddSignal(row,col,ene); - pList->AddNoise(row,col,mod,ene); +void AliITSsimulationSPD::HitToSDigit(AliITSmodule *mod){ + // Does the charge distributions using Gaussian diffusion charge charing. + // Inputs: + // AliITSmodule *mod Pointer to this module + // Output: + // none. + // Return: + // none. + const Double_t kmictocm = 1.0e-4; // convert microns to cm. + const Double_t kBunchLenght = 25e-9; // LHC clock + TObjArray *hits = mod->GetHits(); + Int_t nhits = hits->GetEntriesFast(); + Int_t h,ix,iz,i; + Int_t idtrack; + Double_t x0=0.0,x1=0.0,y0=0.0,y1=0.0,z0=0.0,z1=0.0,de=0.0,ld=0.0; + Double_t x,y,z,t,tp,st,dt=0.2,el,sig,sigx,sigz,fda; + AliITSsegmentationSPD* seg = (AliITSsegmentationSPD*)GetSegmentationModel(0); + AliITSSimuParam *simpar = fDetType->GetSimuParam(); + Double_t thick = 0.5*kmictocm*seg->Dy(); // Half Thickness + simpar->GetSPDSigmaDiffusionAsymmetry(fda); // + + AliDebug(1,Form("(mod=%p) fCoupling=%d",mod,fCoupling)); + if(nhits<=0) return; + for(h=0;h0) { + AliDebug(1,Form("Hits, %d", h)); + cout << *(mod->GetHit(h)) << endl; + } // end if GetDebug + // Check if the hit is inside readout window + if (fStrobe) + if ((mod->GetHit(h)->GetTOF() < fStrobePhase) || + (mod->GetHit(h)->GetTOF() > (fStrobePhase+(Double_t)fStrobeLenght*kBunchLenght))) continue; + if(!mod->LineSegmentL(h,x0,x1,y0,y1,z0,z1,de,idtrack)) continue; + st = TMath::Sqrt(x1*x1+y1*y1+z1*z1); + if(st>0.0){ + st = (Double_t)((Int_t)(st/kmictocm)); // number of microns + if(st<=1.0) st = 1.0; + dt = 1.0/st; + for(t=0.0;t<1.0;t+=dt){ // Integrate over t + tp = t+0.5*dt; + x = x0+x1*tp; + y = y0+y1*tp; + z = z0+z1*tp; + if(!(seg->LocalToDet(x,z,ix,iz))) continue; // outside + //el = res->GeVToCharge((Double_t)(dt*de)); + el = dt * de / simpar->GetGeVToCharge(); + if(GetDebug(1)){ + if(el<=0.0) cout<<"el="<GetNhits(); - for (hitpos=0;hitposGetHits(); + Int_t nhits = hits->GetEntriesFast(); + Int_t h,ix,iz,i; + Int_t idtrack; + Double_t x0=0.0,x1=0.0,y0=0.0,y1=0.0,z0=0.0,z1=0.0,de=0.0,ld=0.0; + Double_t x,y,z,t,st,el,sig,sigx,sigz,fda; + AliITSsegmentationSPD* seg = (AliITSsegmentationSPD*)GetSegmentationModel(0); + AliITSSimuParam* simpar = fDetType->GetSimuParam(); + Double_t thick = 0.5*kmictocm*seg->Dy(); // Half thickness + simpar->GetSPDSigmaDiffusionAsymmetry(fda); +// cout << "Half Thickness " << thick << endl; // dom +// cout << "Diffusion asymm " << fda << endl; // dom + + AliDebug(1,Form("(mod=%p) fCoupling=%d",mod,fCoupling)); + if(nhits<=0) return; + for(h=0;h0) { + AliDebug(1,Form("Hits, %d", h)); + cout << *(mod->GetHit(h)) << endl; + } // end if GetDebug + // Check if the hit is inside readout window + if (fStrobe) + if ((mod->GetHit(h)->GetTOF() < fStrobePhase) || + (mod->GetHit(h)->GetTOF() > (fStrobePhase+(Double_t)fStrobeLenght*kBunchLenght))) continue; + if(!mod->LineSegmentL(h,x0,x1,y0,y1,z0,z1,de,idtrack)) continue; + st = TMath::Sqrt(x1*x1+y1*y1+z1*z1); + if(st>0.0) for(i=0;iLocalToDet(x,z,ix,iz))) continue; // outside + el = kwi[i]*de/simpar->GetGeVToCharge(); + if(GetDebug(1)){ + if(el<=0.0) cout<<"el="<LocalToDet(x,z,ix,iz))) continue; // outside + el = de / simpar->GetGeVToCharge(); + sig = simpar->SigmaDiffusion1D(TMath::Abs(thick + y)); + sigx=sig; + sigz=sig*fda; + if (fLorentz) ld=(y+thick)*fTanLorAng; + SpreadChargeAsym(x,z,ix,iz,el,sigx,sigz,ld,idtrack,h); + } // end if st>0.0 + + } // Loop over all hits h + + // Coupling + switch (fCoupling) { + default: + break; + case 1: // case 3: + for(i=0;iGetEntries();i++) + if(GetMap()->GetpListItem(i)==0) continue; + else{ + GetMap()->GetMapIndex(GetMap()->GetpListItem(i)->GetIndex(),iz,ix); + SetCoupling(iz,ix); + } // end for i + break; + case 2: // case 4: + for(i=0;iGetEntries();i++) + if(GetMap()->GetpListItem(i)==0) continue; + else{ + GetMap()->GetMapIndex(GetMap()->GetpListItem(i)->GetIndex(),iz,ix); + SetCouplingOld(iz,ix); + } // end for i + break; + } // end switch + if(GetDebug(2))Info("HitToSDigit","Finished fCoupling=%d",fCoupling); } //______________________________________________________________________ -void AliITSsimulationSPD::HitToDigit(AliITSmodule *mod,Int_t hitpos, - Int_t *frowpixel,Int_t *fcolpixel, - Double_t *fenepixel,AliITSpList *pList) { - // Steering function to determine the digits associated to a given - // hit (hitpos) - // The digits are created by charge sharing (ChargeSharing) and by - // capacitive coupling (SetCoupling). At all the created digits is - // associated the track number of the hit (ntrack) - Double_t x1l=0.0,y1l=0.0,z1l=0.0,x2l=0.0,y2l=0.0,z2l=0.0; - Int_t r1,r2,c1,c2,row,col,npixel = 0; - Int_t ntrack; - Double_t ene=0.0,etot=0.0; - const Float_t kconv = 10000.; // cm -> microns - const Float_t kconv1= 0.277e9; // GeV -> electrons equivalent - - if(!(mod->LineSegmentL(hitpos,x1l,x2l,y1l,y2l,z1l,z2l,etot,ntrack)))return; - - x2l += x1l; y2l += y1l; z2l += z1l; // Convert to ending coordinate. - // positions shifted and converted in microns - x1l = x1l*kconv + fSegmentation->Dx()/2.; - z1l = z1l*kconv + fSegmentation->Dz()/2.; - // positions shifted and converted in microns - x2l = x2l*kconv + fSegmentation->Dx()/2.; - z2l = z2l*kconv + fSegmentation->Dz()/2.; - etot *= kconv1; // convert from GeV to electrons equivalent. - Int_t module = mod->GetIndex(); - - // to account for the effective sensitive area - // introduced in geometry - if (z1l<0 || z1l>fSegmentation->Dz()) return; - if (z2l<0 || z2l>fSegmentation->Dz()) return; - if (x1l<0 || x1l>fSegmentation->Dx()) return; - if (x2l<0 || x2l>fSegmentation->Dx()) return; - - //Get the col and row number starting from 1 - // the x direction is not inverted for the second layer!!! - fSegmentation->GetPadIxz(x1l, z1l, c1, r1); - fSegmentation->GetPadIxz(x2l, z2l, c2, r2); - - // to account for unexpected equal entrance and - // exit coordinates - if (x1l==x2l) x2l=x2l+x2l*0.1; - if (z1l==z2l) z2l=z2l+z2l*0.1; - - if ((r1==r2) && (c1==c2)){ - // no charge sharing - npixel = 1; - frowpixel[npixel-1] = r1; - fcolpixel[npixel-1] = c1; - fenepixel[npixel-1] = etot; - } else { - // charge sharing - ChargeSharing(x1l,z1l,x2l,z2l,c1,r1,c2,r2,etot, - npixel,frowpixel,fcolpixel,fenepixel); - } // end if r1==r2 && c1==c2. - - for (Int_t npix=0;npixAddSignal(iz0,ix0,t,hi,GetModuleNumber(),el); + if(GetDebug(2)){ + cout << "sig<=0.0=" << sig << endl; + } // end if GetDebug + return; + } // end if + sp = 1.0/(sig*kRoot2); + if(GetDebug(2)){ + cout << "sig=" << sig << " sp=" << sp << endl; + } // end if GetDebug + ixs = TMath::Max(-knx+ix0,0); + ixe = TMath::Min(knx+ix0,seg->Npx()-1); + izs = TMath::Max(-knz+iz0,0); + ize = TMath::Min(knz+iz0,seg->Npz()-1); + for(ix=ixs;ix<=ixe;ix++) for(iz=izs;iz<=ize;iz++){ + seg->DetToLocal(ix,iz,x,z); // pixel center + x1 = x; + z1 = z; + x2 = x1 + 0.5*kmictocm*seg->Dpx(ix); // Upper + x1 -= 0.5*kmictocm*seg->Dpx(ix); // Lower + z2 = z1 + 0.5*kmictocm*seg->Dpz(iz); // Upper + z1 -= 0.5*kmictocm*seg->Dpz(iz); // Lower + x1 -= x0+ld; // Distance from where track traveled (taking into account the Lorentz drift) + x2 -= x0+ld; // Distance from where track traveled (taking into account the Lorentz drift) + z1 -= z0; // Distance from where track traveled + z2 -= z0; // Distance from where track traveled + s = 0.25; // Correction based on definision of Erfc + s *= AliMathBase::ErfcFast(sp*x1) - AliMathBase::ErfcFast(sp*x2); + if(GetDebug(3)){ + cout <<"el="< - -
- - . - - 
-    */
-    //End_Html
-    Float_t xa,za,xb,zb,dx,dz,dtot,dm,refr,refm,refc;
-    Float_t refn=0.;
-    Float_t arefm, arefr, arefn, arefc, azb, az2l, axb, ax2l;
-    Int_t   dirx,dirz,rb,cb;
-    Int_t flag,flagrow,flagcol;
-    Double_t epar;
-
-    npixel = 0;
-    xa     = x1l;
-    za     = z1l;
-//    dx     = x1l-x2l;
-//    dz     = z1l-z2l;
-    dx     = x2l-x1l;
-    dz     = z2l-z1l;
-    dtot   = TMath::Sqrt((dx*dx)+(dz*dz));   
-    if (dtot==0.0) dtot = 0.01;
-    dirx   = (Int_t) TMath::Sign((Float_t)1,dx);
-    dirz   = (Int_t) TMath::Sign((Float_t)1,dz);
-
-    // calculate the x coordinate of  the pixel in the next column    
-    // and the z coordinate of  the pixel in the next row
-    Float_t xpos, zpos;
-
-    fSegmentation->GetPadCxz(c1, r1-1, xpos, zpos); 
-
-    Float_t xsize = fSegmentation->Dpx(0);
-    Float_t zsize = fSegmentation->Dpz(r1-1);
-    
-    if (dirx == 1) refr = xpos+xsize/2.;
-    else refr = xpos-xsize/2.;
-
-    if (dirz == 1) refn = zpos+zsize/2.;
-    else refn = zpos-zsize/2.;
-
-    flag = 0;
-    flagrow = 0;
-    flagcol = 0;
-    do{
-	// calculate the x coordinate of the intersection with the pixel
-	// in the next cell in row  direction
-      if(dz!=0)
-        refm = dx*((refn - z1l)/dz) + x1l;
-      else
-        refm = refr+dirx*xsize;
-   
-	// calculate the z coordinate of the intersection with the pixel
-	// in the next cell in column direction
-      if (dx!=0)
-        refc = dz*((refr - x1l)/dx) + z1l;
-      else
-        refc = refn+dirz*zsize;
-
-	arefm = refm * dirx;
-	arefr = refr * dirx;
-	arefn = refn * dirz;
-	arefc = refc * dirz;
-
-	if ((arefm < arefr) && (arefn < arefc)){
-	    // the track goes in the pixel in the next cell in row direction
-	    xb = refm;
-	    zb = refn;
-	    cb = c1;
-	    rb = r1 + dirz;
-	    azb = zb * dirz;
-	    az2l = z2l * dirz;
-	    if (rb == r2) flagrow=1;
-	    if (azb > az2l) {
-	        zb = z2l;
-	        xb = x2l;
-	    } // end if
-	    // shift to the pixel in the next cell in row direction
-	    Float_t zsizeNext = fSegmentation->Dpz(rb-1);
-	    //to account for cell at the borders of the detector
-	    if(zsizeNext==0) zsizeNext = zsize;
-	    refn += zsizeNext*dirz;
-	}else {
-	    // the track goes in the pixel in the next cell in column direction
-	    xb = refr;
-	    zb = refc;
-	    cb = c1 + dirx;
-	    rb = r1;
-	    axb = xb * dirx;
-	    ax2l = x2l * dirx;
-	    if (cb == c2) flagcol=1;
-	    if (axb > ax2l) {
-	        zb = z2l;
-	        xb = x2l;
-	    } // end ifaxb > ax2l
-
-	    // shift to the pixel in the next cell in column direction
-	    Float_t xsizeNext = fSegmentation->Dpx(cb-1);
-	    //to account for cell at the borders of the detector
-	    if(xsizeNext==0) xsizeNext = xsize;
-	    refr += xsizeNext*dirx;
-	} // end if (arefm < arefr) && (arefn < arefc)
-
-	//calculate the energy lost in the crossed pixel      
-	epar = TMath::Sqrt((xb-xa)*(xb-xa)+(zb-za)*(zb-za)); 
-	epar = etot*(epar/dtot);
-
-	//store row, column and energy lost in the crossed pixel
-	frowpixel[npixel] = r1;
-	fcolpixel[npixel] = c1;
-	fenepixel[npixel] = epar;
-	npixel++;
- 
-	// the exit point of the track is reached
-	if (epar == 0) flag = 1;
-	if ((r1 == r2) && (c1 == c2)) flag = 1;
-	if (flag!=1) {
-	    r1 = rb;
-	    c1 = cb;
-	    xa = xb;
-	    za = zb;
-	} // end if flag!=1
-    } while (flag==0);
+void AliITSsimulationSPD::SpreadChargeAsym(Double_t x0,Double_t z0,
+                                           Int_t ix0,Int_t iz0,
+                                           Double_t el,Double_t sigx,Double_t sigz,
+                                           Double_t ld,Int_t t,Int_t hi){
+   // Spreads the charge over neighboring cells. Assume charge is distributed
+   // as charge(x,z) = (el/2*pi*sigx*sigz)*exp(-arg)
+   // arg=((x-x0)*(x-x0)/2*sigx*sigx)+((z-z0*z-z0)/2*sigz*sigz)
+   // if fLorentz=kTRUE, then x0=x0+ld (Lorentz drift taken into account)
+   // Defined this way, the integral over all x and z is el.
+   // Inputs:
+   //    Double_t x0   x position of point where charge is liberated
+   //    Double_t z0   z position of point where charge is liberated
+   //    Int_t    ix0  row of cell corresponding to point x0
+   //    Int_t    iz0  columb of cell corresponding to point z0
+   //    Double_t el   number of electrons liberated in this step
+   //    Double_t sigx Sigma difusion along x for this step (y0 dependent)
+   //    Double_t sigz Sigma difusion along z for this step (y0 dependent)
+   //    Double_t ld   lorentz drift in x for this stip (y0 dependent)
+   //    Int_t    t    track number
+   //    Int_t    ti   hit track index number
+   //    Int_t    hi   hit "hit" index number
+   // Outputs:
+   //     none.
+   // Return:
+   //     none.
+   const Int_t knx = 3,knz = 2;
+   const Double_t kRoot2 = 1.414213562; // Sqrt(2).
+   const Double_t kmictocm = 1.0e-4; // convert microns to cm.
+   Int_t ix,iz,ixs,ixe,izs,ize;
+   Float_t x,z;
+   Double_t x1,x2,z1,z2,s,spx,spz;
+   AliITSsegmentationSPD* seg = (AliITSsegmentationSPD*)GetSegmentationModel(0);
+
+
+   if(GetDebug(4)) Info("SpreadChargeAsym","(x0=%e,z0=%e,ix0=%d,iz0=%d,el=%e,"
+                        "sigx=%e, sigz=%e, t=%d,i=%d)",x0,z0,ix0,iz0,el,sigx,sigz,t,hi);
+   if(sigx<=0.0 || sigz<=0.0) { // if sig<=0 No diffusion to simulate.
+       GetMap()->AddSignal(iz0,ix0,t,hi,GetModuleNumber(),el);
+       if(GetDebug(2)){
+           cout << "sigx<=0.0=" << sigx << endl;
+           cout << "sigz<=0.0=" << sigz << endl;
+       } // end if GetDebug
+       return;
+   } // end if
+   spx = 1.0/(sigx*kRoot2);     spz = 1.0/(sigz*kRoot2);
+   if(GetDebug(2)){
+       cout << "sigx=" << sigx << " spx=" << spx << endl;
+       cout << "sigz=" << sigz << " spz=" << spz << endl;
+   } // end if GetDebug
+   ixs = TMath::Max(-knx+ix0,0);
+   ixe = TMath::Min(knx+ix0,seg->Npx()-1);
+   izs = TMath::Max(-knz+iz0,0);
+   ize = TMath::Min(knz+iz0,seg->Npz()-1);
+   for(ix=ixs;ix<=ixe;ix++) for(iz=izs;iz<=ize;iz++){
+       seg->DetToLocal(ix,iz,x,z); // pixel center
+       x1  = x;
+       z1  = z;
+       x2  = x1 + 0.5*kmictocm*seg->Dpx(ix); // Upper
+       x1 -= 0.5*kmictocm*seg->Dpx(ix);  // Lower
+       z2  = z1 + 0.5*kmictocm*seg->Dpz(iz); // Upper
+       z1 -= 0.5*kmictocm*seg->Dpz(iz);  // Lower
+       x1 -= x0+ld; // Distance from where track traveled (taking into account the Lorentz drift)
+       x2 -= x0+ld; // Distance from where track traveled (taking into account the Lorentz drift)
+       z1 -= z0; // Distance from where track traveled
+       z2 -= z0; // Distance from where track traveled
+       s   = 0.25; // Correction based on definision of Erfc
+       s  *= AliMathBase::ErfcFast(spx*x1) - AliMathBase::ErfcFast(spx*x2);
+       if(GetDebug(3)){
+           cout <<"el="<
-
-
- - . - - 
-    */
-    //End_Html
-    Int_t j1,j2,flag=0;
-    Double_t pulse1,pulse2;
-    Float_t couplR=0.0,couplC=0.0;
-
-    GetCouplings(couplR,couplC);
-    j1 = row;
-    j2 = col;
-    pulse1 = fMapA2->GetSignal(row,col);
-    pulse2 = pulse1;
-    for (Int_t isign=-1;isign<=1;isign+=2){// loop in row direction
-	do{
-	    j1 += isign;
-	    pulse1 *= couplR;
-	    if ((j1<0) || (j1>GetNPixelsZ()-1) || (pulse1GetSignal(row,col);
-		j1 = row;
-		flag = 1;
-	    }else{
-		UpdateMapSignal(j1,col,ntrack,idhit,module,pulse1,pList);
-		flag = 0;
-	    } // end if
-	} while(flag == 0);
-	// loop in column direction
-      do{
-	  j2 += isign;
-	  pulse2 *= couplC;
-	  if ((j2<0) || (j2>(GetNPixelsX()-1)) || (pulse2GetSignal(row,col);
-	      j2 = col;
-	      flag = 1;
-	  }else{
-	      UpdateMapSignal(row,j2,ntrack,idhit,module,pulse2,pList);
-	      flag = 0;
-	  } // end if
-      } while(flag == 0);
-    } // for isign
+void AliITSsimulationSPD::RemoveDeadPixels(){
+  // Removes dead pixels on each module (ladder)
+  // This should be called before going from sdigits to digits (FrompListToDigits)
+  Int_t mod = GetModuleNumber();
+  AliITSCalibrationSPD* calObj = (AliITSCalibrationSPD*) fDetType->GetCalibrationModel(mod);
+
+  Int_t nrDead = calObj->GetNrBad();
+  for (Int_t i=0; iDeleteHit(calObj->GetBadColAt(i), calObj->GetBadRowAt(i));
+  }
 }
 //______________________________________________________________________
-void AliITSsimulationSPD::CreateDigit(Int_t module,AliITSpList *pList) {
-    // The pixels are fired if the energy deposited inside them is above
-    // the threshold parameter ethr. Fired pixed are interpreted as digits
-    // and stored in the file digitfilename. One also needs to write out
-    // cases when there is only noise (nhits==0).
-
-    static AliITS *aliITS  = (AliITS*)gAlice->GetModule("ITS");
-
-    Int_t digits[3];
-    Int_t tracks[3];
-    Int_t hits[3];
-    Float_t charges[3]; 
-    Int_t j1;
-
-    for (Int_t r=1;r<=GetNPixelsZ();r++) {
-	for (Int_t c=1;c<=GetNPixelsX();c++) {
-	    // check if the deposited energy in a pixel is above the
-	    // threshold 
-	    Float_t signal = (Float_t) fMapA2->GetSignal(r,c);
-	    if ( signal > GetThreshold()) {
-		digits[0] = r-1;  // digits starts from 0
-		digits[1] = c-1;  // digits starts from 0
-		//digits[2] = 1;  
-		digits[2] =  (Int_t) signal;  // the signal is stored in
-                                              //  electrons
-		for(j1=0;j1<3;j1++){
-		    tracks[j1] = pList->GetTrack(r,c,j1);
-		    hits[j1]   = pList->GetHit(r,c,j1);
-		    charges[j1] = 0;
-		} // end for j1
-		Float_t phys = 0;
-		aliITS->AddSimDigit(0,phys,digits,tracks,hits,charges);
-//		cout << " CreateSPDDigit mod=" << fModule << " r,c=" << r;
-//		cout <<","<
-
-
- - . - - 
-    */
-    //End_Html
-    Float_t  thr=0.0,sigm=0.0;
-    Double_t signal,sigma;
-    Int_t iz,ix;
-
-    GetThresholds(thr,sigm);
-    sigma = (Double_t) sigm;
-    for(iz=1;iz<=GetNPixelsZ();iz++){
-	for(ix=1;ix<=GetNPixelsX();ix++){
-	    signal = sigma*gRandom->Gaus(); 
-	    fMapA2->SetHit(iz,ix,signal);
-	    // insert in the label-signal-hit list the pixels fired
-	    // only by noise
-	    pList->AddNoise(iz,ix,module,signal);
-	} // end of loop on pixels
-    } // end of loop on pixels
+void AliITSsimulationSPD::FrompListToDigits(){
+   // add noise and electronics, perform the zero suppression and add the
+   // digit to the list
+   // Inputs:
+   //    none.
+   // Outputs:
+   //    none.
+   // Return:
+   //    none.
+   static AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS");
+   Int_t j,ix,iz;
+   Double_t  electronics;
+   Double_t sig;
+   const Int_t    knmaxtrk=AliITSdigit::GetNTracks();
+   static AliITSdigitSPD dig;
+   AliITSSimuParam *simpar = fDetType->GetSimuParam();
+   if(GetDebug(1)) Info("FrompListToDigits","()");
+   for(iz=0; izRndm() >= eff) continue;
+//        } // end  if 
+// END parametrize the efficiency
+// 
+       electronics = simpar->ApplySPDBaselineAndNoise();
+       UpdateMapNoise(ix,iz,electronics);
+       //
+       // Apply Threshold and write Digits.
+       sig = GetMap()->GetSignalOnly(iz,ix);
+       FillHistograms(ix,iz,sig+electronics);
+       if(GetDebug(3)){
+           cout<threshold("<GetNEntries()) {
+               dig.SetTrack(j,GetMap()->GetTrack(iz,ix,j));
+               dig.SetHit(j,GetMap()->GetHit(iz,ix,j));
+           }else { // Default values
+               dig.SetTrack(j,-3);
+               dig.SetHit(j,-1);
+           } // end if GetMap()
+       } // end for j
+       if(GetDebug(3)){
+           cout<GetpListItem(iz,ix))<AddSimDigit(0,&dig);
+       // simulate fo signal response for this pixel hit:
+       fDetType->ProcessSPDDigitForFastOr(fModule, dig.GetCoord1(), dig.GetCoord2());
+   } //  for ix/iz
 }
 //______________________________________________________________________
-void AliITSsimulationSPD::SetMask() {
-    //  Apply a mask to the SPD module. 1% of the pixel channels are
-    //  masked. When the database will be ready, the masked pixels
-    //  should be read from it.
-    Double_t signal;
-    Int_t iz,ix,im;
-    Float_t totMask;
-    Float_t perc = ((AliITSresponseSPD*)fResponse)->GetFractionDead();
-    // in this way we get the same set of random numbers for all runs.
-    // This is a cluge for now.
-    static TRandom *rnd = new TRandom();
-
-    totMask= perc*GetNPixelsZ()*GetNPixelsX();
-    for(im=1;imRndm()*(GetNPixelsX()-1.) + 1.);
-	} while(ix<=0 || ix>GetNPixelsX());
-	do{
-	    iz=(Int_t)(rnd->Rndm()*(GetNPixelsZ()-1.) + 1.);
-	} while(iz<=0 || iz>GetNPixelsZ());
-	signal = -1.;
-	fMapA2->SetHit(iz,ix,signal);
-    } // end loop on masked pixels
+void AliITSsimulationSPD::CreateHistograms(){
+   // create 1D histograms for tests
+   // Inputs:
+   //    none.
+   // Outputs:
+   //    none.
+   // Return:
+   //     none.
+
+   if(GetDebug(1)) Info("CreateHistograms","create histograms");
+
+   fHis = new TObjArray(GetNPixelsZ());
+   fSPDname="spd_";
+   for(Int_t i=0;iAddAt(new TH1F(fSPDname.Data(),"SPD maps",
+                            GetNPixelsX(),0.,(Double_t)GetNPixelsX()),i);
+   } // end for i
 }
 //______________________________________________________________________
-void AliITSsimulationSPD::CreateHistograms() {
-    // Create Histograms
-    Int_t i;
-
-    fHis=new TObjArray(GetNPixelsZ());
-    for(i=0;iFill(ix,v);
 }
 //______________________________________________________________________
-void AliITSsimulationSPD::ResetHistograms() {
-    // Reset histograms for this detector
-    Int_t i;
+void AliITSsimulationSPD::ResetHistograms(){
+   // Reset histograms for this detector
+   // Inputs:
+   //    none.
+   // Outputs:
+   //    none.
+   // Return:
+   //     none.
+
+   if(!GetHistArray()) return; // Only fill if setup.
+   if(GetDebug(2)) Info("FillHistograms","fill histograms");
+   for ( int i=0;iAt(i))    ((TH1F*)fHis->At(i))->Reset();
+   } // end for i
+}
 
-    for(i=0;iReset();
-    } // end for i
+//______________________________________________________________________
+void AliITSsimulationSPD::SetCoupling(Int_t col, Int_t row) {
+   //  Take into account the coupling between adiacent pixels.
+   //  The parameters probcol and probrow are the probability of the
+   //  signal in one pixel shared in the two adjacent pixels along
+   //  the column and row direction, respectively.
+   //  Note pList is goten via GetMap() and module is not need any more.
+   //  Otherwise it is identical to that coded by Tiziano Virgili (BSN).
+   //Begin_Html
+   /*
+     
+
+
+ + . + + 
+   */
+   //End_Html
+   // Inputs:
+   //    Int_t col            z cell index
+   //    Int_t row            x cell index
+   // Outputs:
+   //    none.
+   // Return:
+   //     none.
+   Int_t j1,j2,flag=0;
+   Double_t pulse1,pulse2;
+   Double_t couplR=0.0,couplC=0.0;
+   Double_t xr=0.;
+
+   GetCouplings(couplC,couplR);
+   if(GetDebug(3)) Info("SetCoupling","(col=%d,row=%d) "
+                        "Calling SetCoupling couplC=%e couplR=%e",
+                        col,row,couplC,couplR);
+   j1 = col;
+   j2 = row;
+   pulse1 = GetMap()->GetSignalOnly(col,row);
+   pulse2 = pulse1;
+   for (Int_t isign=-1;isign<=1;isign+=2){// loop in col direction
+       do{
+           j1 += isign;
+           xr = gRandom->Rndm();
+           if ((j1<0) || (j1>GetNPixelsZ()-1) || (xr>couplC)){
+               j1 = col;
+               flag = 1;
+           }else{
+               UpdateMapNoise(row,j1,pulse1);
+               //  flag = 0;
+               flag = 1; // only first next!!
+           } // end if
+       } while(flag == 0);
+       // loop in row direction
+       do{
+           j2 += isign;
+           xr = gRandom->Rndm();
+           if ((j2<0) || (j2>GetNPixelsX()-1) || (xr>couplR)){
+               j2 = row;
+               flag = 1;
+           }else{
+               UpdateMapNoise(j2,col,pulse2);
+               //  flag = 0;
+               flag = 1; // only first next!!
+           } // end if
+       } while(flag == 0);
+   } // for isign
 }
 //______________________________________________________________________
-void AliITSsimulationSPD::WriteSDigits(AliITSpList *pList){
-    // Fills the Summable digits Tree
-    Int_t i,ni,j,nj;
-    static AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS");
-
-    pList->GetMaxMapIndex(ni,nj);
-    for(i=0;iGetSignalOnly(i,j)>0.0){
-	    aliITS->AddSumDigit(*(pList->GetpListItem(i,j)));
-//	    cout << "pListSPD: " << *(pList->GetpListItem(i,j)) << endl;
-//	    cout << " CreateSPDSDigit mod=" << fModule << " r,c=";
-//	    cout << i  <<","<< j << " sig=" << fpList->GetSignalOnly(i,j);
-//	    cout << " noise=" << fpList->GetNoise(i,j) <
+
+
+ + . + + 
+   */
+   //End_Html
+   // Inputs:
+   //    Int_t col            z cell index
+   //    Int_t row            x cell index
+   //    Int_t module         module number
+   // Outputs:
+   //    none.
+   // Return:
+   //     none.
+   Int_t j1,j2,flag=0;
+   Double_t pulse1,pulse2;
+   Double_t couplR=0.0,couplC=0.0;
+
+   GetCouplings(couplC,couplR);
+
+   //  Debugging ...
+//    cout << "Threshold --> " << GetThreshold() << endl;  // dom
+//    cout << "Couplings --> " << couplC << " " << couplR << endl;  // dom
+
+   if(GetDebug(3)) Info("SetCouplingOld","(col=%d,row=%d) "
+                        "Calling SetCoupling couplC=%e couplR=%e",
+                        col,row,couplC,couplR);
+   for (Int_t isign=-1;isign<=1;isign+=2){// loop in col direction
+   pulse1 = GetMap()->GetSignalOnly(col,row);
+   pulse2 = pulse1;
+   j1 = col;
+   j2 = row;
+       do{
+           j1 += isign;
+           pulse1 *= couplC;
+           if ((j1<0)||(j1>GetNPixelsZ()-1)||(pulse1GetSignalOnly(col,row);
+               j1 = col;
+               flag = 1;
+           }else{
+               UpdateMapNoise(row,j1,pulse1);
+               // flag = 0;
+               flag = 1;  // only first next !!
+           } // end if
+       } while(flag == 0);
+       // loop in row direction
+       do{
+           j2 += isign;
+           pulse2 *= couplR;
+           if((j2<0)||(j2>(GetNPixelsX()-1))||(pulse2GetSignalOnly(col,row);
+               j2 = row;
+               flag = 1;
+           }else{
+               UpdateMapNoise(j2,col,pulse2);
+               // flag = 0;
+               flag = 1; // only first next!!
+           } // end if
+       } while(flag == 0);
+   } // for isign
 }
 //______________________________________________________________________
-void AliITSsimulationSPD::FillMapFrompList(AliITSpList *pList){
-    // Fills fMap2A from the pList of Summable digits
-    Int_t ix,iz;
-
-    for(iz=0;izAddSignal(iz,ix,pList->GetSignal(iz,ix));
-    return;
+void AliITSsimulationSPD::GenerateStrobePhase()
+{
+ // Generate randomly the strobe
+ // phase w.r.t to the LHC clock
+ // Done once per event
+   const Double_t kBunchLenght = 25e-9; // LHC clock
+   fStrobePhase = ((Double_t)gRandom->Integer(fStrobeLenght))*kBunchLenght-
+                  (Double_t)fStrobeLenght*kBunchLenght+
+                  kBunchLenght/2;
 }
+