Hi,
Thank you so much for all of the effort you have put in to help me with this.
In my material definition script I have:
#include "Material.hh"
#include "G4SystemOfUnits.hh"
#include "G4PhysicalConstants.hh"
using namespace std;
Material* Material::me =0;
Material::Material() {
// my material
DefineMaterials();
DefineBoundaries();
DefineProperties();
SetMaterialIndexes();
}
Material::~Material(){}
Material* Material::Get() {
if (!me)
me = new Material();
return me;
}
void Material::DefineMaterials() {
G4String name,
symbol;
G4int z,
ncomponents,
natoms;
G4double a,
density,
fractionmass,
temperature,
pressure;
G4Element* H = new G4Element(name="Hydrogen", symbol="H" , z=1, a=1.00794*g/mole );
G4Element* C = new G4Element(name="Carbon", symbol="C", z=6, a=12.0107*g/mole );
G4Element* O = new G4Element(name="Oxygen", symbol="O", z=8, a=15.9994*g/mole );
G4Element* N = new G4Element(name="Nitrogen", symbol="N", z=7, a=14.00067*g/mole);
...
//Paper (Cellulose C6H10O5)
density=0.8*g/cm3;
fPaper = new G4Material(name="Paper",density,ncomponents=3,kStateSolid);
fPaper->AddElement(C,45*perCent);
fPaper->AddElement(H,6*perCent);
fPaper->AddElement(O,49*perCent);
// Vacuum
density = 1.29e-20*g/cm3;
fVacuum = new G4Material(name="Vacuum", density, ncomponents=2);
fVacuum->AddElement(N, fractionmass=0.7);
fVacuum->AddElement(O, fractionmass=0.3);
...
}
void Material::DefineProperties() {
//---------------------------------------------------------------------------------
// Vacuum
//---------------------------------------------------------------------------------
G4MaterialPropertiesTable *myVacuum = new G4MaterialPropertiesTable();
G4double myVacuumRI[2], myVacuumAbs[2];
myVacuumAbs[0] = 1*km;
myVacuumAbs[1] = 1*km;
myVacuumRI[0] = 1 ;
myVacuumRI[1] = 1 ;
G4double myVacuumEne[2] = {0.1*eV , 20*eV} ;
myVacuum->AddProperty("RINDEX", myVacuumEne , myVacuumRI, 2);
myVacuum->AddProperty("ABSLENGTH", myVacuumEne, myVacuumAbs, 2);
fVacuum->SetMaterialPropertiesTable(myVacuum);
...
//---------------------------------------------------------------------------------
// Paper (Frank 6 December 2021)
//---------------------------------------------------------------------------------
G4MaterialPropertiesTable *myPaper = new G4MaterialPropertiesTable();
G4double myPaperRI[2], myPaperAbs[2];
myPaperAbs[0] = 1*mm;
myPaperAbs[1] = 1*mm;
myPaperRI[0] = 2.27 ;
myPaperRI[1] = 2.27 ;
G4double myPaperEne[2] = {0.1*eV , 20*eV} ;
myPaper->AddProperty("RINDEX", myPaperEne , myPaperRI, 2);
myPaper->AddProperty("ABSLENGTH", myPaperEne, myPaperAbs, 2);
fPaper->SetMaterialPropertiesTable(myPaper);
...
}
And then in my detector construction file I have:
#include "DetectorConstruction.hh"
#include "Material.hh"
#include "IO.hh"
#include "Storage.hh"
#include "G4RunManager.hh"
#include "G4NistManager.hh"
#include "G4Box.hh"
#include "G4Tubs.hh"
#include "G4Polyhedra.hh"
#include "G4Transform3D.hh"
#include "G4Cons.hh"
#include "G4VisAttributes.hh"
#include "G4Colour.hh"
#include "G4Orb.hh"
#include "G4Sphere.hh"
#include "G4Trd.hh"
#include "G4LogicalVolume.hh"
#include "G4PVPlacement.hh"
#include "G4SystemOfUnits.hh"
#include "G4OpticalSurface.hh"
#include "G4MaterialPropertiesTable.hh"
#include "G4LogicalBorderSurface.hh"
#include "G4Element.hh"
#include "G4LogicalSkinSurface.hh"
#include "G4Material.hh"
#include "G4ThreeVector.hh"
#include "G4VSolid.hh"
#include "G4BooleanSolid.hh"
#include "G4SubtractionSolid.hh"
#include "G4ExtrudedSolid.hh"
#include "G4RotationMatrix.hh"
#include "G4ReflectedSolid.hh"
#include "G4ReflectionFactory.hh"
#include "G4Transform3D.hh"
#include <vector>
DetectorConstruction::DetectorConstruction()
: G4VUserDetectorConstruction()
{ }
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
DetectorConstruction::~DetectorConstruction()
{ }
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
G4VPhysicalVolume* DetectorConstruction::Construct()
{
...
//---------------------------------------------------------------------------------
// World
//---------------------------------------------------------------------------------
fSolidWorld = new G4Box("World_Solid", 2.*meter, 2.*meter, 2.*meter);
fLogicWorld = new G4LogicalVolume(fSolidWorld, Material::Get()->GetVacuum(), "World_Logic");
fPhysicWorld = new G4PVPlacement(0,
G4ThreeVector(0,0,0),
"World",
fLogicWorld,
NULL,
false,
0);
//---------------------------------------------------------------------------------
// Reflecting Material
//---------------------------------------------------------------------------------
fSolidReflector = new G4Tubs("SolidReflector", innerRadius, 12.5*mm, motref_h, startAngle, spanningAngle) ;
fLogicReflector = new G4LogicalVolume(fSolidReflector, Material::Get()->GetPaper() , "LogicReflector");
fPhysicReflector = new G4PVPlacement(yRot,
G4ThreeVector(motref_x, 0.*cm, 0.*cm),
"ReflectorLayer",
fLogicReflector,
fPhysicWorld,
false,
0,
myCheckOverlap);
fLogicReflector->SetVisAttributes( G4VisAttributes(myGreen) );;
...
DefineSurfaces();
return fPhysicWorld;
}
void DetectorConstruction::DefineSurfaces(){
G4double energy[2] = {3.1*eV, 1.55*eV};
G4double ref[2] = {0.01,0.01};
G4double tran[2] = {0.99,0.99};
G4double ref_refl[2] = {0.9, 0.9};
G4double tran_refl[2] = {0.1, 0.1};
//---------------------------------------------------------------------------------
// Reflector-Vacuum surface
//---------------------------------------------------------------------------------
G4OpticalSurface *fOpReflectorAirSurface = new G4OpticalSurface("OpAirReflectorSurface", glisur, ground, dielectric_dielectric, 0.0);
new G4LogicalBorderSurface("ReflectorAirSurface", fPhysicReflector, fPhysicWorld, fOpReflectorAirSurface);
G4MaterialPropertiesTable *fReflectorAirSurfProp = new G4MaterialPropertiesTable();
fReflectorAirSurfProp->AddProperty("REFLECTIVITY", energy, ref, 2);
fReflectorAirSurfProp->AddProperty("TRANSMITTANCE", energy, tran, 2);
fOpReflectorAirSurface->SetMaterialPropertiesTable( fReflectorAirSurfProp );
//---------------------------------------------------------------------------------
// Vacuum-Reflector surface
//---------------------------------------------------------------------------------
G4OpticalSurface *fOpAirReflectorSurface = new G4OpticalSurface("OpAirReflectorSurface", glisur, ground, dielectric_dielectric, 0.0);
new G4LogicalBorderSurface("AirReflectorSurface", fPhysicWorld, fPhysicReflector, fOpAirReflectorSurface);
G4MaterialPropertiesTable *fAirReflectorSurfProp = new G4MaterialPropertiesTable();
fAirReflectorSurfProp->AddProperty("REFLECTIVITY", energy, ref_refl, 2);
fAirReflectorSurfProp->AddProperty("TRANSMITTANCE", energy, tran_refl, 2);
fOpAirReflectorSurface->SetMaterialPropertiesTable( fAirReflectorSurfProp );
...
}
The optical photons I use are 3 eV (413.3 nm).
Here is a sample of what the pre and post Step points look like at the 0.99 and 1 REFLECTIVITY settings (TRANSMITTANCE is set as 1-REFLECTIVITY).
For REFLECTIVITY = 1.0, starting in the scintillator and heading perpendicular to the surface of the reflecting material (angle of incidence = 0)
Step pre point: (399.3,0,0)
Material Index: YAP
PDG = 0
TrackID: 1
PID:0
E: 3 eV;
Step post point: (398.3,0,0)
Material Index: YAP
PDG = 0
TrackID: 1
PID:0
E: 3 eV;
Step pre point: (398.3,0,0)
Material Index: Paper
PDG = 0
TrackID: 1
PID:0
E: 3 eV;
Step post point: (398.2,0.0064505,-0.000318219)
Material Index: Paper
PDG = 0
TrackID: 1
PID:0
E: 3 eV;
Step pre point: (398.2,0.0064505,-0.000318219)
Material Index: Vacuum
PDG = 0
TrackID: 1
PID:0
E: 3 eV;
Step post point: (-2000,154.702,-7.63186)
Material Index: Vacuum
PDG = 0
TrackID: 1
PID:0
E: 3 eV;
and for REFLECTIVITY = 0.99, other conditions the same as above:
Step pre point: (399.3,0,0)
Material Index: YAP
PDG = 0
TrackID: 1
PID:0
E: 3 eV;
Step post point: (398.3,0,0)
Material Index: YAP
PDG = 0
TrackID: 1
PID:0
E: 3 eV;
Step pre point: (398.3,0,0)
Material Index: Paper
PDG = 0
TrackID: 1
PID:0
E: 3 eV;
Step post point: (398.3,0,0)
Material Index: Paper
PDG = 0
TrackID: 1
PID:0
E: 3 eV;
Step pre point: (398.3,0,0)
Material Index: YAP
PDG = 0
TrackID: 1
PID:0
E: 3 eV;
Step post point: (400.3,0.971597,-0.186861)
Material Index: YAP
PDG = 0
TrackID: 1
PID:0
E: 3 eV;
I have also tried using the “unified” setting vs “glider” and obtain the same result.
Thank you so much again for all of your time and consideration. Please let me know if there is anything else you would like to see.