BesTofDigitizerEcV2 Class Reference

#include <BesTofDigitizerEcV2.hh>

Inheritance diagram for BesTofDigitizerEcV2:

Public Member Functions
	BesTofDigitizerEcV2 ()
	~BesTofDigitizerEcV2 ()
virtual void	Digitize (ScintSingle *, BesTofDigitsCollection *)
void	ReadData ()
void	TofPmtInit ()
void	TofPmtAccum (BesTofHit *)
void	DirectPh (G4int, G4ThreeVector, G4double &, G4int &)
G4double	Scintillation (G4int)
G4double	TransitTime ()
void	AccuSignal (G4double, G4int)
void	TofPmtRspns (G4int)
Public Member Functions inherited from BesTofDigitizerV
	BesTofDigitizerV ()
	~BesTofDigitizerV ()
void	Initialize ()

Additional Inherited Members
Protected Attributes inherited from BesTofDigitizerV
BesTofDigitsCollection *	m_besTofDigitsCollection
BesTofHitsCollection *	m_THC
ITofCaliSvc *	m_tofCaliSvc
ITofSimSvc *	m_tofSimSvc
ITofQElecSvc *	m_tofQElecSvc
G4double	m_ADC [2]
G4double	m_TDC [2]
G4int	m_trackIndex
G4double	m_globalTime
Static Protected Attributes inherited from BesTofDigitizerV
static bool	m_booked = false
static NTuple::Tuple *	m_tupleTof1 = 0
static NTuple::Item< double >	m_partId
static NTuple::Item< double >	m_scinNb
static NTuple::Item< double >	m_edep
static NTuple::Item< double >	m_nHits
static NTuple::Item< double >	m_time1st0
static NTuple::Item< double >	m_time1st1
static NTuple::Item< double >	m_timelast0
static NTuple::Item< double >	m_timelast1
static NTuple::Item< double >	m_totalPhot0
static NTuple::Item< double >	m_totalPhot1
static NTuple::Item< double >	m_NphAllSteps
static NTuple::Item< double >	m_max0
static NTuple::Item< double >	m_max1
static NTuple::Item< double >	m_tdc0
static NTuple::Item< double >	m_adc0
static NTuple::Item< double >	m_tdc1
static NTuple::Item< double >	m_adc1
static NTuple::Tuple *	m_tupleTof2 = 0
static NTuple::Item< double >	m_eTotal
static NTuple::Item< double >	m_nDigi
static NTuple::Item< double >	m_partIdMPV
static NTuple::Item< double >	m_scinNbMPV
static NTuple::Item< double >	m_edepMPV
static NTuple::Item< double >	m_nDigiOut
static NTuple::Tuple *	m_tupleTof3 = 0
static NTuple::Item< double >	m_forb
static NTuple::Item< double >	m_timeFlight
static NTuple::Item< double >	m_ddT
static NTuple::Item< double >	m_scinSwim
static NTuple::Item< double >	m_scinTime
static NTuple::Item< double >	m_transitTime
static NTuple::Item< double >	m_endTime
static NTuple::Item< double >	m_edepHit

Detailed Description

Definition at line 33 of file BesTofDigitizerEcV2.hh.

Constructor & Destructor Documentation

BesTofDigitizerEcV2()

BesTofDigitizerEcV2::BesTofDigitizerEcV2

(

)

Definition at line 28 of file BesTofDigitizerEcV2.cc.

{
    ReadData();
    m_timeBinSize=0.005;
 
    //retrieve G4Svc
    ISvcLocator* svcLocator = Gaudi::svcLocator();
    IG4Svc* tmpSvc;
    StatusCode sc = svcLocator->service("G4Svc", tmpSvc);
    m_G4Svc = dynamic_cast<G4Svc *>(tmpSvc);
 
}

~BesTofDigitizerEcV2()

BesTofDigitizerEcV2::~BesTofDigitizerEcV2

(

)

Definition at line 70 of file BesTofDigitizerEcV2.cc.

{
    ;
}

Member Function Documentation

AccuSignal()

void BesTofDigitizerEcV2::AccuSignal	(	G4double	endTime,
		G4int	forb )

Definition at line 302 of file BesTofDigitizerEcV2.cc.

{
    G4int ihst;
    ihst=G4int(endTime/m_timeBinSize);
    if (ihst>0 &&ihst<m_profBinNEcV2)
    {
        m_nPhot[ihst][forb]=m_nPhot[ihst][forb]+1;
        m_totalPhot[forb]=m_totalPhot[forb]+1;
    }
}

Referenced by TofPmtAccum().

Digitize()

void BesTofDigitizerEcV2::Digitize ( ScintSingle * scint, BesTofDigitsCollection * DC )

virtual

Reimplemented from BesTofDigitizerV.

Definition at line 75 of file BesTofDigitizerEcV2.cc.

{
    m_beamTime = m_G4Svc->GetBeamTime() * ns;
    m_besTofDigitsCollection = DC;
 
    G4DigiManager* digiManager = G4DigiManager::GetDMpointer();
 
    G4int THCID = digiManager->GetHitsCollectionID("BesTofHitsCollection");
    m_THC = (BesTofHitsCollection*) (digiManager->GetHitsCollection(THCID));
 
    if (m_THC)
    {
        //for each digi, compute TDC and ADC
        G4int partId, scinNb, nHits;
        G4double edep;
        BesTofHit* hit;
        partId=scint->GetPartId();
        scinNb=scint->GetScinNb();
        edep = scint->GetEdep();
        nHits=scint->GetHitIndexes()->size();
 
        TofPmtInit();
 
        if (edep>0.01)
        {
            for (G4int j=0;j<nHits;j++)
            {
                hit= (*m_THC)[( *(scint->GetHitIndexes()) )[j]];
                TofPmtAccum(hit);
            }
 
            //get final tdc and adc
            TofPmtRspns(partId);
 
            G4double temp0 = m_ADC[0]+m_TDC[0];
            G4double temp1 = m_ADC[1]+m_TDC[1];
            if ( (partId==1&&temp0>0&&temp1>0) || ((partId!=1)&&temp0>0) )
            {
                if (m_ADC[0]>1000) m_ADC[0]=1000;
                if (m_ADC[1]>1000) m_ADC[1]=1000;
                BesTofDigi* digi = new BesTofDigi;
                digi->SetTrackIndex(m_trackIndex);
                digi->SetPartId(partId);
                digi->SetScinNb(scinNb);
                digi->SetForwADC( m_ADC[0]) ;
                digi->SetBackADC( m_ADC[1]) ;
                if (m_TDC[0]!=-999)
                    m_TDC[0] = m_TDC[0]+m_beamTime;
                if (m_TDC[1]!=-999)
                    m_TDC[1] = m_TDC[1]+m_beamTime;
                digi->SetForwTDC( m_TDC[0]) ;
                digi->SetBackTDC( m_TDC[1]) ;
                m_besTofDigitsCollection->insert(digi);
            }
        }
    }
}

Referenced by BesTofDigitizer::Digitize().

DirectPh()

void BesTofDigitizerEcV2::DirectPh	(	G4int	partId,
		G4ThreeVector	emtPos,
		G4double &	pathL,
		G4int &	forb )

Definition at line 224 of file BesTofDigitizerEcV2.cc.

{
    //generation photon have random direction
    //optical parameters of scintillation simulation
    //G4double cos_spanEc = 1-1/m_refIndex;
    G4double cos_spanEc = 1;
    G4double ran;
    //dcos=cos_span*(ran*2.0-1.0);
    ran=G4UniformRand();
    //assuming uniform phi distribution, simulate cos distr only
    G4double dcosEc = cos_spanEc*(ran*2.0-1.0);
    G4double r2=sqrt(emtPos.x()*emtPos.x()+emtPos.y()*emtPos.y());
    if (dcosEc > 0.0 && dcosEc != 1)
    {
        forb=0;  //forward
        pathL = (r2-m_ecR1)/(1.0-dcosEc);
    }
    else if (dcosEc < 0 && dcosEc != -1)
    {
        forb=0;
        pathL=(2*838-r2-m_ecR1)/(1.0+dcosEc);
    }
}

Referenced by TofPmtAccum().

ReadData()

void BesTofDigitizerEcV2::ReadData

(

)

Definition at line 41 of file BesTofDigitizerEcV2.cc.

{
    BesTofGeoParameter* tofPara = BesTofGeoParameter::GetInstance();
 
    m_ecR1         = tofPara->GetEcR1();
    m_tau1Ec       = tofPara->GetTau1Ec();
    m_tau2Ec       = tofPara->GetTau2Ec();
    m_tau3Ec       = tofPara->GetTau3Ec();
    m_tauRatioEc   = tofPara->GetTauRatioEc();
    m_refIndexEc   = tofPara->GetRefIndexEc();
    m_phNConstEc   = tofPara->GetPhNConstEc();
    m_Cpe2pmtEc    = tofPara->GetCpe2pmtEc();
    m_rAngleEc     = tofPara->GetRAngleEc();
    m_QEEc         = tofPara->GetQEEc();
    m_CEEc         = tofPara->GetCEEc();
    m_peCorFacEc   = tofPara->GetPeCorFacEc();
    m_attenEc      = tofPara->GetAttenEc();
 
    m_ttsMeanEc    = tofPara->GetTTSmeanEc();
    m_ttsSigmaEc   = tofPara->GetTTSsigmaEc();
    m_PMTgainEc    = tofPara->GetPMTgainEc();
    m_CeEc         = tofPara->GetCeEc();
    m_riseTimeEc   = tofPara->GetRiseTimeEc();
    m_LLthreshEc   = tofPara->GetLLthreshEc();
    m_HLthreshEc   = tofPara->GetHLthreshEc();
    m_preGainEc    = tofPara->GetPreGainEc();
    m_noiseSigmaEc = tofPara->GetNoiseSigmaEc();
}

Referenced by BesTofDigitizerEcV2().

Scintillation()

G4double BesTofDigitizerEcV2::Scintillation

(

G4int

partId

)

Definition at line 248 of file BesTofDigitizerEcV2.cc.

{
    G4double tmp_tauRatio,tmp_tau1,tmp_tau2,tmp_tau3;
    tmp_tauRatio = m_tauRatioEc;
    tmp_tau1 = m_tau1Ec;
    tmp_tau2 = m_tau2Ec;
    tmp_tau3 = m_tau3Ec;
 
    G4double UniformR = tmp_tauRatio/(1+tmp_tauRatio);
    G4double EmissionTime;
    if (G4UniformRand()>UniformR) {
        while (1) {
            EmissionTime = -tmp_tau2*log( G4UniformRand() );
            if (G4UniformRand()-exp(EmissionTime/tmp_tau2-EmissionTime/tmp_tau1)>1.E-8)
                break;
        }
    }
    else EmissionTime = -tmp_tau3*log( G4UniformRand() );
    return EmissionTime;
}

Referenced by TofPmtAccum().

TofPmtAccum()

void BesTofDigitizerEcV2::TofPmtAccum

(

BesTofHit *

hit

)

Definition at line 149 of file BesTofDigitizerEcV2.cc.

{
    G4double cvelScint=c_light/m_refIndexEc;
    //Get information of this step
    G4ThreeVector pos=hit->GetPos();
    G4int trackIndex = hit->GetTrackIndex();
    G4int partId =hit->GetPartId();
    G4double edep=hit->GetEdep();
    G4double stepL=hit->GetStepL();
    G4double deltaT=hit->GetDeltaT();
    G4double timeFlight=hit->GetTime()-m_beamTime;
    if (timeFlight < m_globalTime)
    {
        m_globalTime = timeFlight;
        m_trackIndex = trackIndex;
    }
 
    G4ThreeVector pDirection=hit->GetPDirection();
    G4double nx=pDirection.x();
    G4double ny=pDirection.y();
    G4double nz=pDirection.z();
 
    //Cpe2pmt=cathode area/scint area
    //peCorFac=correction factor for eff. of available Npe
 
    //number of photons generated in this step
    G4int NphStep;
    NphStep=G4int(edep*m_phNConstEc*m_Cpe2pmtEc*0.66*m_QEEc*m_CEEc*m_peCorFacEc);
 
    G4double ddS, ddT;
    if (NphStep>0)
    {
        for (G4int i=0;i<NphStep;i++)
        {
            //uniform scintilation in each step
            ddS=stepL*G4UniformRand();
            ddT=deltaT*G4UniformRand();
            G4ThreeVector emtPos;
            emtPos.setX(pos.x() + nx*ddS);
            emtPos.setY(pos.y() + ny*ddS);
            emtPos.setZ(pos.z() + nz*ddS);
 
            //check scinillation light whether to hit the pmt or not
            //forb=0/1 for forward(z>0, east) and backward(z<0, west)
            G4double pathL=0;
            G4int forb;
            DirectPh(partId, emtPos, pathL, forb);
 
            //check if photon can reach PMT or not, after attenuation
            G4double ran=G4UniformRand();
            if (pathL>0 && exp(-pathL/m_attenEc) > ran)
            {
                //propagation time in scintillator
                G4double scinSwim=pathL/cvelScint;
                //scintillation timing
                //G4double scinTime=GenPhoton(partId);
                G4double scinTime=Scintillation(partId);
 
                //PMT transit time
                G4double transitTime=TransitTime();
                //sum up all time components
                G4double endTime= timeFlight + ddT + scinSwim + scinTime + transitTime;
 
                //store timing into binned buffer
                AccuSignal(endTime, forb);
 
                //update 1st and last timings here
                if (m_t1st[forb]>endTime)   m_t1st[forb] = endTime;
                if (m_tLast[forb]<endTime)  m_tLast[forb]= endTime;
            }
        }
    }
}

Referenced by Digitize().

TofPmtInit()

void BesTofDigitizerEcV2::TofPmtInit

(

)

Definition at line 133 of file BesTofDigitizerEcV2.cc.

{
    m_trackIndex=-999;
    m_globalTime = 9999;
    m_t1st[0]=100;
    m_t1st[1]=100;
    m_tLast[0]=0.;
    m_tLast[1]=0;
    m_totalPhot[0]=0;
    m_totalPhot[1]=0;
    for (G4int i=0;i<2;i++)
        for (G4int j=0;j<m_profBinNEcV2;j++)
            m_nPhot[j][i]=0;
 
}

Referenced by Digitize().

TofPmtRspns()

void BesTofDigitizerEcV2::TofPmtRspns

(

G4int

partId

)

Definition at line 313 of file BesTofDigitizerEcV2.cc.

{
    //to generate PMT signal shape for single photoelectron.
    //only one time for a job.
    static G4double snpe[m_snpeBinNEcV2];
    static G4int istore_snpe=-1;
 
    //Model: f(t)=Gain*mv_1pe* t**2 * exp-(t/tau)**2/normal-const
    //normalization const =sqrt(pi)*tau*tau*tau/4.0
    G4double tau = m_riseTimeEc;
    G4double norma_const=sqrt(M_PI)*tau*tau*tau/4.0;    //in unit of ns**3
    G4double echarge=1.6e-7;   //in unit of PC
 
    //time profile of pmt signals for Back and Forward PMT.
    G4double profPmt[m_profBinNEcV2][2];
 
    G4double t;
    G4int n1, n2, ii;
    G4int phtn;
 
    if (istore_snpe<0)
    {
        istore_snpe = 1;
        for (G4int i=0;i<m_snpeBinNEcV2;i++)
        {
            t=(i+1)*m_timeBinSize;
            snpe[i]=m_PMTgainEc*m_CeEc*t*t*exp(- (t/tau) * (t/tau) )/norma_const;
        }
    }
    //for barrel and endcap tof: fb=2 or 1
    G4int fb=1;
 
    for (G4int j=0;j<2;j++)
    {
        m_TDC[j]=-999.0;
        m_ADC[j]=-999.0;
    }
    for (G4int j=0; j<fb; j++)
    {
        if (m_totalPhot[j] > 0)
        {
            n1=G4int(m_t1st[j]/m_timeBinSize);
            n2=G4int(m_tLast[j]/m_timeBinSize);
 
            for (G4int i=0;i<m_profBinNEcV2;i++)
                profPmt[i][j]=0.0;
 
            //generate PMT pulse
            n2 = n2<m_profBinNEcV2 ? n2:m_profBinNEcV2;
            for (G4int i=n1;i<n2;i++)
            {
                phtn=m_nPhot[i][j];
                if (phtn>0)
                {
                    for (G4int ihst=0; ihst<m_snpeBinNEcV2; ihst++)
                    {
                        ii=i+ihst;
                        if (ii<m_profBinNEcV2)
                            profPmt[ii][j] += phtn*snpe[ihst];
                    }
                }
            }
 
            //add preamplifier and noise
            for (int i=0;i<m_profBinNEcV2;i++)
            {
                if (profPmt[i][j]>0)
                    profPmt[i][j] = m_preGainEc*profPmt[i][j]+G4RandGauss::shoot(0,m_noiseSigmaEc);
            }
 
            //get pulse height
            G4double max=0;
            for (int i=0;i<m_profBinNEcV2;i++)
            {
                if (profPmt[i][j]>max)
                    max=profPmt[i][j];
            }
 
            G4double tmp_HLthresh, tmp_LLthresh;
            if (partId==1)  {
                tmp_HLthresh = m_HLthreshEc;
                tmp_LLthresh = m_LLthreshEc;
            }
            else {
                tmp_HLthresh = m_HLthreshEc;
                tmp_LLthresh = m_LLthreshEc;
            }
 
            //get final tdc and adc
            if (max>=tmp_HLthresh)
            {
                for (int i=0;i<m_profBinNEcV2;i++)
                {
                    if ( profPmt[i][j] >= tmp_LLthresh )
                    {
                        m_TDC[j] = i*m_timeBinSize;
                        m_ADC[j] = m_preGainEc*m_totalPhot[j]*echarge*m_PMTgainEc;
                        break;
                    }
                }
            }
        }
    }
}

Referenced by Digitize().

TransitTime()

G4double BesTofDigitizerEcV2::TransitTime

(

)

Definition at line 296 of file BesTofDigitizerEcV2.cc.

{
    //get time transit spread
    return G4RandGauss::shoot(m_ttsMeanEc,m_ttsSigmaEc);
}

Referenced by TofPmtAccum().

---

The documentation for this class was generated from the following files:

• BesTofDigitizerEcV2.hh
• BesTofDigitizerEcV2.cc