/search.css" rel="stylesheet" type="text/css"/> /search.js">
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

In This Package:

Public Member Functions | Private Attributes
DsPhysConsOptical Class Reference

Construct Optical Processes. More...

#include <DsPhysConsOptical.h>

List of all members.

Public Member Functions

 DsPhysConsOptical (const std::string &type, const std::string &name, const IInterface *parent)
virtual ~DsPhysConsOptical ()
void ConstructParticle ()
void ConstructProcess ()
StatusCode initialize ()

Private Attributes

bool m_useCerenkov
 Property: UseCerenkov, UseScintillation, UseRayleigh, UseAbsorption Turn on/off optical processes.
bool m_useScintillation
bool m_useRayleigh
bool m_useAbsorption
bool m_applyWaterQe
double m_cerenPhotonScaleWeight
 Property: CerenPhotonScaleWeight Number (>= 1.0) used to scale down the mean number of optical photons produced.
int m_cerenMaxPhotonPerStep
 Property: CerenMaxPhotonsPerStep Maximum number of photons per step to limit step size.
double m_scintPhotonScaleWeight
 ScintPhotonScaleWeight: Scale down number of produced scintillation photons by this much.
double m_ScintillationYieldFactor
 ScintillationYieldFactor: scale the number of produced scintillation photons per MeV by this much.
bool m_useFastMu300nsTrick
double m_birksConstant1
 Birks constants C1 and C2.
double m_birksConstant2
double m_gammaSlowerTime
double m_gammaSlowerRatio
double m_neutronSlowerTime
double m_neutronSlowerRatio
double m_alphaSlowerTime
double m_alphaSlowerRatio
bool m_doReemission
 ScintDoReemission: Do reemission in scintilator.
bool m_doScintAndCeren
 ScintDoScintAndCeren: Do both scintillation and Cerenkov in scintilator.

Detailed Description

Construct Optical Processes.

bv@bnl.gov Tue Jan 29 15:25:22 2008

Definition at line 20 of file DsPhysConsOptical.h.


Constructor & Destructor Documentation

DsPhysConsOptical::DsPhysConsOptical ( const std::string &  type,
const std::string &  name,
const IInterface *  parent 
)

Definition at line 32 of file DsPhysConsOptical.cc.

    : GiGaPhysConstructorBase(type,name,parent)
{
    declareProperty("CerenMaxPhotonsPerStep",m_cerenMaxPhotonPerStep = 300,
                    "Limit step to at most this many (unscaled) Cerenkov photons.");

    declareProperty("ScintDoReemission",m_doReemission = true,
                    "Do reemission in scintilator.");
    declareProperty("ScintDoScintAndCeren",m_doScintAndCeren = true,
                    "Do both scintillation and Cerenkov in scintilator.");

    declareProperty("UseCerenkov", m_useCerenkov=true, 
                    "Use the Cerenkov process?");
    declareProperty("ApplyWaterQe", m_applyWaterQe=true,
                    "Apply QE for water cerenkov process when OP is created?"
                    "If it is true the CerenPhotonScaleWeight will be disabled in water,"
                    "but it still works for AD and others");  
                    // wz: Maybe we can set the weight of a OP to >1 in future.

    declareProperty("UseScintillation",m_useScintillation=true,
                    "Use the Scintillation process?");
    declareProperty("UseRayleigh", m_useRayleigh=true,
                    "Use the Rayleigh scattering process?");
    declareProperty("UseAbsorption", m_useAbsorption=true,
                    "Use light absorption process?");
    declareProperty("UseFastMu300nsTrick", m_useFastMu300nsTrick=false,
                    "Use Fast muon simulation?");
    declareProperty("ScintillationYieldFactor",m_ScintillationYieldFactor = 1.0,
                    "Scale the number of scintillation photons per MeV by this much.");
   
    declareProperty("BirksConstant1", m_birksConstant1 = 6.5e-3*g/cm2/MeV, 
                    "Birks constant C1");
    declareProperty("BirksConstant2", m_birksConstant2 = 3.0e-6*(g/cm2/MeV)*(g/cm2/MeV), 
                   "Birks constant C2");

    declareProperty("GammaSlowerTime", m_gammaSlowerTime = 149*ns, 
                    "Gamma Slower time constant");
    declareProperty("GammaSlowerRatio", m_gammaSlowerRatio = 0.338, 
                   "Gamma Slower time ratio");

    declareProperty("NeutronSlowerTime", m_neutronSlowerTime = 220*ns, 
                    "Neutron Slower time constant");
    declareProperty("NeutronSlowerRatio", m_neutronSlowerRatio = 0.34, 
                   "Neutron Slower time ratio");

    declareProperty("AlphaSlowerTime", m_alphaSlowerTime = 220*ns, 
                    "Alpha Slower time constant");
    declareProperty("AlphaSlowerRatio", m_alphaSlowerRatio = 0.35, 
                   "Alpha Slower time ratio");
     
    declareProperty("CerenPhotonScaleWeight",m_cerenPhotonScaleWeight = 3.125,
                    "Scale down number of produced Cerenkov photons by this much."); 
    declareProperty("ScintPhotonScaleWeight",m_scintPhotonScaleWeight = 3.125,           
                    "Scale down number of produced scintillation photons by this much."); 
}
DsPhysConsOptical::~DsPhysConsOptical ( ) [virtual]

Definition at line 102 of file DsPhysConsOptical.cc.

{
}

Member Function Documentation

void DsPhysConsOptical::ConstructParticle ( )

Definition at line 106 of file DsPhysConsOptical.cc.

{
}
void DsPhysConsOptical::ConstructProcess ( )

Definition at line 110 of file DsPhysConsOptical.cc.

{
#ifdef USE_CUSTOM_CERENKOV
    
    info () << "Using customized DsG4Cerenkov." << endreq;
    DsG4Cerenkov* cerenkov = 0;
    if (m_useCerenkov) {
        cerenkov = new DsG4Cerenkov();
        cerenkov->SetMaxNumPhotonsPerStep(m_cerenMaxPhotonPerStep);
        cerenkov->SetApplyPreQE(m_cerenPhotonScaleWeight>1.);
        cerenkov->SetPreQE(1./m_cerenPhotonScaleWeight);
        
        // wangzhe   Give user a handle to control it.   
        cerenkov->SetApplyWaterQe(m_applyWaterQe);
        // wz
        cerenkov->SetTrackSecondariesFirst(true);
    }
#else
    info () << "Using standard G4Cerenkov." << endreq;
    G4Cerenkov* cerenkov = 0;
    if (m_useCerenkov) {
        cerenkov = new G4Cerenkov();
        cerenkov->SetMaxNumPhotonsPerStep(m_cerenMaxPhotonPerStep);
        cerenkov->SetTrackSecondariesFirst(true);
    }
#endif

#ifdef USE_CUSTOM_SCINTILLATION
    DsG4Scintillation* scint = 0;
    info() << "Using customized DsG4Scintillation." << endreq;
    scint = new DsG4Scintillation();
    scint->SetBirksConstant1(m_birksConstant1);
    scint->SetBirksConstant2(m_birksConstant2);
    scint->SetGammaSlowerTimeConstant(m_gammaSlowerTime);
    scint->SetGammaSlowerRatio(m_gammaSlowerRatio);
    scint->SetNeutronSlowerTimeConstant(m_neutronSlowerTime);
    scint->SetNeutronSlowerRatio(m_neutronSlowerRatio);
    scint->SetAlphaSlowerTimeConstant(m_alphaSlowerTime);
    scint->SetAlphaSlowerRatio(m_alphaSlowerRatio);
    scint->SetDoReemission(m_doReemission);
    scint->SetDoBothProcess(m_doScintAndCeren);
    scint->SetApplyPreQE(m_scintPhotonScaleWeight>1.);
    scint->SetPreQE(1./m_scintPhotonScaleWeight);
    scint->SetScintillationYieldFactor(m_ScintillationYieldFactor); //1.);
    scint->SetUseFastMu300nsTrick(m_useFastMu300nsTrick);
    scint->SetTrackSecondariesFirst(true);
    if (!m_useScintillation) {
        scint->SetNoOp();
    }
#else  // standard G4 scint
    G4Scintillation* scint = 0;
    if (m_useScintillation) {
        info() << "Using standard G4Scintillation." << endreq;
        scint = new G4Scintillation();
        scint->SetScintillationYieldFactor(m_ScintillationYieldFactor); // 1.);
        scint->SetTrackSecondariesFirst(true);
    }
#endif

    G4OpAbsorption* absorb = 0;
    if (m_useAbsorption) {
        absorb = new G4OpAbsorption();
    }

    DsG4OpRayleigh* rayleigh = 0;
    if (m_useRayleigh) {
        rayleigh = new DsG4OpRayleigh();
        //        rayleigh->SetVerboseLevel(2);
    }

    //G4OpBoundaryProcess* boundproc = new G4OpBoundaryProcess();
    DsG4OpBoundaryProcess* boundproc = new DsG4OpBoundaryProcess();
    boundproc->SetModel(unified);

    G4FastSimulationManagerProcess* fast_sim_man
        = new G4FastSimulationManagerProcess("fast_sim_man");
    
    theParticleIterator->reset();
    while( (*theParticleIterator)() ) {

        G4ParticleDefinition* particle = theParticleIterator->value();
        G4ProcessManager* pmanager = particle->GetProcessManager();
    
        // Caution: as of G4.9, Cerenkov becomes a Discrete Process.
        // This code assumes a version of G4Cerenkov from before this version.

        if(cerenkov && cerenkov->IsApplicable(*particle)) {
            pmanager->AddProcess(cerenkov);
            pmanager->SetProcessOrdering(cerenkov, idxPostStep);
            debug() << "Process: adding Cherenkov to " 
                    << particle->GetParticleName() << endreq;
        }

        if(scint && scint->IsApplicable(*particle)) {
            pmanager->AddProcess(scint);
            pmanager->SetProcessOrderingToLast(scint, idxAtRest);
            pmanager->SetProcessOrderingToLast(scint, idxPostStep);
            debug() << "Process: adding Scintillation to "
                    << particle->GetParticleName() << endreq;
        }

        if (particle == G4OpticalPhoton::Definition()) {
            if (absorb)
                pmanager->AddDiscreteProcess(absorb);
            if (rayleigh)
                pmanager->AddDiscreteProcess(rayleigh);
            pmanager->AddDiscreteProcess(boundproc);
            //pmanager->AddDiscreteProcess(pee);
            pmanager->AddDiscreteProcess(fast_sim_man);
        }
    }
}
StatusCode DsPhysConsOptical::initialize ( )

Definition at line 90 of file DsPhysConsOptical.cc.

{
    info()<<"Photons prescaling is "<<( m_cerenPhotonScaleWeight>1.?"on":"off" )
          <<" for Cerenkov. Preliminary applied efficiency is "
          <<1./m_cerenPhotonScaleWeight<<" (weight="<<m_cerenPhotonScaleWeight<<")"<<endreq;
    info()<<"Photons prescaling is "<<( m_scintPhotonScaleWeight>1.?"on":"off" )
          <<" for Scintillation. Preliminary applied efficiency is "
          <<1./m_scintPhotonScaleWeight<<" (weight="<<m_scintPhotonScaleWeight<<")"<<endreq;
    info()<<"WaterQE is turned "<<(m_applyWaterQe?"on":"off")<<" for Cerenkov."<<endreq;
    return this->GiGaPhysConstructorBase::initialize();
}

Member Data Documentation

Property: UseCerenkov, UseScintillation, UseRayleigh, UseAbsorption Turn on/off optical processes.

Default, all are on.

Definition at line 37 of file DsPhysConsOptical.h.

Definition at line 37 of file DsPhysConsOptical.h.

Definition at line 37 of file DsPhysConsOptical.h.

Definition at line 37 of file DsPhysConsOptical.h.

Definition at line 41 of file DsPhysConsOptical.h.

Property: CerenPhotonScaleWeight Number (>= 1.0) used to scale down the mean number of optical photons produced.

For each actual secondary optical photon track produced, it will be given a weight equal to this scale for scaling up if detected later. Default is 1.0.

Definition at line 48 of file DsPhysConsOptical.h.

Property: CerenMaxPhotonsPerStep Maximum number of photons per step to limit step size.

This value is independent from PhotonScaleWeight. Default is 300.

Definition at line 53 of file DsPhysConsOptical.h.

ScintPhotonScaleWeight: Scale down number of produced scintillation photons by this much.

Definition at line 57 of file DsPhysConsOptical.h.

ScintillationYieldFactor: scale the number of produced scintillation photons per MeV by this much.

This controls the primary yield of scintillation photons per MeV of deposited energy.

Definition at line 64 of file DsPhysConsOptical.h.

Definition at line 65 of file DsPhysConsOptical.h.

Birks constants C1 and C2.

Definition at line 68 of file DsPhysConsOptical.h.

Definition at line 69 of file DsPhysConsOptical.h.

Definition at line 71 of file DsPhysConsOptical.h.

Definition at line 72 of file DsPhysConsOptical.h.

Definition at line 73 of file DsPhysConsOptical.h.

Definition at line 74 of file DsPhysConsOptical.h.

Definition at line 75 of file DsPhysConsOptical.h.

Definition at line 76 of file DsPhysConsOptical.h.

ScintDoReemission: Do reemission in scintilator.

Definition at line 80 of file DsPhysConsOptical.h.

ScintDoScintAndCeren: Do both scintillation and Cerenkov in scintilator.

Definition at line 83 of file DsPhysConsOptical.h.


The documentation for this class was generated from the following files:
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

Generated on Fri May 16 2014 10:17:58 for DetSim by doxygen 1.7.4