Message: Scintillation by protons doesn't work Not Logged In (login)
 Next-in-Thread Next-in-Thread
 Next-in-Forum Next-in-Forum

Question Scintillation by protons doesn't work 

Forum: Processes Involving Optical Photons
Date: 20 Aug, 2018
From: Kyrre Skjerdal <Kyrre Skjerdal>

Dear experts I am trying to simulate a scintillation detector, and I am interesting in the response from protons and neutrons.

I don't get any scintillation in the detector. When I send an electron, I get a lot of optical photons, but I have found out that is from the Cerenkov process. When I send a proton through my detector I get no scintillation. (The protons may deflect or other particles may be emitted, but no optical photons.)
I have done this in my physics list:
fScintProcess = new G4Scintillation();
fScintProcess->SetScintillationYieldFactor(1.);
fScintProcess->SetTrackSecondariesFirst(true);
fScintProcess->SetScintillationExcitationRatio(0.0);
fScintProcess->SetScintillationByParticleType(true);

G4ProcessManager* pManager = G4OpticalPhoton::OpticalPhoton()->GetProcessManager();
pManager->AddDiscreteProcess(fAbsorptionProcess);

pManager->AddDiscreteProcess(fRayleighScattering);
pManager->AddDiscreteProcess(fMieHGScatteringProcess);

pManager->AddDiscreteProcess(fBoundaryProcess);

auto particleIterator=GetParticleIterator();
particleIterator->reset();
while ( (*particleIterator)() ){

G4ParticleDefinition* particle = particleIterator->value();
G4String particleName = particle->GetParticleName();
pManager = particle->GetProcessManager();
if (!pManager) {
std::ostringstream o;
o << "Particle " << particleName << "without a Process Manager"; G4Exception("WLSOpticalPhysics::ConstructProcess()","", FatalException,o.str().c_str()); }
if(fCerenkovProcess->IsApplicable(*particle)){
pManager->AddProcess(fCerenkovProcess);

pManager->SetProcessOrdering(fCerenkovProcess,idxPostStep);
}
if(fScintProcess->IsApplicable(*particle)){
pManager->AddProcess(fScintProcess);
pManager->SetProcessOrderingToLast(fScintProcess,idxAtRest);
pManager->SetProcessOrderingToLast(fScintProcess,idxPostStep);
}

}

And here is my materials properties table, from DetectorConstruction:
fStilbeneMPT = new G4MaterialPropertiesTable();
fStilbeneMPT->AddProperty("RINDEX", scintEnergy, scintRIND, scintNum);
fStilbeneMPT->AddProperty("ABSLENGTH", scintEnergy, scintABSL, scintNum);
G4double pEnergy[] = {1.0*MeV, 100.0*MeV};
G4double pYIELD[] = {14000., 14000.};
assert(sizeof(pYIELD) == sizeof(pEnergy));
const G4int pNum = sizeof(pEnergy)/sizeof(G4double);
fStilbeneMPT->AddProperty("PROTONSCINTILLATIONYIELD", pEnergy, pYIELD, pNum);
fStilbeneMPT->AddProperty("ELECTRONSCINTILLATIONYIELD", pEnergy, pYIELD, pNum);
fStilbene->SetMaterialPropertiesTable(fStilbeneMPT);

I have also tried with a constant scintillation yield.

Any suggestions are welcome!

Inline Depth:
 1 1
 All All
Outline Depth:
 1 1
 2 2
 All All
Add message: (add)

1 None: Re: Scintillation by protons doesn't work   (Daren Sawkey - 21 Aug, 2018)
(_ None: Re: Scintillation by protons doesn't work   (Kyrre Skjerdal - 21 Aug, 2018)
(_ None: Re: Scintillation by protons doesn't work   (Daren Sawkey - 21 Aug, 2018)
(_ Ok: Re: Scintillation by protons doesn't work   (Kyrre Skjerdal - 22 Aug, 2018)
(_ None: Re: Scintillation by protons doesn't work   (Daren Sawkey - 22 Aug, 2018)
(_ Ok: Re: Scintillation by protons doesn't work   (Kyrre Skjerdal - 23 Aug, 2018)
 Add Message Add Message
to: "Scintillation by protons doesn't work"

 Subscribe Subscribe

This site runs SLAC HyperNews version 1.11-slac-98, derived from the original HyperNews