-/*
- questions to experts: why RemoveDeadPixels should be called before FrompListToDigits ?
-
-
-*/
/**************************************************************************
* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
* provided "as is" without express or implied warranty. *
**************************************************************************/
-#include <Riostream.h>
#include <TGeoGlobalMagField.h>
#include <TH1.h>
+#include <TH2.h>
#include <TString.h>
#include "AliITSU.h"
#include "AliITSUDigitPix.h"
#include "AliITSUHit.h"
-#include "AliITSUModule.h"
+#include "AliITSUChip.h"
#include "AliITSUSensMap.h"
#include "AliITSUCalibrationPix.h"
#include "AliITSUSegmentationPix.h"
#include "AliMathBase.h"
#include "AliITSUSimuParam.h"
#include "AliITSUSDigit.h"
+#include "AliITSUParamList.h"
+
+using std::cout;
+using std::endl;
+using namespace TMath;
ClassImp(AliITSUSimulationPix)
////////////////////////////////////////////////////////////////////////
//
// AliITSUSimulationPix is to do the simulation of pixels
//
+// 2013 Feb: Added MonoPix response and nois calculation al la MIMOSA32 (levente.molnar@cern.ch)
+//
////////////////////////////////////////////////////////////////////////
//______________________________________________________________________
AliITSUSimulationPix::AliITSUSimulationPix()
: fTanLorAng(0)
- ,fStrobe(kTRUE)
- ,fStrobeLenght(4)
- ,fStrobePhase(-12.5e-9)
+ ,fGlobalChargeScale(1.0)
+ ,fSpread2DHisto(0)
+ ,fSpreadFun(0)
+ ,fROTimeFun(0)
{
// Default constructor.
+ SetUniqueID(AliITSUGeomTGeo::kChipTypePix);
}
//______________________________________________________________________
AliITSUSimulationPix::AliITSUSimulationPix(AliITSUSimuParam* sim,AliITSUSensMap* map)
:AliITSUSimulation(sim,map)
,fTanLorAng(0)
- ,fStrobe(kTRUE)
- ,fStrobeLenght(4)
- ,fStrobePhase(-12.5e-9)
+ ,fGlobalChargeScale(1.0)
+ ,fSpread2DHisto(0)
+ ,fSpreadFun(0)
+ ,fROTimeFun(0)
{
// standard constructor
+ SetUniqueID(AliITSUGeomTGeo::kChipTypePix);
Init();
}
AliITSUSimulationPix::AliITSUSimulationPix(const AliITSUSimulationPix &s)
:AliITSUSimulation(s)
,fTanLorAng(s.fTanLorAng)
- ,fStrobe(s.fStrobe)
- ,fStrobeLenght(s.fStrobeLenght)
- ,fStrobePhase(s.fStrobePhase)
+ ,fGlobalChargeScale(s.fGlobalChargeScale)
+ ,fSpread2DHisto(s.fSpread2DHisto)
+ ,fSpreadFun(s.fSpreadFun)
+ ,fROTimeFun(s.fROTimeFun)
{
// Copy Constructor
}
// Assignment operator
if (&s == this) return *this;
AliITSUSimulation::operator=(s);
- fStrobe = s.fStrobe;
- fStrobeLenght = s.fStrobeLenght;
- fStrobePhase = s.fStrobePhase;
+ fSpread2DHisto = s.fSpread2DHisto;
+ //
+ fGlobalChargeScale = s.fGlobalChargeScale;
+ fSpreadFun = s.fSpreadFun;
+ fROTimeFun = s.fROTimeFun;
+ //
return *this;
}
AliMagF* fld = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
if (!fld) AliFatal("The field is not initialized");
Double_t bz = fld->SolenoidField();
- fTanLorAng = TMath::Tan(weightHole*fSimuParam->LorentzAngleHole(bz) +
+ fTanLorAng = Tan(weightHole*fSimuParam->LorentzAngleHole(bz) +
weightEle*fSimuParam->LorentzAngleElectron(bz));
return kTRUE;
}
//_____________________________________________________________________
-void AliITSUSimulationPix::SDigitiseModule(AliITSUModule *mod, Int_t /*mask*/,Int_t event, AliITSsegmentation* seg)
+void AliITSUSimulationPix::SDigitiseChip()
{
// This function begins the work of creating S-Digits.
- if (!(mod->GetNHits())) {
- AliDebug(1,Form("In event %d module %d there are %d hits returning.",
- event, mod->GetIndex(),mod->GetNHits()));
- return;
- }
- //
- if (fStrobe) if (event != fEvent) GenerateStrobePhase();
- InitSimulationModule(mod->GetIndex(), event, seg );
- // Hits2SDigits(mod);
- Hits2SDigitsFast(mod);
- if (fSimuParam->GetPixAddNoisyFlag()) AddNoisyPixels();
- if (fSimuParam->GetPixRemoveDeadFlag()) RemoveDeadPixels();
+
+ AliDebug(10,Form("In event %d chip %d there are %d hits", fEvent, fChip->GetIndex(),fChip->GetNHits()));
+ if (fChip->GetNHits()) Hits2SDigitsFast();
+ if (!fSensMap->GetEntries()) return;
WriteSDigits();
ClearMap();
+ //
}
//______________________________________________________________________
}
//______________________________________________________________________
-void AliITSUSimulationPix::FinishSDigitiseModule()
+void AliITSUSimulationPix::FinishSDigitiseChip()
{
// This function calls SDigitsToDigits which creates Digits from SDigits
FrompListToDigits();
}
//______________________________________________________________________
-void AliITSUSimulationPix::DigitiseModule(AliITSUModule *mod,Int_t /*mask*/, Int_t event, AliITSsegmentation* seg)
+void AliITSUSimulationPix::DigitiseChip()
{
// 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 Hits2SDigits
//
- if (fStrobe) if (event != fEvent) GenerateStrobePhase();
- InitSimulationModule( mod->GetIndex(), event, seg );
- // Hits2SDigits(mod);
- Hits2SDigitsFast(mod);
- //
- if (fSimuParam->GetPixAddNoisyFlag()) AddNoisyPixels();
- if (fSimuParam->GetPixRemoveDeadFlag()) RemoveDeadPixels();
- FrompListToDigits();
- ClearMap();
+ // pick charge spread function
+ Hits2SDigitsFast();
+ FinishSDigitiseChip();
}
//______________________________________________________________________
-void AliITSUSimulationPix::Hits2SDigits(AliITSUModule *mod)
+void AliITSUSimulationPix::Hits2SDigits()
{
// Does the charge distributions using Gaussian diffusion charge charing.
- const Double_t kmictocm = 1.0e-4; // convert microns to cm.
- const Double_t kcmtomic = 1.e4;
- const Double_t kBunchLenght = 25e-9; // LHC clock
- Int_t nhits = mod->GetNHits();
+ Int_t nhits = fChip->GetNHits();
if (!nhits) return;
//
Int_t h,ix,iz,i;
Int_t idtrack;
Float_t x,y,z; // keep coordinates float (required by AliSegmentation)
- 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 t,tp,st,dt=0.2,el,sig,sigx,sigz,fda;
- Double_t thick = 0.5*kmictocm*fSeg->Dy(); // Half Thickness
- fSimuParam->GetPixSigmaDiffusionAsymmetry(fda);
+ Double_t tof,x0=0.0,x1=0.0,y0=0.0,y1=0.0,z0=0.0,z1=0.0,de=0.0;
+ Double_t t,tp,st,dt=0.2,el;
+ Double_t thick = 0.5*fSeg->Dy(); // Half Thickness
+
//
for (h=0;h<nhits;h++) {
- if (fStrobe &&
- ((mod->GetHit(h)->GetTOF()<fStrobePhase) ||
- (mod->GetHit(h)->GetTOF()>(fStrobePhase+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 (!fChip->LineSegmentL(h,x0,x1,y0,y1,z0,z1,de,tof,idtrack)) continue;
+ st = Sqrt(x1*x1+y1*y1+z1*z1);
if (st>0.0) {
- st = (Double_t)((Int_t)(st*kcmtomic)); // number of microns
+ st = (Double_t)((Int_t)(st*1e4)); // number of microns
if (st<=1.0) st = 1.0;
- dt = 1.0/st;
+ dt = 1.0/st; // RS TODO: do we need 1 micron steps?
+ double dy = dt*thick;
+ y = -0.5*dy;
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;
+ y += dy;
z = z0+z1*tp;
if (!(fSeg->LocalToDet(x,z,ix,iz))) continue; // outside
- el = dt * de / fSimuParam->GetGeVToCharge();
- //
- if (GetDebug(1)) if (el<=0.0) cout<<"el="<<el<<" dt="<<dt<<" de="<<de<<endl; // end if GetDebug
+ el = fGlobalChargeScale * dt * de / fSimuParam->GetGeVToCharge();
//
- sig = fSimuParam->SigmaDiffusion1D(TMath::Abs(thick + y));
- sigx=sig;
- sigz=sig*fda;
- if (fSimuParam->GetPixLorentzDrift()) ld=(y+thick)*fTanLorAng;
- SpreadChargeAsym(x,z,ix,iz,el,sigx,sigz,ld,idtrack,h);
+ if (fSimuParam->GetPixLorentzDrift()) x += y*fTanLorAng;
+ SpreadCharge2D(x,z,y,ix,iz,el,tof,idtrack,h);
} // end for t
} else { // st == 0.0 deposit it at this point
x = x0;
- y = y0;
+ y = y0 + 0.5*thick;
z = z0;
if (!(fSeg->LocalToDet(x,z,ix,iz))) continue; // outside
- el = de / fSimuParam->GetGeVToCharge();
- sig = fSimuParam->SigmaDiffusion1D(TMath::Abs(thick + y));
- sigx = sig;
- sigz = sig*fda;
- if (fSimuParam->GetPixLorentzDrift()) ld=(y+thick)*fTanLorAng;
- SpreadChargeAsym(x,z,ix,iz,el,sigx,sigz,ld,idtrack,h);
+ el = fGlobalChargeScale * de / fSimuParam->GetGeVToCharge();
+ if (fSimuParam->GetPixLorentzDrift()) x += y*fTanLorAng;
+ SpreadCharge2D(x,z,y,ix,iz,el,tof,idtrack,h);
} // end if st>0.0
} // Loop over all hits h
//
int nd = fSensMap->GetEntriesUnsorted(); // use unsorted access when possible, since it is faster
AliITSUSDigit* dg = 0;
switch (fSimuParam->GetPixCouplingOption()) {
+ case AliITSUSimuParam::kNoCouplingPix :
+ break;
case AliITSUSimuParam::kNewCouplingPix :
for (i=nd;i--;) {
dg = (AliITSUSDigit*)fSensMap->AtUnsorted(i);
if (fSensMap->IsDisabled(dg)) continue;
- SetCoupling(dg,idtrack,h);
+ SetCoupling(dg);
}
break;
case AliITSUSimuParam::kOldCouplingPix:
for (i=nd;i--;) {
dg = (AliITSUSDigit*)fSensMap->AtUnsorted(i);
if (fSensMap->IsDisabled(dg)) continue;
- SetCouplingOld(dg,idtrack,h);
+ SetCouplingOld(dg);
}
break;
default:
}
//______________________________________________________________________
-void AliITSUSimulationPix::Hits2SDigitsFast(AliITSUModule *mod)
+void AliITSUSimulationPix::Hits2SDigitsFast()
{
// Does the charge distributions using Gaussian diffusion charge charing. // Inputs:
- // AliITSUModule *mod Pointer to this module
- // Output:
- // none.
- // Return:
- // none.
- const Double_t kmictocm = 1.0e-4; // convert microns to cm.
- const Int_t kn10=10;
- const Double_t kti[kn10]={7.443716945e-3,2.166976971e-1,3.397047841e-1,
- 4.325316833e-1,4.869532643e-1,5.130467358e-1,
- 5.674683167e-1,6.602952159e-1,7.833023029e-1,
- 9.255628306e-1};
- const Double_t kwi[kn10]={1.477621124e-1,1.346333597e-1,1.095431813e-1,
- 7.472567455e-2,3.333567215e-2,3.333567215e-2,
- 7.472567455e-2,1.095431813e-1,1.346333597e-1,
- 1.477621124e-1};
- const Double_t kBunchLenght = 25e-9; // LHC clock
- TObjArray *hits = mod->GetHits();
+ // AliITSUChip *mod Pointer to this chip
+ //
+ TObjArray *hits = fChip->GetHits();
Int_t nhits = hits->GetEntriesFast();
if (nhits<=0) return;
//
Int_t h,ix,iz,i;
Int_t idtrack;
Float_t x,y,z; // keep coordinates float (required by AliSegmentation)
- Double_t x0=0.0,x1=0.0,y0=0.0,y1=0.0,z0=0.0,z1=0.0;
- Double_t t,st,el,sig,sigx,sigz,fda,de=0.0,ld=0.0;
- Double_t thick = 0.5*kmictocm*fSeg->Dy(); // Half thickness
- fSimuParam->GetPixSigmaDiffusionAsymmetry(fda);
+ Double_t tof,x0=0.0,x1=0.0,y0=0.0,y1=0.0,z0=0.0,z1=0.0;
+ Double_t step,st,el,de=0.0;
+ Double_t minDim = Min(fSeg->Dpx(1),fSeg->Dpz(1)); // RStmp: smallest pitch
+ Double_t thick = fSeg->Dy();
//
for (h=0;h<nhits;h++) {
- // Check if the hit is inside readout window
- if (fStrobe &&
- ((mod->GetHit(h)->GetTOF()<fStrobePhase) ||
- (mod->GetHit(h)->GetTOF()>(fStrobePhase+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;i<kn10;i++) { // Integrate over t
- t = kti[i];
- x = x0+x1*t;
- y = y0+y1*t;
- z = z0+z1*t;
+ if (!fChip->LineSegmentL(h,x0,x1,y0,y1,z0,z1,de,tof,idtrack)) continue;
+ //
+ st = Sqrt(x1*x1+y1*y1+z1*z1);
+ if (st>0.0) {
+ int np = int(1.5*st/minDim); //RStmp: inject the points in such a way that there is ~1.5 point per cell
+ np = TMath::Max(1.0*np,fResponseParam->GetParameter(kSpreadFunMinSteps));
+ AliDebug(10,Form(" Number of charge injection steps is set to %d ",np));
+ double dstep = 1./np;
+ double dy = dstep*thick;
+ y = -0.5*dy;
+ step = -0.5*dstep;
+ for (i=0;i<np;i++) { //RStmp Integrate over t
+ // for (i=0;i<kn10;i++) { // Integrate over t
+ step += dstep; // RStmp kti[i];
+ x = x0+x1*step;
+ y += dy;
+ z = z0+z1*step;
if (!(fSeg->LocalToDet(x,z,ix,iz))) continue; // outside
- el = kwi[i]*de/fSimuParam->GetGeVToCharge();
- sig = fSimuParam->SigmaDiffusion1D(TMath::Abs(thick + y));
- sigx=sig;
- sigz=sig*fda;
- if (fSimuParam->GetPixLorentzDrift()) ld=(y+thick)*fTanLorAng;
- SpreadChargeAsym(x,z,ix,iz,el,sigx,sigz,ld,idtrack,h);
- // cout << "sigx sigz " << sigx << " " << sigz << endl; // dom
+ el = fGlobalChargeScale * dstep * de/fSimuParam->GetGeVToCharge();
+ if (fSimuParam->GetPixLorentzDrift()) x += y*fTanLorAng;
+ SpreadCharge2D(x,z,y,ix,iz,el,tof,idtrack,h);
} // end for i // End Integrate over t
+ }
else { // st == 0.0 deposit it at this point
x = x0;
- y = y0;
+ y = y0+0.5*thick;
z = z0;
if (!(fSeg->LocalToDet(x,z,ix,iz))) continue; // outside
- el = de / fSimuParam->GetGeVToCharge();
- sig = fSimuParam->SigmaDiffusion1D(TMath::Abs(thick + y));
- sigx=sig;
- sigz=sig*fda;
- if (fSimuParam->GetPixLorentzDrift()) ld=(y+thick)*fTanLorAng;
- SpreadChargeAsym(x,z,ix,iz,el,sigx,sigz,ld,idtrack,h);
+ el = fGlobalChargeScale * de / fSimuParam->GetGeVToCharge();
+ if (fSimuParam->GetPixLorentzDrift()) x += y*fTanLorAng;
+ SpreadCharge2D(x,z,y,ix,iz,el,tof,idtrack,h);
} // end if st>0.0
} // Loop over all hits h
int nd = fSensMap->GetEntriesUnsorted(); // use unsorted access when possible, since it is faster
AliITSUSDigit* dg = 0;
switch (fSimuParam->GetPixCouplingOption()) {
+ case AliITSUSimuParam::kNoCouplingPix :
+ break;
case AliITSUSimuParam::kNewCouplingPix :
for (i=nd;i--;) {
dg = (AliITSUSDigit*)fSensMap->AtUnsorted(i);
if (fSensMap->IsDisabled(dg)) continue;
- SetCoupling(dg,idtrack,h);
+ SetCoupling(dg);
}
case AliITSUSimuParam::kOldCouplingPix:
for (i=nd;i--;) {
dg = (AliITSUSDigit*)fSensMap->AtUnsorted(i);
if (fSensMap->IsDisabled(dg)) continue;
- SetCouplingOld(dg,idtrack,h);
+ SetCouplingOld(dg);
}
break;
default:
}
//______________________________________________________________________
-void AliITSUSimulationPix::SpreadCharge(Double_t x0,Double_t z0,
- Int_t ix0,Int_t iz0,
- Double_t el,Double_t sig,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*sig*sig)*exp(-arg)
- // arg=((x-x0)*(x-x0)/2*sig*sig)+((z-z0*z-z0)/2*sig*sig)
- // if fSimuParam->GetPixLorentzDrift()=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 sig Sigma difusion for this step (y0 dependent)
- // Double_t ld lorentz drift in x for this step (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; // keep coordinates float (required by AliSegmentation)
- Double_t s,sp,x1,x2,z1,z2;
- //
- if (GetDebug(2)) AliInfo(Form("(x0=%e,z0=%e,ix0=%d,iz0=%d,el=%e,sig=%e,t=%d,i=%d)",x0,z0,ix0,iz0,el,sig,t,hi));
- if (sig<=0.0) { // if sig<=0 No diffusion to simulate.
- UpdateMapSignal(iz0,ix0,t,hi,el);
- return;
- } // end if
- sp = 1.0/(sig*kRoot2);
- ixs = TMath::Max(-knx+ix0,0);
- ixe = TMath::Min(knx+ix0,fSeg->Npx()-1);
- izs = TMath::Max(-knz+iz0,0);
- ize = TMath::Min(knz+iz0,fSeg->Npz()-1);
- for (ix=ixs;ix<=ixe;ix++)
- for (iz=izs;iz<=ize;iz++) {
- fSeg->DetToLocal(ix,iz,x,z); // pixel center
- double dxi = 0.5*kmictocm*fSeg->Dpx(ix);
- double dzi = 0.5*kmictocm*fSeg->Dpz(iz);
- x1 = x;
- z1 = z;
- x2 = x1 + dxi; // Upper
- x1 -= dxi; // Lower
- z2 = z1 + dzi; // Upper
- z1 -= dzi; // 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 = el*0.25; // Correction based on definision of Erfc
- s *= AliMathBase::ErfcFast(sp*x1) - AliMathBase::ErfcFast(sp*x2);
- if (GetDebug(3)) {
- cout <<"el="<<el<<" ix0="<<ix0<<" ix="<<ix<<" x0="<<x<<
- " iz0="<<iz0<<" iz="<<iz<<" z0="<<z<<
- " sp*x1="<<sp*x1<<" sp*x2="<<sp*x2<<" s="<<s;
- } // end if GetDebug
- s *= AliMathBase::ErfcFast(sp*z1) - AliMathBase::ErfcFast(sp*z2);
- if (GetDebug(3)) cout<<" sp*z1="<<sp*z1<<" sp*z2="<<sp*z2<<" s="<<s<< endl; // end if GetDebug
- if (s>fSimuParam->GetPixMinElToAdd()) UpdateMapSignal(iz,ix,t,hi,s);
- } // end for ix, iz
- //
-}
-
-//______________________________________________________________________
-void AliITSUSimulationPix::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)
+void AliITSUSimulationPix::SpreadCharge2D(Double_t x0,Double_t z0, Double_t dy, Int_t ix0,Int_t iz0,
+ Double_t el, Double_t tof, Int_t tID, Int_t hID)
{
// 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 fSimuParam->GetPixLorentzDrift()=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
+ // Double_t x0 x position of point where charge is liberated (local)
+ // Double_t z0 z position of point where charge is liberated (local)
+ // Double_t dy distance from the entrance surface (diffusion sigma may depend on it)
// 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 (z0 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 = 3; // RS: TO TUNE
- const Double_t kRoot2 = 1.414213562; // Sqrt(2).
- const Double_t kmictocm = 1.0e-4; // convert microns to cm.
+ // Int_t tID track number
+ // Int_t hID hit "hit" index number
+ //
Int_t ix,iz,ixs,ixe,izs,ize;
Float_t x,z; // keep coordinates float (required by AliSegmentation)
- Double_t s,spx,spz,x1,x2,z1,z2;
- //
- if (GetDebug(2)) AliInfo(Form("SpreadChargeAsym","(x0=%e,z0=%e,ix0=%d,iz0=%d,el=%e,sigx=%e, sigz=%e, t=%d,i=%d,ld=%e)",x0,z0,ix0,iz0,el,sigx,sigz,t,hi,ld));
- if (sigx<=0.0 || sigz<=0.0) { // if sig<=0 No diffusion to simulate.
- UpdateMapSignal(iz0,ix0,t,hi,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,fSeg->Npx()-1);
- izs = TMath::Max(-knz+iz0,0);
- ize = TMath::Min(knz+iz0,fSeg->Npz()-1);
+ Float_t xdioshift = 0 , zdioshift = 0 ;
+ Double_t s,dtIn[kNDtSpread]; // data transfered to spread function for integral calculation
+ //
+ if (GetDebug(2)) AliInfo(Form("(x0=%e,z0=%e,dy=%e, ix0=%d,iz0=%d,el=%e,tID=%d,hID=%d)",x0,z0,dy,ix0,iz0,el,tID,hID));
+ //
+ Double_t &x1 = dtIn[kCellX1];
+ Double_t &x2 = dtIn[kCellX2];
+ Double_t &z1 = dtIn[kCellZ1];
+ Double_t &z2 = dtIn[kCellZ2];
+ //
+ int nx = GetResponseParam()->GetParameter(kSpreadFunParamNXoffs);
+ int nz = GetResponseParam()->GetParameter(kSpreadFunParamNZoffs);
+ //
+ dtIn[kCellYDepth] = dy;
+ ixs = Max(-nx+ix0,0);
+ ixe = Min( nx+ix0,fSeg->Npx()-1);
+ izs = Max(-nz+iz0,0);
+ ize = Min( nz+iz0,fSeg->Npz()-1);
for (ix=ixs;ix<=ixe;ix++)
for (iz=izs;iz<=ize;iz++) {
+ //
+ // Check if the hit is inside readout window
+ int cycleRO = (((AliITSUSimulationPix*)this)->*AliITSUSimulationPix::fROTimeFun)(ix,iz,tof);
+ if (Abs(cycleRO)>kMaxROCycleAccept) continue;
+ //
fSeg->DetToLocal(ix,iz,x,z); // pixel center
- double dxi = 0.5*kmictocm*fSeg->Dpx(ix);
- double dzi = 0.5*kmictocm*fSeg->Dpz(iz);
- x1 = x;
- z1 = z;
+ xdioshift = zdioshift = 0;
+ double dxi = fSeg->Dpx(ix);
+ double dzi = fSeg->Dpz(iz);
+ CalcDiodeShiftInPixel(ix,iz,xdioshift,zdioshift); // Check and apply diode shift if needed
+ xdioshift *= dxi;
+ zdioshift *= dzi;
+ dxi *= 0.5;
+ dzi *= 0.5;
+ // printf("DShift: %d %d -> %.4f %.4f\n",ix,iz,xdioshift,zdioshift);
+ x1 = (x + xdioshift) - x0; // calculate distance of cell boundaries from injection center
+ z1 = (z + zdioshift) - z0;
x2 = x1 + dxi; // Upper
x1 -= dxi; // Lower
z2 = z1 + dzi; // Upper
z1 -= dzi; // 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 = el*0.25; // Correction based on definision of Erfc
- s *= AliMathBase::ErfcFast(spx*x1) - AliMathBase::ErfcFast(spx*x2);
- if (GetDebug(3)) cout <<"el="<<el<<" ix0="<<ix0<<" ix="<<ix<<" x0="<<x<<" iz0="<<iz0<<" iz="<<iz<<" z0="<<z<<" spx*x1="<<spx*x1<<" spx*x2="<<spx*x2<<" s="<<s; // end if GetDebug
- s *= AliMathBase::ErfcFast(spz*z1) - AliMathBase::ErfcFast(spz*z2);
- if (GetDebug(3)) cout<<" spz*z1="<<spz*z1<<" spz*z2="<<spz*z2<<" s="<<s<< endl; // end if GetDebug
- if (s>fSimuParam->GetPixMinElToAdd()) UpdateMapSignal(iz,ix,t,hi,s);
+ s = el* (((AliITSUSimulationPix*)this)->*AliITSUSimulationPix::fSpreadFun)(dtIn);
+ if (s>fSimuParam->GetPixMinElToAdd()) UpdateMapSignal(iz,ix,tID,hID,s,cycleRO);
} // end for ix, iz
//
}
+//______________________________________________________________________
+Double_t AliITSUSimulationPix::SpreadFunDoubleGauss2D(const Double_t *dtIn)
+{
+ // calculate integral of charge in the cell with boundaries at X=dtIn[kCellX1]:dtIn[kCellX2]
+ // and Z=dtIn[kCellZ1]:dtIn[kCellZ2]
+ // The spread function is assumed to be double gaussian in 2D
+ // Parameters should be: mean0,sigma0, mean1,sigma1, relative strenght of 2nd gaussian wrt 1st one
+ //
+ // 1st gaussian
+ double intg1 = GausInt2D(fResponseParam->GetParameter(kG2SigX0), // sigX
+ dtIn[kCellX1]-fResponseParam->GetParameter(kG2MeanX0), // x1-xmean
+ dtIn[kCellX2]-fResponseParam->GetParameter(kG2MeanX0), // x2-xmean
+ fResponseParam->GetParameter(kG2SigZ0), // sigZ
+ dtIn[kCellZ1]-fResponseParam->GetParameter(kG2MeanZ0), // z1-zmean
+ dtIn[kCellZ2]-fResponseParam->GetParameter(kG2MeanZ0)); // z2-zmean
+ // 2nd gaussian
+ double intg2 = GausInt2D(fResponseParam->GetParameter(kG2SigX1), // sigX
+ dtIn[kCellX1]-fResponseParam->GetParameter(kG2MeanX1), // x1-xmean
+ dtIn[kCellX2]-fResponseParam->GetParameter(kG2MeanX1), // x2-xmean
+ fResponseParam->GetParameter(kG2SigZ1), // sigZ
+ dtIn[kCellZ1]-fResponseParam->GetParameter(kG2MeanZ1), // z1-zmean
+ dtIn[kCellZ2]-fResponseParam->GetParameter(kG2MeanZ1)); // z2-zmean
+ double scl = fResponseParam->GetParameter(kG2ScaleG2);
+ return (intg1+intg2*scl)/(1+scl);
+ //
+}
+
+
+//______________________________________________________________________
+Double_t AliITSUSimulationPix::SpreadFrom2DHisto(const Double_t *dtIn)
+{
+ // calculate integral of charge in the cell with boundaries at X=dtIn[kCellX1]:dtIn[kCellX2]
+ // and Z=dtIn[kCellZ1]:dtIn[kCellZ2]
+ // The spread function integral is taken from fSpread2DHisto extracted from the sensor response parameters
+ // list in the method SetResponseParam. The histo must return the fraction of charge integrates in the
+ // cell whose center is on the distance X=(dtIn[kCellX1]+dtIn[kCellX2])/2 and Z=(dtIn[kCellZ1]+dtIn[kCellZ2])/2
+ // from the injection point.
+ //
+ Double_t qpixfrac = 0;
+ Double_t xintp = 1e4*(dtIn[kCellX1]+dtIn[kCellX2])/2.0;
+ Double_t zintp = 1e4*(dtIn[kCellZ1]+dtIn[kCellZ2])/2.0;
+ //
+ qpixfrac = fSpread2DHisto->Interpolate(xintp,zintp); //the PSF map is given in um but the dtIn is in cm so we need to convert it
+ //
+ return qpixfrac;
+}
+
+//______________________________________________________________________
+Double_t AliITSUSimulationPix::SpreadFunGauss2D(const Double_t *dtIn)
+{
+ // calculate integral of charge in the cell with boundaries at X=dtIn[kCellX1]:dtIn[kCellX2]
+ // and Z=dtIn[kCellZ1]:dtIn[kCellZ2]
+ // The spread function is assumed to be gaussian in 2D
+ // Parameters should be: mean0,sigma0
+ return GausInt2D(fResponseParam->GetParameter(kG1SigX), // sigX
+ dtIn[kCellX1]-fResponseParam->GetParameter(kG1MeanX), // x1-xmean
+ dtIn[kCellX2]-fResponseParam->GetParameter(kG1MeanX), // x2-xmean
+ fResponseParam->GetParameter(kG1SigZ), // sigZ
+ dtIn[kCellZ1]-fResponseParam->GetParameter(kG1MeanZ), // z1-zmean
+ dtIn[kCellZ2]-fResponseParam->GetParameter(kG1MeanZ)); // z2-zmean
+ //
+}
+
//______________________________________________________________________
void AliITSUSimulationPix::RemoveDeadPixels()
{
- // Removes dead pixels on each module (ladder)
- // This should be called before going from sdigits to digits (FrompListToDigits)
+ // Removes dead pixels on each chip (ladder)
+ // This should be called before going from sdigits to digits (i.e. from FrompListToDigits)
AliITSUCalibrationPix* calObj = (AliITSUCalibrationPix*) GetCalibDead();
if (!calObj) return;
//
- if (calObj->IsBad()) {ClearMap(); return;} // whole module is masked
+ if (calObj->IsBad()) {ClearMap(); return;} // whole chip is masked
//
+ // prepare the list of r/o cycles seen
+ Char_t cyclesSeen[2*kMaxROCycleAccept+1];
+ int ncycles = 0;
+ for (int i=(2*kMaxROCycleAccept+1);i--;) if (fCyclesID[i]) cyclesSeen[ncycles++]=i-kMaxROCycleAccept;
+
// remove single bad pixels one by one
int nsingle = calObj->GetNrBadSingle();
UInt_t col,row;
+ Int_t cycle;
for (int i=nsingle;i--;) {
calObj->GetBadPixelSingle(i,row,col);
- fSensMap->DeleteItem(col,row);
+ for (int icl=ncycles;icl--;) fSensMap->DeleteItem(col,row,cyclesSeen[icl]);
}
int nsd = fSensMap->GetEntriesUnsorted();
for (int isd=nsd;isd--;) {
AliITSUSDigit* sd = (AliITSUSDigit*)fSensMap->AtUnsorted(isd);
if (fSensMap->IsDisabled(sd)) continue;
- fSensMap->GetMapIndex(sd->GetUniqueID(),col,row);
+ fSensMap->GetMapIndex(sd->GetUniqueID(),col,row,cycle);
int chip = fSeg->GetChipFromChannel(0,col);
// if (calObj->IsChipMarkedBad(chip)) fSensMap->Disable(sd); // this will simple mark the hit as bad
if (calObj->IsChipMarkedBad(chip)) fSensMap->DeleteItem(sd); // this will suppress hit in the sorted list
//______________________________________________________________________
void AliITSUSimulationPix::AddNoisyPixels()
{
- // Adds noisy pixels on each module (ladder)
- // This should be called before going from sdigits to digits (FrompListToDigits)
+ // Adds noisy pixels on each chip (ladder)
+ // This should be called before going from sdigits to digits (i.e. FrompListToDigits)
AliITSUCalibrationPix* calObj = (AliITSUCalibrationPix*) GetCalibNoisy();
- if (!calObj) return;
+ if (!calObj) { AliDebug(10,Form(" No Calib Object for Noise!!! ")); return;}
for (Int_t i=calObj->GetNrBad(); i--;) UpdateMapNoise(calObj->GetBadColAt(i), calObj->GetBadRowAt(i),
- 10*fSimuParam->GetPixThreshold(fModule));
+ 10*fSimuParam->GetPixThreshold(fChip->GetIndex()));
//
}
//______________________________________________________________________
void AliITSUSimulationPix::FrompListToDigits()
{
- // add noise and electronics, perform the zero suppression and add the
- // digit to the list
- static AliITSU *aliITS = (AliITSU*)gAlice->GetModule("ITS");
- UInt_t ix,iz;
- Double_t sig;
- const Int_t knmaxtrk=AliITSdigit::GetNTracks();
- static AliITSUDigitPix dig;
- // RS: in principle:
- // 1) for every pixel w/o hit we have to generate a noise and activate the pixel if the noise exceeds the threshold.
- // 2) for every pixel with hit we should add random noise and check if the total signal exceeds the threshold
- // With many channels this will be too time consuming, hence I do the following
- // 1) Precalculate the probability that the nois alone will exceed the threshold.
- // 2) Chose randomly empty pixels according to this probability and apply the noise above threshold.
- // 3) For pixels having a hits apply the usual noise and compare with threshold
+ // add noise and electronics, perform the zero suppression and add the digits to the list
//
// RS may use for ordered random sample generation dl.acm.org/ft_gateway.cfm?id=356313&type=pdf
//
- int maxInd = fSensMap->GetMaxIndex();
- double minProb = 0.1/maxInd;
+ int nsd = fSensMap->GetEntriesUnsorted(); // sdigits added from the signal
//
- int nsd = fSensMap->GetEntries();
- Int_t prevID=0,curID=0;
- TArrayI ordSampleInd(100),ordSample(100);
+ // add different kinds of noise.
+ Bool_t addNoisy = fSimuParam->GetPixAddNoisyFlag() && (nsd>0 || fSimuParam->GetPixNoiseInAllMod()); // do we generate noise?
+ if (addNoisy) {
+ AddNoisyPixels(); // constantly noisy channels
+ AddRandomNoisePixels(0.0); // random noise: at the moment generate noise only for instance 0
+ nsd = fSensMap->GetEntriesUnsorted();
+ }
//
- double probNoisy,noiseSig,noiseMean,thresh = fSimuParam->GetPixThreshold(fModule);
- fSimuParam->GetPixNoise(fModule, noiseSig, noiseMean);
- probNoisy = AliITSUSimuParam::CalcProbNoiseOverThreshold(noiseMean,noiseSig,thresh); // prob. to have noise above threshold
+ if (nsd && fSimuParam->GetPixRemoveDeadFlag()) {
+ RemoveDeadPixels();
+ // note that here we shall use GetEntries instead of GetEntriesUnsorted since the
+ // later operates on the array where the elements are not removed by flagged
+ nsd = fSensMap->GetEntries();
+ }
+ if (!nsd) return; // nothing to digitize
+ //
+ UInt_t row,col;
+ Int_t iCycle,modId = fChip->GetIndex();
+ Double_t sig;
+ const Int_t knmaxtrk=AliITSdigit::GetNTracks();
+ static AliITSU *aliITS = (AliITSU*)gAlice->GetModule("ITS");
+ static AliITSUDigitPix dig;
//
for (int i=0;i<nsd;i++) {
AliITSUSDigit* sd = (AliITSUSDigit*)fSensMap->At(i); // ordered in index
if (fSensMap->IsDisabled(sd)) continue;
- curID = (int)sd->GetUniqueID();
- //
- if (probNoisy>minProb) { // generate randomly noisy pixels above the threshold, with ID's between previous hit and current
- CreateNoisyDigits(prevID,curID,probNoisy, noiseSig, noiseMean);
- prevID = curID+1;
- }
//
- if ((sig=sd->GetSumSignal())<=fSimuParam->GetPixThreshold(fModule)) continue;
- if (TMath::Abs(sig)>2147483647.0) { //RS?
+ if ((sig=sd->GetSumSignal())<=fSimuParam->GetPixThreshold(modId)) continue;
+ if (Abs(sig)>2147483647.0) { //RS?
//PH 2147483647 is the max. integer
//PH This apparently is a problem which needs investigation
AliWarning(Form("Too big or too small signal value %f",sig));
}
- fSensMap->GetMapIndex(sd->GetUniqueID(),iz,ix);
- dig.SetCoord1(iz);
- dig.SetCoord2(ix);
- dig.SetSignal(1);
+ fSensMap->GetMapIndex(sd->GetUniqueID(),col,row,iCycle);
+ dig.SetCoord1(col);
+ dig.SetCoord2(row);
+ dig.SetROCycle(iCycle);
+ dig.SetSignal((Int_t)sig);
dig.SetSignalPix((Int_t)sig);
int ntr = sd->GetNTracks();
for (int j=0;j<ntr;j++) {
dig.SetTrack(j,-3);
dig.SetHit(j,-1);
}
- aliITS->AddSimDigit(AliITSUGeomTGeo::kDetTypePix, &dig);
+ aliITS->AddSimDigit(AliITSUGeomTGeo::kChipTypePix, &dig);
}
- // if needed, add noisy pixels with id from last real hit to maxID
- if (probNoisy>minProb) CreateNoisyDigits(prevID,maxInd,probNoisy, noiseSig, noiseMean);
//
}
//______________________________________________________________________
-Int_t AliITSUSimulationPix::CreateNoisyDigits(Int_t minID,Int_t maxID,double probNoisy, double noise, double base)
+Int_t AliITSUSimulationPix::AddRandomNoisePixels(Double_t tof)
{
- // create random noisy digits above threshold within id range [minID,maxID[
- // see FrompListToDigits for details
+ // create random noisy sdigits above threshold
//
- static AliITSU *aliITS = (AliITSU*)gAlice->GetModule("ITS");
- UInt_t ix,iz;
- static AliITSUDigitPix dig;
+ int modId = fChip->GetIndex();
+ int npix = fSeg->GetNPads();
+ int ncand = gRandom->Poisson( npix*fSimuParam->GetPixFakeRate() );
+ if (ncand<1) return 0;
+ //
+ double probNoisy,noiseSig,noiseMean,thresh = fSimuParam->GetPixThreshold(modId);
+ fSimuParam->GetPixNoise(modId, noiseSig, noiseMean);
+ probNoisy = AliITSUSimuParam::CalcProbNoiseOverThreshold(noiseMean,noiseSig,thresh); // prob. to have noise above threshold
+ UInt_t row,col;
+ Int_t iCycle;
static TArrayI ordSampleInd(100),ordSample(100); //RS!!! static is not thread-safe!!!
- const Int_t knmaxtrk=AliITSdigit::GetNTracks();
//
- Int_t ncand = 0;
- int npix = maxID-minID;
- if (npix<1 || (ncand=gRandom->Poisson(npix*probNoisy))<1) return ncand; // decide how many noisy pixels will be added
ncand = GenOrderedSample(npix,ncand,ordSample,ordSampleInd);
int* ordV = ordSample.GetArray();
int* ordI = ordSampleInd.GetArray();
for (int j=0;j<ncand;j++) {
- fSensMap->GetMapIndex((UInt_t)ordV[ordI[j]],iz,ix); // create noisy digit
- dig.SetCoord1(iz);
- dig.SetCoord2(ix);
- dig.SetSignal(1);
- dig.SetSignalPix((Int_t)AliITSUSimuParam::GenerateNoiseQFunction(probNoisy,base,noise));
- for (int k=knmaxtrk;k--;) {
- dig.SetTrack(k,-3);
- dig.SetHit(k,-1);
- }
- aliITS->AddSimDigit(AliITSUGeomTGeo::kDetTypePix,&dig);
- if (GetDebug(2)) AliInfo(Form("Add noisy pixel %d(%d/%d) Noise=%d",ordV[ordI[j]],iz,ix,dig.GetSignalPix()));
+ fSensMap->GetMapIndex((UInt_t)ordV[ordI[j]],col,row,iCycle); // create noisy digit
+ iCycle = (((AliITSUSimulationPix*)this)->*AliITSUSimulationPix::fROTimeFun)(row,col,tof);
+ UpdateMapNoise(col,row,AliITSUSimuParam::GenerateNoiseQFunction(probNoisy,noiseMean,noiseSig), iCycle);
}
return ncand;
}
+
//______________________________________________________________________
-void AliITSUSimulationPix::SetCoupling(AliITSUSDigit* old, Int_t ntrack, Int_t idhit)
+void AliITSUSimulationPix::SetCoupling(AliITSUSDigit* old)
{
// 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.
+ // Note pList is goten via GetMap() and chip is not need any more.
// Otherwise it is identical to that coded by Tiziano Virgili (BSN).
- //Begin_Html
- /*
- <img src="picts/ITS/barimodel_3.gif">
- </pre>
- <br clear=left>
- <font size=+2 color=red>
- <a href="mailto:tiziano.virgili@cern.ch"></a>.
- </font>
- <pre>
- */
- //End_Html
- // Inputs:
- // old existing AliITSUSDigit
- // Int_t ntrack track incex number
- // Int_t idhit hit index number
UInt_t col,row;
+ Int_t iCycle;
Double_t pulse1,pulse2;
Double_t couplR=0.0,couplC=0.0;
Double_t xr=0.;
//
- fSensMap->GetMapIndex(old->GetUniqueID(),col,row);
+ fSensMap->GetMapIndex(old->GetUniqueID(),col,row,iCycle);
+ int cycle = iCycle;
fSimuParam->GetPixCouplingParam(couplC,couplR);
- if (GetDebug(2)) AliInfo(Form("(col=%d,row=%d,ntrack=%d,idhit=%d) couplC=%e couplR=%e",
- col,row,ntrack,idhit,couplC,couplR));
+ if (GetDebug(2)) AliInfo(Form("(col=%d,row=%d,couplC=%e couplR=%e",
+ col,row,couplC,couplR));
pulse2 = pulse1 = old->GetSignal();
if (pulse1<fSimuParam->GetPixMinElToAdd()) return; // too small signal
for (Int_t isign=-1;isign<=1;isign+=2) {
// loop in col direction
int j1 = int(col) + isign;
xr = gRandom->Rndm();
- if ( !((j1<0) || (j1>fSeg->Npz()-1) || (xr>couplC)) ) UpdateMapSignal(UInt_t(j1),row,ntrack,idhit,pulse1);
+ if ( !((j1<0) || (j1>fSeg->Npz()-1) || (xr>couplC)) ) UpdateMapSignal(UInt_t(j1),row,old->GetTrack(0),old->GetHit(0),pulse1,cycle);
//
// loop in row direction
int j2 = int(row) + isign;
xr = gRandom->Rndm();
- if ( !((j2<0) || (j2>fSeg->Npx()-1) || (xr>couplR)) ) UpdateMapSignal(col,UInt_t(j2),ntrack,idhit,pulse2);
+ if ( !((j2<0) || (j2>fSeg->Npx()-1) || (xr>couplR)) ) UpdateMapSignal(col,UInt_t(j2),old->GetTrack(0),old->GetHit(0),pulse2,cycle);
}
//
}
//______________________________________________________________________
-void AliITSUSimulationPix::SetCouplingOld(AliITSUSDigit* old, Int_t ntrack,Int_t idhit)
+void AliITSUSimulationPix::SetCouplingOld(AliITSUSDigit* old)
{
// Take into account the coupling between adiacent pixels.
// The parameters probcol and probrow are the fractions of the
// signal in one pixel shared in the two adjacent pixels along
// the column and row direction, respectively.
- //Begin_Html
- /*
- <img src="picts/ITS/barimodel_3.gif">
- </pre>
- <br clear=left>
- <font size=+2 color=red>
- <a href="mailto:Rocco.Caliandro@ba.infn.it"></a>.
- </font>
- <pre>
- */
- //End_Html
// Inputs:
// old existing AliITSUSDigit
// ntrack track incex number
// idhit hit index number
- // module module number
+ // chip chip number
//
UInt_t col,row;
+ Int_t cycle;
+ Int_t modId = fChip->GetIndex();
Double_t pulse1,pulse2;
Double_t couplR=0.0,couplC=0.0;
//
- fSensMap->GetMapIndex(old->GetUniqueID(),col,row);
+ fSensMap->GetMapIndex(old->GetUniqueID(),col,row,cycle);
fSimuParam->GetPixCouplingParam(couplC,couplR);
- if (GetDebug(3)) AliInfo(Form("(col=%d,row=%d,ntrack=%d,idhit=%d) couplC=%e couplR=%e",
- col,row,ntrack,idhit,couplC,couplR));
+ if (GetDebug(3)) AliInfo(Form("(col=%d,row=%d,roCycle=%d) couplC=%e couplR=%e",col,row,cycle,couplC,couplR));
//
if (old->GetSignal()<fSimuParam->GetPixMinElToAdd()) return; // too small signal
for (Int_t isign=-1;isign<=1;isign+=2) {// loop in col direction
//
int j1 = int(col)+isign;
pulse1 *= couplC;
- if ((j1<0)||(j1>fSeg->Npz()-1)||(pulse1<fSimuParam->GetPixThreshold(fModule))) pulse1 = old->GetSignal();
- else UpdateMapSignal(UInt_t(j1),row,ntrack,idhit,pulse1);
+ if ((j1<0)||(j1>fSeg->Npz()-1)||(pulse1<fSimuParam->GetPixThreshold(modId))) pulse1 = old->GetSignal();
+ else UpdateMapSignal(UInt_t(j1),row,old->GetTrack(0),old->GetHit(0),pulse1,cycle);
// loop in row direction
int j2 = int(row) + isign;
pulse2 *= couplR;
- if ((j2<0)||(j2>(fSeg->Npx()-1))||(pulse2<fSimuParam->GetPixThreshold(fModule))) pulse2 = old->GetSignal();
- else UpdateMapSignal(col,UInt_t(j2),ntrack,idhit,pulse2);
+ if ((j2<0)||(j2>(fSeg->Npx()-1))||(pulse2<fSimuParam->GetPixThreshold(modId))) pulse2 = old->GetSignal();
+ else UpdateMapSignal(col,UInt_t(j2),old->GetTrack(0),old->GetHit(0),pulse2,cycle);
} // for isign
}
//______________________________________________________________________
-void AliITSUSimulationPix::GenerateStrobePhase()
+void AliITSUSimulationPix::SetResponseParam(AliITSUParamList* resp)
+{
+ // attach response parameterisation data
+ fResponseParam = resp;
+ //
+ int spreadID = Nint(fResponseParam->GetParameter(AliITSUSimulationPix::kChargeSpreadType));
+ const char* hname = 0;
+ fSpread2DHisto = 0;
+ //
+ switch (spreadID) {
+ //
+ case kSpreadFunHisto:
+ fSpreadFun = &AliITSUSimulationPix::SpreadFrom2DHisto;
+ hname = fResponseParam->GetParName(AliITSUSimulationPix::kChargeSpreadType);
+ if (!(fSpread2DHisto=(TH2*)fResponseParam->GetParamObject(hname)))
+ AliFatal(Form("Did not find 2D histo %s for charge spread parameterization",hname));
+ break;
+ //
+ case kSpreadFunDoubleGauss2D:
+ fSpreadFun = &AliITSUSimulationPix::SpreadFunDoubleGauss2D;
+ break;
+ //
+ case kSpreadFunGauss2D:
+ fSpreadFun = &AliITSUSimulationPix::SpreadFunGauss2D;
+ break;
+ //
+ default: AliFatal(Form("Did not find requested spread function id=%d",spreadID));
+ }
+ //
+ int readoutType = Nint(fResponseParam->GetParameter(kReadOutSchemeType));
+ switch (readoutType) {
+ case kReadOutStrobe:
+ fROTimeFun = &AliITSUSimulationPix::GetReadOutCycle;
+ break;
+ case kReadOutRollingShutter:
+ fROTimeFun = &AliITSUSimulationPix::GetReadOutCycleRollingShutter;
+ break;
+ default: AliFatal(Form("Did not find requested readout time type id=%d",readoutType));
+ }
+ //___ Set the Rolling Shutter read-out window
+ fReadOutCycleLength = fResponseParam->GetParameter(kReadOutCycleLength);
+ //___ Pixel discrimination threshold, and the S/N cut
+ fSimuParam->SetPixThreshold(fResponseParam->GetParameter(kPixNoiseMPV) *fResponseParam->GetParameter(kPixSNDisrcCut) , fResponseParam->GetParameter(kPixSNDisrcCut),-1); //for all chips
+ //___ Minimum number of electrons to add
+ fSimuParam->SetPixMinElToAdd(fResponseParam->GetParameter(kPixMinElToAdd));
+ //___ Set the Pixel Noise MPV and Sigma (the noise distribution is Landau not Gauss due to RTN)
+ fSimuParam->SetPixNoise( fResponseParam->GetParameter(kPixNoiseMPV), fResponseParam->GetParameter(kPixNoiseSigma), -1); //for all chips
+ //___ Pixel fake hit rate
+ fSimuParam->SetPixFakeRate( fResponseParam->GetParameter(kPixFakeRate) );
+ //___ To apply the noise or not
+ if ( fResponseParam->GetParameter(kPixNoiseIsOn) > 0.01) fSimuParam->SetPixAddNoisyFlag(kTRUE);
+ else fSimuParam->SetPixAddNoisyFlag(kFALSE);
+ //
+ if(fResponseParam->GetParameter(kPixNoiseInAllMod) > 0.01 ) fSimuParam->SetPixNoiseInAllMod(kTRUE);
+ else fSimuParam->SetPixNoiseInAllMod(kFALSE);
+ //
+ // Double_t vGeVToQ = fSimuParam->GetGeVToCharge();
+ fGlobalChargeScale = fResponseParam->GetParameter(kSpreadFunGlobalQScale);
+
+ AliDebug(10,Form("=============== Setting the response start ============================"));
+ AliDebug(10,Form("=============== RO type: %d",readoutType));
+ AliDebug(10,Form("=============== RO cycle lenght: %lf",fReadOutCycleLength));
+ AliDebug(10,Form("=============== Noise MPV: %lf",fResponseParam->GetParameter(kPixNoiseMPV)));
+ AliDebug(10,Form("=============== Noise Sigma: %lf",fResponseParam->GetParameter(kPixNoiseSigma)));
+ AliDebug(10,Form("=============== Fake rate: %lf",fResponseParam->GetParameter(kPixFakeRate)));
+ AliDebug(10,Form("=============== Noise On/Off: %d",fSimuParam->GetPixAddNoisyFlag()));
+ AliDebug(10,Form("=============== Noise in all mod on/off: %d",fSimuParam->GetPixNoiseInAllMod()));
+ AliDebug(10,Form("=============== Global Charge scale: %lf",fGlobalChargeScale));
+ AliDebug(10,Form("=============== Setting the response done ============================"));
+
+
+
+}
+
+//______________________________________________________________________
+Int_t AliITSUSimulationPix::GetReadOutCycleRollingShutter(Int_t row, Int_t col, Double_t hitTime)
{
- // 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;
+ //
+ // Get the read-out cycle of the hit in the given column/row of the sensor.
+ // hitTime is the time of the subhit (hit is divided to nstep charge deposit) in seconds
+ // globalPhaseShift gives the start of the RO for the cycle in pixel wrt the LHC clock
+ // GetRollingShutterWindow give the with of the rolling shutter read out window
+ //
+ double tmin = fReadOutCycleOffset + fReadOutCycleLength*(double(row)/fSeg->Npx()-1.);
+ int cycle = Nint( (hitTime-tmin)/fReadOutCycleLength - 0.5 );
+ AliDebug(3,Form("Rolling shutter at row%d/col%d: particle time: %e, tmin: %e : tmax %e -> cycle:%d",row,col,hitTime,tmin,
+ tmin+fReadOutCycleLength,cycle));
+ return cycle;
+ //
}
+//______________________________________________________________________
+Int_t AliITSUSimulationPix::GetReadOutCycle(Int_t row, Int_t col, Double_t hitTime)
+{
+ //
+ // Check whether the hit is in the read out window of the given column/row of the sensor
+ //
+ AliDebug(3,Form("Strobe readout: row%d/col%d: particle time: %e, tmin: %e, tmax %e",
+ row,col,hitTime,fReadOutCycleOffset,fReadOutCycleOffset+fReadOutCycleLength));
+ hitTime -= fReadOutCycleOffset+0.5*fReadOutCycleLength;
+ return (hitTime<0 || hitTime>fReadOutCycleLength) ? kMaxROCycleAccept+1 : 0;
+ //
+}
+
+//_______________________________________________________________________
+void AliITSUSimulationPix::CalcDiodeShiftInPixel(Int_t xrow, Int_t zcol, Float_t &x, Float_t &z)
+{
+ //
+ // Calculates the shift of the diode wrt the geometric center of the pixel.
+ // It is needed for staggerred pixel layout or double diode pixels with assymetric center
+ // The shift can depend on the column or line or both...
+ // The x and z are passed in cm
+ //
+ ((AliITSUSegmentationPix*)fSeg)->GetDiodShift(xrow,zcol,x,z);
+ //
+}
+//_______________________________________________________________________