aboutsummaryrefslogtreecommitdiffstats
path: root/src
diff options
context:
space:
mode:
Diffstat (limited to 'src')
-rw-r--r--src/#G4BeamTestSiHit.cxx#115
-rw-r--r--src/#G4TestEMPhysics.cxx#0
-rw-r--r--src/G4BeamTestDetectorConstruction.cxx200
-rw-r--r--src/G4BeamTestEventAction.cxx56
-rw-r--r--src/G4BeamTestHadronPhysics.cxx.backup412
-rw-r--r--src/G4BeamTestIonPhysics.cxx.backup100
-rw-r--r--src/G4BeamTestRunManager.cxx.backup131
-rw-r--r--src/G4BeamTestSC4SD.cxx90
-rw-r--r--src/G4BeamTestSiSD.cxx125
-rw-r--r--src/G4BeamTestTank.cxx197
-rw-r--r--src/G4BeamTestUserSteppingAction.cxx2
-rw-r--r--src/G4BeamTestUserTrackingAction.cxx57
-rw-r--r--src/G4TankIceSD.cxx166
13 files changed, 385 insertions, 1266 deletions
diff --git a/src/#G4BeamTestSiHit.cxx# b/src/#G4BeamTestSiHit.cxx#
deleted file mode 100644
index 663a03f..0000000
--- a/src/#G4BeamTestSiHit.cxx#
+++ /dev/null
@@ -1,115 +0,0 @@
-#include <iostream>
-#include <iomanip>
-
-#include "G4BeamTestSiHit.h"
-#include "G4UnitsTable.hh"
-#include "G4VVisManager.hh"
-#include "G4Circle.hh"
-#include "G4Colour.hh"
-#include "G4VisAttributes.hh"
-
-G4ThreadLocal G4Allocator<G4BeamTestSiHit>* G4BeamTestSiHitAllocator=0;
-static std::fstream testnew("./testnew.txt", std::ofstream::out);
-
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-
-G4BeamTestSiHit::G4BeamTestSiHit()
- : G4VHit(),
- fTrackID(-1),
-
- fEdep(0.),
- fTime(0.),
- fPos(G4ThreeVector())
-{
- /* G4cout << "opening file" << G4endl; */
- /* testnew.open(testnew_out); */
-}
-
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-
-G4BeamTestSiHit ::~G4BeamTestSiHit() {
- /* G4cout << "closing file" << G4endl; */
- /* testnew.close(); */
-}
-
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-
-G4BeamTestSiHit::G4BeamTestSiHit(const G4BeamTestSiHit& right)
- : G4VHit()
-{
- fTrackID = right.fTrackID;
- // fChamberNb = right.fChamberNb;
- fEdep = right.fEdep;
- fPos = right.fPos;
- fTime = right.fTime;
-}
-
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-
-const G4BeamTestSiHit& G4BeamTestSiHit::operator=(const G4BeamTestSiHit& right)
-{
- fTrackID = right.fTrackID;
- // fChamberNb = right.fChamberNb;
- fEdep = right.fEdep;
- fPos = right.fPos;
- fTime = right.fTime;
-
- return *this;
-}
-
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-
-G4int G4BeamTestSiHit::operator==(const G4BeamTestSiHit& right) const
-{
- return ( this == &right ) ? 1 : 0;
-}
-
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-
-void G4BeamTestSiHit::Draw()
-{
- /* G4VVisManager* pVVisManager = G4VVisManager::GetConcreteInstance(); */
- /* if(pVVisManager) */
- /* { */
- /* G4Circle circle(fPos); */
- /* circle.SetScreenSize(4.); */
- /* circle.SetFillStyle(G4Circle::filled); */
- /* G4Colour colour(1.,0.,0.); */
- /* G4VisAttributes attribs(colour); */
- /* circle.SetVisAttributes(attribs); */
- /* pVVisManager->Draw(circle); */
- /* } */
-}
-
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-
-
-void G4BeamTestSiHit::Print()
-{
- G4cout
- << " trackID: " << fTrackID
- << "Edep: "
- << std::setw(7) << G4BestUnit(fEdep,"Energy")
- << " Position: "
- << std::setw(7) << G4BestUnit(fPos ,"Length")
- << "Time: "
- << std::setw(7) << G4BestUnit(fTime,"Time")
- << G4endl;
-}
-
-void G4BeamTestSiHit::Dataout()
-{
-testnew
- << " trackID: " << fTrackID
- << "Edep: "
- << std::setw(7) << G4BestUnit(fEdep,"Energy")
- << " Position: "
- << std::setw(7) << G4BestUnit( fPos ,"Length")
- << "Time: "
- << std::setw(7) << G4BestUnit( fTime,"Time")
- << G4endl;
-
- }
-
-
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
diff --git a/src/#G4TestEMPhysics.cxx# b/src/#G4TestEMPhysics.cxx#
deleted file mode 100644
index e69de29..0000000
--- a/src/#G4TestEMPhysics.cxx#
+++ /dev/null
diff --git a/src/G4BeamTestDetectorConstruction.cxx b/src/G4BeamTestDetectorConstruction.cxx
index 495330a..2942ec6 100644
--- a/src/G4BeamTestDetectorConstruction.cxx
+++ b/src/G4BeamTestDetectorConstruction.cxx
@@ -13,6 +13,7 @@
#include "G4BeamTestDetectorConstruction.h"
#include "G4BeamTestTank.h"
+#include "G4BeamTestSC4SD.h"
G4BeamTestDetectorConstruction::G4BeamTestDetectorConstruction():
origin_(0, 0, 0), verboseLevel_(0)/* , tankList_(0) */
@@ -27,22 +28,18 @@ G4BeamTestDetectorConstruction::~G4BeamTestDetectorConstruction()
G4VPhysicalVolume* G4BeamTestDetectorConstruction::Construct()
{
- /* if(tankList_.empty()) return NULL; */
-
CreateMaterials();
- /* // World origin in IceCube coordinates */
- /* origin_.set(delaunay.GetOrigin().x(), delaunay.GetOrigin().y(), zSnowBottom + zHalfLength); */
-
// Determine World dimensions
+ // G4double xWorld = 20.0 * CLHEP::m;
G4double xWorld = 4.0 * CLHEP::m;
G4double yWorld = 4.0 * CLHEP::m;
G4double zWorld = 4.0 * CLHEP::m;
// SC4 dimensions
- G4double scinHeight_ = 1 * 2.54 * CLHEP::cm;
- G4double scinWidth_ = 1 * 2.54 * CLHEP::cm;
- G4double scinThickness_ = 1 * 2.54 * CLHEP::cm;
+ G4double scinX_ = 0.25 * 2.54 * CLHEP::cm;
+ G4double scinY_ = 1 * 2.54 * CLHEP::cm;
+ G4double scinZ_ = 1 * 2.54 * CLHEP::cm;
// Create world volume
G4Box* world_box = new G4Box("solid_world", xWorld*0.5, yWorld*0.5, zWorld*0.5);
@@ -51,39 +48,61 @@ G4VPhysicalVolume* G4BeamTestDetectorConstruction::Construct()
G4VPhysicalVolume* worldPhys =
new G4PVPlacement(0, G4ThreeVector(), worldLog, "world", 0, false, 0);
- /* // Instantiation of a set of visualization attributes with cyan colour */
- /* G4VisAttributes * snowVisAtt = new G4VisAttributes(G4Colour(0., 1., 1.)); */
- /* // Assignment of the visualization attributes to the logical volume */
- /* snowLog->SetVisAttributes(snowVisAtt); */
-
- // Install tank
- /* BOOST_FOREACH(G4BeamTestTank* tank, tankList_) */
- /* { */
- tank_->InstallTank(worldPhys, origin_);
- // }
+ // Location calculations
+ G4ThreeVector mwpc1Vec = G4ThreeVector(198.8*CLHEP::cm-xWorld*0.5, 0, 0);
+ G4ThreeVector sc1Vec = mwpc1Vec + G4ThreeVector(-50.0*CLHEP::cm, 0, 0);
+ G4ThreeVector sc2Vec = mwpc1Vec + G4ThreeVector( 80.0*CLHEP::cm, 0, 0);
+ G4ThreeVector sc3Vec = mwpc1Vec + G4ThreeVector(180.0*CLHEP::cm, 0, 0);
+ G4ThreeVector mwpc2Vec = mwpc1Vec + G4ThreeVector(275.3*CLHEP::cm, 0, 0);
+ G4ThreeVector mwpc3Vec = mwpc2Vec + G4ThreeVector(284.6*CLHEP::cm, 0, 0);
+ G4ThreeVector mwpc4Vec = mwpc3Vec + G4ThreeVector(771.8*CLHEP::cm, 0, 0);
+ // G4ThreeVector tankVec = mwpc4Vec + G4ThreeVector(300.0*CLHEP::cm, 0, 0);
+ G4ThreeVector tankVec = G4ThreeVector(0, 0, 0);
+ G4ThreeVector sc4Vec = tankVec + G4ThreeVector(120.0*CLHEP::cm, 0, 0);
+
+ tank_->InstallTank(worldPhys, tankVec);
+
+ // // Define SC1
+ // G4Box* sc1_box = new G4Box("sc1",scinX_*0.5, scinY_*2.0, scinZ_*2.0);
+ // G4LogicalVolume* sc1Log =
+ // new G4LogicalVolume(sc1_box, G4Material::GetMaterial("SC4"), "log_sc1", 0, 0, 0);
+ // G4VPhysicalVolume* sc1Phys =
+ // new G4PVPlacement(0, sc1Vec, sc1Log, "sc1", worldLog, false, 0);
+
+ // // Define SC2
+ // G4Box* sc2_box = new G4Box("sc2",scinX_*0.5, scinY_*2.0, scinZ_*2.0);
+ // G4LogicalVolume* sc2Log =
+ // new G4LogicalVolume(sc2_box, G4Material::GetMaterial("SC4"), "log_sc2", 0, 0, 0);
+ // G4VPhysicalVolume* sc2Phys =
+ // new G4PVPlacement(0, sc2Vec, sc2Log, "sc2", worldLog, false, 0);
+
+ // // Define SC3
+ // G4Box* sc3_box = new G4Box("sc3",scinX_*0.5, scinY_*2.0, scinZ_*2.0);
+ // G4LogicalVolume* sc3Log =
+ // new G4LogicalVolume(sc3_box, G4Material::GetMaterial("SC4"), "log_sc3", 0, 0, 0);
+ // G4VPhysicalVolume* sc3Phys =
+ // new G4PVPlacement(0, sc3Vec, sc3Log, "sc3", worldLog, false, 0);
// Define SC4
- G4Box* sc4_box = new G4Box("sc4",scinHeight_*0.5, scinWidth_*0.5, scinThickness_*0.5);
+ G4Box* sc4_box = new G4Box("sc4",scinX_, scinY_*0.5, scinZ_*0.5);
G4LogicalVolume* sc4Log =
new G4LogicalVolume(sc4_box, G4Material::GetMaterial("SC4"), "log_sc4", 0, 0, 0);
G4VPhysicalVolume* sc4Phys =
- new G4PVPlacement(0, G4ThreeVector(1.2*CLHEP::m,0,0), sc4Log, "sc4", worldLog, false, 0);
+ new G4PVPlacement(0, sc4Vec, sc4Log, "sc4", worldLog, false, 0);
// User limits (energy cutoffs)
- // Do not create photons or electrons below cherenkov threshold
- // See also corresponding UserSpecialCuts in Physicslist !!!!
G4UserLimits* energyLimit = new G4UserLimits();
- /* energyLimit->SetUserMinEkine(280.0 * CLHEP::keV); // Cherenkov threshold of electrons in ice */
energyLimit->SetUserMinEkine(2.26 * CLHEP::eV); // Lower threshold of PMT - 550nm
- // energyLimit->SetUserMaxEkine(3.55 * CLHEP::eV); //upper threshold of PMT - 350nm
worldLog->SetUserLimits(energyLimit);
- /* snowLog->SetUserLimits(energyLimit); */
+ // sc1Log->SetUserLimits(energyLimit);
+ // sc2Log->SetUserLimits(energyLimit);
+ // sc3Log->SetUserLimits(energyLimit);
sc4Log->SetUserLimits(energyLimit);
-
- // G4SDManager* sdManager = G4SDManager::GetSDMpointer();
- // sc4SD_ = new G4BeamTestSC4SD("sc4_SD_", "HitsCollection");
- // sdManager->AddNewDetector(sc4SD_);
- // sc4Log->SetSensitiveDetector(sc4SD_);
+
+ G4SDManager* sdManager = G4SDManager::GetSDMpointer();
+ sc4SD_ = new G4BeamTestSC4SD("sc4_SD_", "SDHitsCollection");
+ sdManager->AddNewDetector(sc4SD_);
+ sc4Log->SetSensitiveDetector(sc4SD_);
return worldPhys;
}
@@ -93,16 +112,11 @@ G4VPhysicalVolume* G4BeamTestDetectorConstruction::Construct()
void G4BeamTestDetectorConstruction::CreateMaterials()
{
CreateAir();
- /* CreateIce(); */
- /* CreateSnow(); */
CreateWater();
CreatePlastic();
- /* CreatePerlite(); */
CreateGlassSphere();
CreateEffectiveDOMMaterial();
CreateSC4();
-
- //if(verboseLevel_>0) G4cout << *G4Material::GetMaterialTable() << G4endl;
}
/*****************************************************************/
@@ -122,16 +136,86 @@ void G4BeamTestDetectorConstruction::CreateWater()
water->AddElement(nistManager->FindOrBuildElement("H"), 2);
water->AddElement(nistManager->FindOrBuildElement("O"), 1);
- // pmt spectral response 300-650nm
- // https://docushare.icecube.wisc.edu/dsweb/Get/Document-6637/R7081-02%20data%20sheet.pdf<Paste>
- // TODO(shivesh): add more properties?
- const G4int water_bins = 2;
- // G4double water_ephot[water_bins] = {1.91 * CLHEP::eV, 4.13 * CLHEP::eV};
- G4double water_ephot[water_bins] = {0.1 * CLHEP::eV, 10 * CLHEP::eV};
- G4double water_refr[water_bins] = {1.33, 1.33};
+ // const G4int water_bins = 2;
+ // G4double water_ephot[water_bins] = {0.1 * CLHEP::eV, 10 * CLHEP::eV};
+ // G4double water_refr[water_bins] = {1.33, 1.33};
+
+ //From SFDETSIM water absorption
+ const G4int NUMENTRIES_water=60;
+
+ G4double ENERGY_water[NUMENTRIES_water] =
+ { 1.56962e-09*GeV, 1.58974e-09*GeV, 1.61039e-09*GeV, 1.63157e-09*GeV,
+ 1.65333e-09*GeV, 1.67567e-09*GeV, 1.69863e-09*GeV, 1.72222e-09*GeV,
+ 1.74647e-09*GeV, 1.77142e-09*GeV,1.7971e-09*GeV, 1.82352e-09*GeV,
+ 1.85074e-09*GeV, 1.87878e-09*GeV, 1.90769e-09*GeV, 1.93749e-09*GeV,
+ 1.96825e-09*GeV, 1.99999e-09*GeV, 2.03278e-09*GeV, 2.06666e-09*GeV,
+ 2.10169e-09*GeV, 2.13793e-09*GeV, 2.17543e-09*GeV, 2.21428e-09*GeV,
+ 2.25454e-09*GeV, 2.29629e-09*GeV, 2.33962e-09*GeV, 2.38461e-09*GeV,
+ 2.43137e-09*GeV, 2.47999e-09*GeV, 2.53061e-09*GeV, 2.58333e-09*GeV,
+ 2.63829e-09*GeV, 2.69565e-09*GeV, 2.75555e-09*GeV, 2.81817e-09*GeV,
+ 2.88371e-09*GeV, 2.95237e-09*GeV, 3.02438e-09*GeV, 3.09999e-09*GeV,
+ 3.17948e-09*GeV, 3.26315e-09*GeV, 3.35134e-09*GeV, 3.44444e-09*GeV,
+ 3.54285e-09*GeV, 3.64705e-09*GeV, 3.75757e-09*GeV, 3.87499e-09*GeV,
+ 3.99999e-09*GeV, 4.13332e-09*GeV, 4.27585e-09*GeV, 4.42856e-09*GeV,
+ 4.59258e-09*GeV, 4.76922e-09*GeV, 4.95999e-09*GeV, 5.16665e-09*GeV,
+ 5.39129e-09*GeV, 5.63635e-09*GeV, 5.90475e-09*GeV, 6.19998e-09*GeV };
+
+ // M Fechner : new ; define the water refraction index using refsg.F
+ //from skdetsim using the whole range.
+ G4double RINDEX1[NUMENTRIES_water] =
+ {1.32885, 1.32906, 1.32927, 1.32948, 1.3297, 1.32992, 1.33014,
+ 1.33037, 1.3306, 1.33084, 1.33109, 1.33134, 1.3316, 1.33186, 1.33213,
+ 1.33241, 1.3327, 1.33299, 1.33329, 1.33361, 1.33393, 1.33427, 1.33462,
+ 1.33498, 1.33536, 1.33576, 1.33617, 1.3366, 1.33705, 1.33753, 1.33803,
+ 1.33855, 1.33911, 1.3397, 1.34033, 1.341, 1.34172, 1.34248, 1.34331,
+ 1.34419, 1.34515, 1.3462, 1.34733, 1.34858, 1.34994, 1.35145, 1.35312,
+ 1.35498, 1.35707, 1.35943, 1.36211, 1.36518, 1.36872, 1.37287, 1.37776,
+ 1.38362, 1.39074, 1.39956, 1.41075, 1.42535};
+
+ G4double ABWFF = 1.0;
+ //T. Akiri: Values from Skdetsim
+ G4double ABSORPTION_water[NUMENTRIES_water] =
+ {
+ 16.1419*cm*ABWFF, 18.278*cm*ABWFF, 21.0657*cm*ABWFF, 24.8568*cm*ABWFF, 30.3117*cm*ABWFF,
+ 38.8341*cm*ABWFF, 54.0231*cm*ABWFF, 81.2306*cm*ABWFF, 120.909*cm*ABWFF, 160.238*cm*ABWFF,
+ 193.771*cm*ABWFF, 215.017*cm*ABWFF, 227.747*cm*ABWFF, 243.85*cm*ABWFF, 294.036*cm*ABWFF,
+ 321.647*cm*ABWFF, 342.81*cm*ABWFF, 362.827*cm*ABWFF, 378.041*cm*ABWFF, 449.378*cm*ABWFF,
+ 739.434*cm*ABWFF, 1114.23*cm*ABWFF, 1435.56*cm*ABWFF, 1611.06*cm*ABWFF, 1764.18*cm*ABWFF,
+ 2100.95*cm*ABWFF, 2292.9*cm*ABWFF, 2431.33*cm*ABWFF, 3053.6*cm*ABWFF, 4838.23*cm*ABWFF,
+ 6539.65*cm*ABWFF, 7682.63*cm*ABWFF, 9137.28*cm*ABWFF, 12220.9*cm*ABWFF, 15270.7*cm*ABWFF,
+ 19051.5*cm*ABWFF, 23671.3*cm*ABWFF, 29191.1*cm*ABWFF, 35567.9*cm*ABWFF, 42583*cm*ABWFF,
+ 49779.6*cm*ABWFF, 56465.3*cm*ABWFF, 61830*cm*ABWFF, 65174.6*cm*ABWFF, 66143.7*cm*ABWFF,
+ 64820*cm*ABWFF, 61635*cm*ABWFF, 57176.2*cm*ABWFF, 52012.1*cm*ABWFF, 46595.7*cm*ABWFF,
+ 41242.1*cm*ABWFF, 36146.3*cm*ABWFF, 31415.4*cm*ABWFF, 27097.8*cm*ABWFF, 23205.7*cm*ABWFF,
+ 19730.3*cm*ABWFF, 16651.6*cm*ABWFF, 13943.6*cm*ABWFF, 11578.1*cm*ABWFF, 9526.13*cm*ABWFF
+ };
+
+ G4double RAYFF = 0.625;
+ //T. Akiri: Values from Skdetsim
+ G4double RAYLEIGH_water[NUMENTRIES_water] = {
+ 386929*cm*RAYFF, 366249*cm*RAYFF, 346398*cm*RAYFF, 327355*cm*RAYFF, 309097*cm*RAYFF,
+ 291603*cm*RAYFF, 274853*cm*RAYFF, 258825*cm*RAYFF, 243500*cm*RAYFF, 228856*cm*RAYFF,
+ 214873*cm*RAYFF, 201533*cm*RAYFF, 188816*cm*RAYFF, 176702*cm*RAYFF, 165173*cm*RAYFF,
+ 154210*cm*RAYFF, 143795*cm*RAYFF, 133910*cm*RAYFF, 124537*cm*RAYFF, 115659*cm*RAYFF,
+ 107258*cm*RAYFF, 99318.2*cm*RAYFF, 91822.2*cm*RAYFF, 84754*cm*RAYFF, 78097.3*cm*RAYFF,
+ 71836.5*cm*RAYFF, 65956*cm*RAYFF, 60440.6*cm*RAYFF, 55275.4*cm*RAYFF, 50445.6*cm*RAYFF,
+ 45937*cm*RAYFF, 41735.2*cm*RAYFF, 37826.6*cm*RAYFF, 34197.6*cm*RAYFF, 30834.9*cm*RAYFF,
+ 27725.4*cm*RAYFF, 24856.6*cm*RAYFF, 22215.9*cm*RAYFF, 19791.3*cm*RAYFF, 17570.9*cm*RAYFF,
+ 15543*cm*RAYFF, 13696.6*cm*RAYFF, 12020.5*cm*RAYFF, 10504.1*cm*RAYFF, 9137.15*cm*RAYFF,
+ 7909.45*cm*RAYFF, 6811.3*cm*RAYFF, 5833.25*cm*RAYFF, 4966.2*cm*RAYFF, 4201.36*cm*RAYFF,
+ 3530.28*cm*RAYFF, 2944.84*cm*RAYFF, 2437.28*cm*RAYFF, 2000.18*cm*RAYFF, 1626.5*cm*RAYFF,
+ 1309.55*cm*RAYFF, 1043.03*cm*RAYFF, 821.016*cm*RAYFF, 637.97*cm*RAYFF, 488.754*cm*RAYFF
+ };
+
+ // G4MaterialPropertiesTable *mpt_water = new G4MaterialPropertiesTable ();
+ // mpt_water->AddProperty ("RINDEX", water_ephot, water_refr, water_bins);
+ // water->SetMaterialPropertiesTable(mpt_water);
G4MaterialPropertiesTable *mpt_water = new G4MaterialPropertiesTable ();
- mpt_water->AddProperty ("RINDEX", water_ephot, water_refr, water_bins);
+ mpt_water->AddProperty ("RINDEX", ENERGY_water, RINDEX1, NUMENTRIES_water);
+ mpt_water->AddProperty("ABSLENGTH",ENERGY_water, ABSORPTION_water, NUMENTRIES_water);
+ mpt_water->AddProperty("RAYLEIGH",ENERGY_water,RAYLEIGH_water,NUMENTRIES_water);
water->SetMaterialPropertiesTable(mpt_water);
+
}
/*****************************************************************/
@@ -139,52 +223,32 @@ void G4BeamTestDetectorConstruction::CreateWater()
void G4BeamTestDetectorConstruction::CreatePlastic()
{
G4NistManager* nistManager = G4NistManager::Instance();
- /*
- G4Material* plastic = new G4Material("Plastic", 1.5 * CLHEP::g / CLHEP::cm3, 3, kStateSolid);
- plastic->AddElement(nistManager->FindOrBuildElement("H"), 20);
- plastic->AddElement(nistManager->FindOrBuildElement("C"), 10);
- plastic->AddElement(nistManager->FindOrBuildElement("O"), 5);
- */
-
- // POM
G4Material* plastic = new G4Material("Plastic", 1.425 * CLHEP::g / CLHEP::cm3, 3, kStateSolid);
plastic->AddElement(nistManager->FindOrBuildElement("H"), 2);
plastic->AddElement(nistManager->FindOrBuildElement("C"), 1);
plastic->AddElement(nistManager->FindOrBuildElement("O"), 1);
-
- //nistManager->FindOrBuildMaterial("G4_POLYOXYMETHYLENE");
-
}
/*****************************************************************/
void G4BeamTestDetectorConstruction::CreateGlassSphere()
{
- // Elemental composition not exact for detailed line-up look at
- // docushare collection 690 -> XRF and ICP Analysis of Pressure Sphere Glass
-
- // 20 lbs. weight = 9072 g
- // 6.5" outer radius & 0.5" thickness = 4024 cm3
- G4NistManager* nistManager = G4NistManager::Instance();
- // G4Material* glass = new G4Material("Glass", 2.254 * CLHEP::g / CLHEP::cm3, 2, kStateSolid);
- // glass->AddElement(nistManager->FindOrBuildElement("Si"), 1);
- // glass->AddElement(nistManager->FindOrBuildElement("O"), 2);
-
// http://hypernews.slac.stanford.edu/HyperNews/geant4/get/AUX/2013/03/11/12.39-85121-chDetectorConstruction.cc
+ G4NistManager* nistManager = G4NistManager::Instance();
// Define elements for all materials not found in the NIST database
G4Element* Si = nistManager->FindOrBuildElement("Si");
G4Element* B = nistManager->FindOrBuildElement("B");
G4Element* O = nistManager->FindOrBuildElement("O");
G4Element* Na = nistManager->FindOrBuildElement("Na");
G4Element* Al = nistManager->FindOrBuildElement("Al");
- G4Element* K = nistManager->FindOrBuildElement("K");
+ G4Element* K = nistManager->FindOrBuildElement("K");
G4double density;
G4int ncomponents;
G4double fractionmass;
G4Material* glass = new G4Material("Glass", density= 2.23*CLHEP::g/CLHEP::cm3, ncomponents=6);
glass->AddElement(B, fractionmass=0.040064);
- glass->AddElement(O, fractionmass=0.539562);
+ glass->AddElement(O, fractionmass=0.539562);
glass->AddElement(Na, fractionmass=0.028191);
glass->AddElement(Al, fractionmass=0.011644);
glass->AddElement(Si, fractionmass=0.377220);
@@ -192,7 +256,6 @@ void G4BeamTestDetectorConstruction::CreateGlassSphere()
// pmt spectral response 300-650nm
// https://docushare.icecube.wisc.edu/dsweb/Get/Document-6637/R7081-02%20data%20sheet.pdf<Paste>
- // TODO(shivesh): add more properties?
const G4int glass_bins = 2;
// G4double glass_ephot[glass_bins] = {1.91 * CLHEP::eV, 4.13 * CLHEP::eV};
G4double glass_ephot[glass_bins] = {0.1 * CLHEP::eV, 10 * CLHEP::eV};
@@ -216,7 +279,7 @@ void G4BeamTestDetectorConstruction::CreateEffectiveDOMMaterial()
// const G4int glass_bins = 2;
// G4double glass_ephot[glass_bins] = {1.91 * CLHEP::eV, 4.13 * CLHEP::eV};
- // glass
+ // glass
// G4double glass_ephot[glass_bins] = {0.1 * CLHEP::eV, 10 * CLHEP::eV};
// G4double glass_refr[glass_bins] = {1.47, 1.47};
// G4MaterialPropertiesTable *mpt_glass = new G4MaterialPropertiesTable ();
@@ -227,7 +290,6 @@ void G4BeamTestDetectorConstruction::CreateEffectiveDOMMaterial()
void G4BeamTestDetectorConstruction::CreateSC4()
{
G4NistManager* nistManager = G4NistManager::Instance();
- // POM
G4Material* plastic = new G4Material("SC4", 1.425 * CLHEP::g / CLHEP::cm3, 3, kStateSolid);
plastic->AddElement(nistManager->FindOrBuildElement("H"), 2);
plastic->AddElement(nistManager->FindOrBuildElement("C"), 1);
diff --git a/src/G4BeamTestEventAction.cxx b/src/G4BeamTestEventAction.cxx
index 53402c9..bf6b09b 100644
--- a/src/G4BeamTestEventAction.cxx
+++ b/src/G4BeamTestEventAction.cxx
@@ -9,6 +9,7 @@
#include "G4Event.hh"
#include "G4RunManager.hh"
+#include "G4Interface.h"
#include "G4BeamTestSiSD.h"
#include "G4BeamTestSiHit.h"
#include "G4BeamTestEventAction.h"
@@ -16,10 +17,11 @@
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
G4BeamTestEventAction::G4BeamTestEventAction()
: G4UserEventAction(),
-
- SiCollID(0)
+
+ SiCollID(0),
+ SC4CollID(0)
{
-}
+}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
@@ -29,75 +31,67 @@ G4BeamTestEventAction::~G4BeamTestEventAction()
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
void G4BeamTestEventAction::BeginOfEventAction(const G4Event* event )
-{
-
-
+{
G4SDManager * SDman = G4SDManager::GetSDMpointer();
SDman->ListTree();
- if(SiCollID<0)
- {
- G4String colNam;
- SiCollID = SDman->GetCollectionID(colNam="G4BeamTestSiSDCollection");
-
- }
-
-
+ if(SiCollID<0) SiCollID = SDman->GetCollectionID("ice_SD_");
+ if(SC4CollID<0) SC4CollID = SDman->GetCollectionID("sc4_SD_");
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
void G4BeamTestEventAction::EndOfEventAction(const G4Event* event)
-{
+{
G4cout << ">>> Summary of Event " << event->GetEventID() << G4endl;
// testnew << ">>> Summary of Event " << event->GetEventID() << G4endl;
G4cout << SiCollID << G4endl;
+ G4cout << SC4CollID << G4endl;
if(SiCollID<0) return;
+ if(SC4CollID<0) return;
G4HCofThisEvent* HCE = event->GetHCofThisEvent();
G4BeamTestSiHitsCollection* SiHC = 0;
+ G4BeamTestSiHitsCollection* SC4HC = 0;
- if(HCE)
- {
- SiHC = (G4BeamTestSiHitsCollection*)(HCE->GetHC(SiCollID));
-
- }
-
+ G4cout << "# collections = " << HCE->GetNumberOfCollections() << G4endl;
+ if(HCE) {
+ SiHC = (G4BeamTestSiHitsCollection*)(HCE->GetHC(SiCollID));
+ SC4HC = (G4BeamTestSiHitsCollection*)(HCE->GetHC(1));
+ }
- if(SiHC)
+ if(SiHC && SC4HC)
{
+ // SC4 veto
+ int sc4_hit = SC4HC->entries();
+ if (sc4_hit < 1) return;
+
int n_hit = SiHC->entries();
- // testnew << std::flush;
G4cout << G4endl;
// G4cout << "Si hits " <<
// "--------------------------------------------------------------"
// << G4endl;
- G4cout << n_hit << " hits are stored in G4BeamTestSiHitsCollection."
+ G4cout << n_hit << " hits are stored in ice_SD_"
<< G4endl;
/* G4cout << "List of hits in tracker" << G4endl; */
- // testnew << G4endl;
- // testnew << "Si hits " <<
// "--------------------------------------------------------------"
- testnew << n_hit << " hits are stored in G4BeamTestSiHitsCollection."
+ testnew << n_hit << " hits are stored in ice_SD_"
<< G4endl;
- /* testnew << "List of hits in tracker" << G4endl; */
for(int i=0;i<n_hit;i++)
{
/* (*SiHC)[i]->Print(); */
(*SiHC)[i]->Dataout();
}
-
// G4cout << "sid + " << SiCollID << G4endl;
// testnew << "sid + " << SiCollID << G4endl;
testnew << std::flush;
- }
+ }
}
-
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
diff --git a/src/G4BeamTestHadronPhysics.cxx.backup b/src/G4BeamTestHadronPhysics.cxx.backup
deleted file mode 100644
index 1e37ed5..0000000
--- a/src/G4BeamTestHadronPhysics.cxx.backup
+++ /dev/null
@@ -1,412 +0,0 @@
-#include <iomanip>
-
-#include <globals.hh>
-#include <G4ios.hh>
-#include <G4ParticleDefinition.hh>
-#include <G4ParticleTable.hh>
-#include <G4MesonConstructor.hh>
-#include <G4BaryonConstructor.hh>
-#include <G4ShortLivedConstructor.hh>
-#include <G4ProcessManager.hh>
-
-#include "G4BeamTestHadronPhysics.h"
-
-G4BeamTestHadronPhysics::G4BeamTestHadronPhysics()
- : G4VPhysicsConstructor("hadron") {
-}
-
-/********************************************************************/
-
-G4BeamTestHadronPhysics::~G4BeamTestHadronPhysics()
-{
- delete stringDecay_;
-}
-
-/********************************************************************/
-
-void G4BeamTestHadronPhysics::ConstructParticle()
-{
-
- // Construct all mesons
- G4MesonConstructor mesonConstructor;
- mesonConstructor.ConstructParticle();
-
- // Construct all barions
- G4BaryonConstructor baryonConstructor;
- baryonConstructor.ConstructParticle();
-
- // Construct resonaces and quarks
- G4ShortLivedConstructor shortLivedConstructor;
- shortLivedConstructor.ConstructParticle();
-
-}
-
-/********************************************************************/
-
-void G4BeamTestHadronPhysics::ConstructProcess()
-{
- G4ProcessManager *pManager = 0;
-
- // G4cout << "" << G4endl;
- // G4cout << "You are using the G4BeamTestHadronPhysics" << G4endl;
- // G4cout << " - Note that this hadronic physics list is not optimized for any particular usage" << G4endl;
- // G4cout << " - If you wish to have a starting point tailored for a particular area of work," << G4endl;
- // G4cout << " please use one of the available physics lists by use-case." << G4endl;
- // G4cout << " More information can also be found from the Geant4 HyperNews." << G4endl;
- // G4cout << "" << G4endl;
-
- // Elastic Process
- elasticModel_ = new G4LElastic();
- elasticProcess_.RegisterMe(elasticModel_);
-
- // pi+ and pi-
- preEquilib_ = new G4PreCompoundModel(&excitationHandler_);
- cascade_.SetDeExcitation(preEquilib_);
- theoModel_.SetTransport(&cascade_);
- theoModel_.SetHighEnergyGenerator(&stringModel_);
- stringDecay_ = new G4ExcitedStringDecay(&fragmentation_);
- stringModel_.SetFragmentationModel(stringDecay_);
- theoModel_.SetMinEnergy(15*CLHEP::GeV);
- theoModel_.SetMaxEnergy(100*CLHEP::TeV);
-
- // PionPlus
- pManager = G4PionPlus::PionPlus()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- lePionPlusModel_ = new G4LEPionPlusInelastic();
- pionPlusInelastic_.RegisterMe(lePionPlusModel_);
- pionPlusInelastic_.RegisterMe(&theoModel_);
- pManager->AddDiscreteProcess(&pionPlusInelastic_);
-
- pManager->AddProcess(&pionPlusIonisation_, ordInActive,2, 2);
-
- pManager->AddProcess(&pionPlusMult_);
- pManager->SetProcessOrdering(&pionPlusMult_, idxAlongStep, 1);
- pManager->SetProcessOrdering(&pionPlusMult_, idxPostStep, 1);
-
- // PionMinus
- pManager = G4PionMinus::PionMinus()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- lePionMinusModel_ = new G4LEPionMinusInelastic();
- pionMinusInelastic_.RegisterMe(lePionMinusModel_);
- pionMinusInelastic_.RegisterMe(&theoModel_);
- pManager->AddDiscreteProcess(&pionMinusInelastic_);
-
- pManager->AddProcess(&pionMinusIonisation_, ordInActive,2, 2);
-
- pManager->AddProcess(&pionMinusMult_);
- pManager->SetProcessOrdering(&pionMinusMult_, idxAlongStep, 1);
- pManager->SetProcessOrdering(&pionMinusMult_, idxPostStep, 1);
-
- pManager->AddRestProcess(&pionMinusAbsorption_, ordDefault);
-
- // KaonPlus
- pManager = G4KaonPlus::KaonPlus()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leKaonPlusModel_ = new G4LEKaonPlusInelastic();
- heKaonPlusModel_ = new G4HEKaonPlusInelastic();
- kaonPlusInelastic_.RegisterMe(leKaonPlusModel_);
- kaonPlusInelastic_.RegisterMe(&theoModel_);
- pManager->AddDiscreteProcess(&kaonPlusInelastic_);
-
- pManager->AddProcess(&kaonPlusIonisation_, ordInActive,2, 2);
-
- pManager->AddProcess(&kaonPlusMult_);
- pManager->SetProcessOrdering(&kaonPlusMult_, idxAlongStep, 1);
- pManager->SetProcessOrdering(&kaonPlusMult_, idxPostStep, 1);
-
- // KaonMinus
- pManager = G4KaonMinus::KaonMinus()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leKaonMinusModel_ = new G4LEKaonMinusInelastic();
- heKaonMinusModel_ = new G4HEKaonMinusInelastic();
- kaonMinusInelastic_.RegisterMe(leKaonMinusModel_);
- kaonMinusInelastic_.RegisterMe(heKaonMinusModel_);
- pManager->AddDiscreteProcess(&kaonMinusInelastic_);
-
- pManager->AddProcess(&kaonMinusIonisation_, ordInActive,2, 2);
-
- pManager->AddProcess(&kaonMinusMult_);
- pManager->SetProcessOrdering(&kaonMinusMult_, idxAlongStep, 1);
- pManager->SetProcessOrdering(&kaonMinusMult_, idxPostStep, 1);
-
- pManager->AddRestProcess(&kaonMinusAbsorption_, ordDefault);
-
- // KaonZeroL
- pManager = G4KaonZeroLong::KaonZeroLong()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leKaonZeroLModel_ = new G4LEKaonZeroLInelastic();
- heKaonZeroLModel_ = new G4HEKaonZeroInelastic();
- kaonZeroLInelastic_.RegisterMe(leKaonZeroLModel_);
- kaonZeroLInelastic_.RegisterMe(heKaonZeroLModel_);
- pManager->AddDiscreteProcess(&kaonZeroLInelastic_);
-
- // KaonZeroS
- pManager = G4KaonZeroShort::KaonZeroShort()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leKaonZeroSModel_ = new G4LEKaonZeroSInelastic();
- heKaonZeroSModel_ = new G4HEKaonZeroInelastic();
- kaonZeroSInelastic_.RegisterMe(leKaonZeroSModel_);
- kaonZeroSInelastic_.RegisterMe(heKaonZeroSModel_);
- pManager->AddDiscreteProcess(&kaonZeroSInelastic_);
-
- // Proton
- pManager = G4Proton::Proton()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leProtonModel_ = new G4LEProtonInelastic();
- heProtonModel_ = new G4HEProtonInelastic();
- protonInelastic_.RegisterMe(leProtonModel_);
- protonInelastic_.RegisterMe(&theoModel_);
- pManager->AddDiscreteProcess(&protonInelastic_);
-
- pManager->AddProcess(&protonIonisation_, ordInActive,2, 2);
-
- pManager->AddProcess(&protonMult_);
- pManager->SetProcessOrdering(&protonMult_, idxAlongStep, 1);
- pManager->SetProcessOrdering(&protonMult_, idxPostStep, 1);
-
- // anti-Proton
- pManager = G4AntiProton::AntiProton()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leAntiProtonModel_ = new G4LEAntiProtonInelastic();
- heAntiProtonModel_ = new G4HEAntiProtonInelastic();
- antiProtonInelastic_.RegisterMe(leAntiProtonModel_);
- antiProtonInelastic_.RegisterMe(heAntiProtonModel_);
- pManager->AddDiscreteProcess(&antiProtonInelastic_);
-
- pManager->AddProcess(&antiProtonIonisation_, ordInActive,2, 2);
-
- pManager->AddProcess(&antiProtonMult_);
- pManager->SetProcessOrdering(&antiProtonMult_, idxAlongStep, 1);
- pManager->SetProcessOrdering(&antiProtonMult_, idxPostStep, 1);
-
- pManager->AddRestProcess(&antiProtonAnnihilation_);
-
- // Neutron
- pManager = G4Neutron::Neutron()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leNeutronModel_ = new G4LENeutronInelastic();
- heNeutronModel_ = new G4HENeutronInelastic();
- neutronInelastic_.RegisterMe(leNeutronModel_);
- neutronInelastic_.RegisterMe(&theoModel_);
- pManager->AddDiscreteProcess(&neutronInelastic_);
-
- neutronFissionModel_ = new G4LFission();
- neutronFission_.RegisterMe(neutronFissionModel_);
- pManager->AddDiscreteProcess(&neutronFission_);
-
- neutronCaptureModel_ = new G4LCapture();
- neutronCapture_.RegisterMe(neutronCaptureModel_);
- pManager->AddDiscreteProcess(&neutronCapture_);
-
- // AntiNeutron
- pManager = G4AntiNeutron::AntiNeutron()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leAntiNeutronModel_ = new G4LEAntiNeutronInelastic();
- heAntiNeutronModel_ = new G4HEAntiNeutronInelastic();
- antiNeutronInelastic_.RegisterMe(leAntiNeutronModel_);
- antiNeutronInelastic_.RegisterMe(heAntiNeutronModel_);
- pManager->AddDiscreteProcess(&antiNeutronInelastic_);
-
- pManager->AddRestProcess(&antiNeutronAnnihilation_);
-
- // Lambda
- pManager = G4Lambda::Lambda()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leLambdaModel_ = new G4LELambdaInelastic();
- heLambdaModel_ = new G4HELambdaInelastic();
- lambdaInelastic_.RegisterMe(leLambdaModel_);
- lambdaInelastic_.RegisterMe(heLambdaModel_);
- pManager->AddDiscreteProcess(&lambdaInelastic_);
-
- // AntiLambda
- pManager = G4AntiLambda::AntiLambda()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leAntiLambdaModel_ = new G4LEAntiLambdaInelastic();
- heAntiLambdaModel_ = new G4HEAntiLambdaInelastic();
- antiLambdaInelastic_.RegisterMe(leAntiLambdaModel_);
- antiLambdaInelastic_.RegisterMe(heAntiLambdaModel_);
- pManager->AddDiscreteProcess(&antiLambdaInelastic_);
-
- // SigmaMinus
- pManager = G4SigmaMinus::SigmaMinus()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leSigmaMinusModel_ = new G4LESigmaMinusInelastic();
- heSigmaMinusModel_ = new G4HESigmaMinusInelastic();
- sigmaMinusInelastic_.RegisterMe(leSigmaMinusModel_);
- sigmaMinusInelastic_.RegisterMe(heSigmaMinusModel_);
- pManager->AddDiscreteProcess(&sigmaMinusInelastic_);
-
- pManager->AddProcess(&sigmaMinusIonisation_, ordInActive,2, 2);
-
- pManager->AddProcess(&sigmaMinusMult_);
- pManager->SetProcessOrdering(&sigmaMinusMult_, idxAlongStep, 1);
- pManager->SetProcessOrdering(&sigmaMinusMult_, idxPostStep, 1);
-
- // anti-SigmaMinus
- pManager = G4AntiSigmaMinus::AntiSigmaMinus()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leAntiSigmaMinusModel_ = new G4LEAntiSigmaMinusInelastic();
- heAntiSigmaMinusModel_ = new G4HEAntiSigmaMinusInelastic();
- antiSigmaMinusInelastic_.RegisterMe(leAntiSigmaMinusModel_);
- antiSigmaMinusInelastic_.RegisterMe(heAntiSigmaMinusModel_);
- pManager->AddDiscreteProcess(&antiSigmaMinusInelastic_);
-
- pManager->AddProcess(&antiSigmaMinusIonisation_, ordInActive,2, 2);
-
- pManager->AddProcess(&antiSigmaMinusMult_);
- pManager->SetProcessOrdering(&antiSigmaMinusMult_, idxAlongStep, 1);
- pManager->SetProcessOrdering(&antiSigmaMinusMult_, idxPostStep, 1);
-
- // SigmaPlus
- pManager = G4SigmaPlus::SigmaPlus()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leSigmaPlusModel_ = new G4LESigmaPlusInelastic();
- heSigmaPlusModel_ = new G4HESigmaPlusInelastic();
- sigmaPlusInelastic_.RegisterMe(leSigmaPlusModel_);
- sigmaPlusInelastic_.RegisterMe(heSigmaPlusModel_);
- pManager->AddDiscreteProcess(&sigmaPlusInelastic_);
-
- pManager->AddProcess(&sigmaPlusIonisation_, ordInActive,2, 2);
-
- pManager->AddProcess(&sigmaPlusMult_);
- pManager->SetProcessOrdering(&sigmaPlusMult_, idxAlongStep, 1);
- pManager->SetProcessOrdering(&sigmaPlusMult_, idxPostStep, 1);
-
- // anti-SigmaPlus
- pManager = G4AntiSigmaPlus::AntiSigmaPlus()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leAntiSigmaPlusModel_ = new G4LEAntiSigmaPlusInelastic();
- heAntiSigmaPlusModel_ = new G4HEAntiSigmaPlusInelastic();
- antiSigmaPlusInelastic_.RegisterMe(leAntiSigmaPlusModel_);
- antiSigmaPlusInelastic_.RegisterMe(heAntiSigmaPlusModel_);
- pManager->AddDiscreteProcess(&antiSigmaPlusInelastic_);
-
- pManager->AddProcess(&antiSigmaPlusIonisation_, ordInActive,2, 2);
-
- pManager->AddProcess(&antiSigmaPlusMult_);
- pManager->SetProcessOrdering(&antiSigmaPlusMult_, idxAlongStep, 1);
- pManager->SetProcessOrdering(&antiSigmaPlusMult_, idxPostStep, 1);
-
- // XiMinus
- pManager = G4XiMinus::XiMinus()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leXiMinusModel_ = new G4LEXiMinusInelastic();
- heXiMinusModel_ = new G4HEXiMinusInelastic();
- xiMinusInelastic_.RegisterMe(leXiMinusModel_);
- xiMinusInelastic_.RegisterMe(heXiMinusModel_);
- pManager->AddDiscreteProcess(&xiMinusInelastic_);
-
- pManager->AddProcess(&xiMinusIonisation_, ordInActive,2, 2);
-
- pManager->AddProcess(&xiMinusMult_);
- pManager->SetProcessOrdering(&xiMinusMult_, idxAlongStep, 1);
- pManager->SetProcessOrdering(&xiMinusMult_, idxPostStep, 1);
-
- // anti-XiMinus
- pManager = G4AntiXiMinus::AntiXiMinus()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leAntiXiMinusModel_ = new G4LEAntiXiMinusInelastic();
- heAntiXiMinusModel_ = new G4HEAntiXiMinusInelastic();
- antiXiMinusInelastic_.RegisterMe(leAntiXiMinusModel_);
- antiXiMinusInelastic_.RegisterMe(heAntiXiMinusModel_);
- pManager->AddDiscreteProcess(&antiXiMinusInelastic_);
-
- pManager->AddProcess(&antiXiMinusIonisation_, ordInActive,2, 2);
-
- pManager->AddProcess(&antiXiMinusMult_);
- pManager->SetProcessOrdering(&antiXiMinusMult_, idxAlongStep, 1);
- pManager->SetProcessOrdering(&antiXiMinusMult_, idxPostStep, 1);
-
- // XiZero
- pManager = G4XiZero::XiZero()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leXiZeroModel_ = new G4LEXiZeroInelastic();
- heXiZeroModel_ = new G4HEXiZeroInelastic();
- xiZeroInelastic_.RegisterMe(leXiZeroModel_);
- xiZeroInelastic_.RegisterMe(heXiZeroModel_);
- pManager->AddDiscreteProcess(&xiZeroInelastic_);
-
- // anti-XiZero
- pManager = G4AntiXiZero::AntiXiZero()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leAntiXiZeroModel_ = new G4LEAntiXiZeroInelastic();
- heAntiXiZeroModel_ = new G4HEAntiXiZeroInelastic();
- antiXiZeroInelastic_.RegisterMe(leAntiXiZeroModel_);
- antiXiZeroInelastic_.RegisterMe(heAntiXiZeroModel_);
- pManager->AddDiscreteProcess(&antiXiZeroInelastic_);
-
- // OmegaMinus
- pManager = G4OmegaMinus::OmegaMinus()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leOmegaMinusModel_ = new G4LEOmegaMinusInelastic();
- heOmegaMinusModel_ = new G4HEOmegaMinusInelastic();
- omegaMinusInelastic_.RegisterMe(leOmegaMinusModel_);
- omegaMinusInelastic_.RegisterMe(heOmegaMinusModel_);
- pManager->AddDiscreteProcess(&omegaMinusInelastic_);
-
- pManager->AddProcess(&omegaMinusIonisation_, ordInActive,2, 2);
-
- pManager->AddProcess(&omegaMinusMult_);
- pManager->SetProcessOrdering(&omegaMinusMult_, idxAlongStep, 1);
- pManager->SetProcessOrdering(&omegaMinusMult_, idxPostStep, 1);
-
- // anti-OmegaMinus
- pManager = G4AntiOmegaMinus::AntiOmegaMinus()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- leAntiOmegaMinusModel_ = new G4LEAntiOmegaMinusInelastic();
- heAntiOmegaMinusModel_ = new G4HEAntiOmegaMinusInelastic();
- antiOmegaMinusInelastic_.RegisterMe(leAntiOmegaMinusModel_);
- antiOmegaMinusInelastic_.RegisterMe(heAntiOmegaMinusModel_);
- pManager->AddDiscreteProcess(&antiOmegaMinusInelastic_);
-
- pManager->AddProcess(&antiOmegaMinusIonisation_, ordInActive,2, 2);
-
- pManager->AddProcess(&antiOmegaMinusMult_);
- pManager->SetProcessOrdering(&antiOmegaMinusMult_, idxAlongStep, 1);
- pManager->SetProcessOrdering(&antiOmegaMinusMult_, idxPostStep, 1);
-
-}
diff --git a/src/G4BeamTestIonPhysics.cxx.backup b/src/G4BeamTestIonPhysics.cxx.backup
deleted file mode 100644
index a8cc573..0000000
--- a/src/G4BeamTestIonPhysics.cxx.backup
+++ /dev/null
@@ -1,100 +0,0 @@
-#include <iomanip>
-
-#include <globals.hh>
-#include <G4ios.hh>
-#include <G4ParticleDefinition.hh>
-#include <G4ParticleTable.hh>
-#include <G4IonConstructor.hh>
-#include <G4ProcessManager.hh>
-
-#include "G4BeamTestIonPhysics.h"
-
-G4BeamTestIonPhysics::G4BeamTestIonPhysics()
- : G4VPhysicsConstructor ("ion")
-{}
-
-/********************************************************************/
-
-G4BeamTestIonPhysics::~G4BeamTestIonPhysics()
-{}
-
-
-/********************************************************************/
-
-void G4BeamTestIonPhysics::ConstructParticle()
-{
- // Construct light ions
- G4IonConstructor pConstructor;
- pConstructor.ConstructParticle();
-}
-
-/********************************************************************/
-
-
-
-
-void G4BeamTestIonPhysics::ConstructProcess()
-{
- G4ProcessManager *pManager = 0;
-
- // Elastic Process
- elasticModel_ = new G4LElastic;
- elasticProcess_.RegisterMe(elasticModel_);
-
- // Generic Ion
- pManager = G4GenericIon::GenericIon()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- pManager->AddProcess(&ionMultipleScattering_, -1, 1, 1);
- pManager->AddProcess(&ionIonisation_, -1, 2, 2);
-
- // Deuteron
- pManager = G4Deuteron::Deuteron()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- deuteronModel_ = new G4LEDeuteronInelastic;
- deuteronProcess_.RegisterMe(deuteronModel_);
- pManager->AddDiscreteProcess(&deuteronProcess_);
-
- pManager->AddProcess(&deuteronMultipleScattering_, -1, 1, 1);
- pManager->AddProcess(&deuteronIonisation_, -1, 2, 2);
-
- // Triton
- pManager = G4Triton::Triton()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- tritonModel_ = new G4LETritonInelastic;
- tritonProcess_.RegisterMe(tritonModel_);
- pManager->AddDiscreteProcess(&tritonProcess_);
-
- pManager->AddProcess(&tritonMultipleScattering_, -1, 1, 1);
- pManager->AddProcess(&tritonIonisation_, -1, 2, 2);
-
- // Alpha
- pManager = G4Alpha::Alpha()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- alphaModel_ = new G4LEAlphaInelastic;
- alphaProcess_.RegisterMe(alphaModel_);
- pManager->AddDiscreteProcess(&alphaProcess_);
-
- pManager->AddProcess(&alphaMultipleScattering_, -1, 1, 1);
- pManager->AddProcess(&alphaIonisation_, -1, 2, 2);
-
- // He3
- pManager = G4He3::He3()->GetProcessManager();
- // add process
- pManager->AddDiscreteProcess(&elasticProcess_);
-
- pManager->AddProcess(&he3MultipleScattering_, -1, 1, 1);
- pManager->AddProcess(&he3Ionisation_, -1, 2, 2);
-
-
-}
-
-
-
diff --git a/src/G4BeamTestRunManager.cxx.backup b/src/G4BeamTestRunManager.cxx.backup
deleted file mode 100644
index 5ed2e05..0000000
--- a/src/G4BeamTestRunManager.cxx.backup
+++ /dev/null
@@ -1,131 +0,0 @@
-// On Sun, to prevent conflict with ObjectSpace, G4Timer.hh has to be
-// loaded *before* globals.hh...
-#include "G4Timer.hh"
-
-#include "G4BeamTestRunManager.h"
-#include <G4ParticleGun.hh>
-#include <G4Run.hh>
-
-
-G4BeamTestRunManager::G4BeamTestRunManager(): G4RunManager()
-{
-
-}
-
-
-void G4BeamTestRunManager::BeamOn(G4int n_event,const char* macroFile,G4int n_select)
-{
- G4String text = "BeamOn method is not supported in G4BeamTestRunManager!";
- G4Exception("G4BeamTestRunManager::BeamOn()", "G4BeamTestRunManager001", FatalException, text);
-}
-
-
-void G4BeamTestRunManager::InitializeRun()
-{
- G4bool cond = ConfirmBeamOnCondition();
- if(cond)
- {
- // Reset the event counter
- numberOfEventToBeProcessed = 0;
- ConstructScoringWorlds();
- RunInitialization();
-
- if(verboseLevel>0) timer->Start();
- }
-}
-
-
-void G4BeamTestRunManager::InjectParticle(G4ParticleGun* particleGun)
-{
- if(!currentRun)
- {
- G4String text = "Run needs to be initialized before injecting a particle.";
- G4Exception("G4BeamTestRunManager::InjectParticle()", "G4BeamTestRunManager002", FatalException, text);
- }
- assert(currentRun); // the G4Exception() above calls abort(). This assert() silences the clang static analyzer
-
- numberOfEventToBeProcessed++;
- currentRun->SetNumberOfEventToBeProcessed(numberOfEventToBeProcessed);
-
- currentEvent = GenerateEvent(numberOfEventToBeProcessed);
- particleGun->GeneratePrimaryVertex(currentEvent);
-
- eventManager->ProcessOneEvent(currentEvent);
- AnalyzeEvent(currentEvent);
- Update_Scoring();
- StackPreviousEvent(currentEvent);
- currentEvent = 0;
-
- if(runAborted) TerminateRun();
-}
-
-
-G4Event* G4BeamTestRunManager::GenerateEvent(G4int i_event)
-{
- G4Event* anEvent = new G4Event(i_event);
-
- if(storeRandomNumberStatusToG4Event==1 || storeRandomNumberStatusToG4Event==3)
- {
- std::ostringstream oss;
- CLHEP::HepRandom::saveFullState(oss);
- randomNumberStatusForThisEvent = oss.str();
- anEvent->SetRandomNumberStatus(randomNumberStatusForThisEvent);
- }
-
- if(storeRandomNumberStatus)
- {
- G4String fileN = randomNumberStatusDir + "currentEvent.rndm";
- CLHEP::HepRandom::saveEngineStatus(fileN);
- }
-
- return anEvent;
-}
-
-
-void G4BeamTestRunManager::TerminateRun()
-{
- if(verboseLevel>0)
- {
- timer->Stop();
- G4cout << "Run terminated." << G4endl;
- G4cout << "Run Summary" << G4endl;
- if(runAborted)
- {
- G4cout << " Run Aborted after " << numberOfEventToBeProcessed << " events processed." << G4endl;
- }
- else
- {
- G4cout << " Number of events processed : " << numberOfEventToBeProcessed << G4endl;
- }
- G4cout << " " << *timer << G4endl;
- }
-
- RunTermination();
-}
-
-
-//
-// The following method is an exact copy of
-// UpdateScoring which is private in the G4RunManager
-//
-
-#include <G4ScoringManager.hh>
-#include <G4HCofThisEvent.hh>
-#include <G4VHitsCollection.hh>
-
-void G4BeamTestRunManager::Update_Scoring()
-{
- G4ScoringManager* ScM = G4ScoringManager::GetScoringManagerIfExist();
- if(!ScM) return;
- G4int nPar = ScM->GetNumberOfMesh();
- if(nPar<1) return;
-
- G4HCofThisEvent* HCE = currentEvent->GetHCofThisEvent();
- if(!HCE) return;
- G4int nColl = HCE->GetCapacity();
- for(G4int i=0;i<nColl;i++)
- {
- G4VHitsCollection* HC = HCE->GetHC(i);
- if(HC) ScM->Accumulate(HC);
- }
-}
diff --git a/src/G4BeamTestSC4SD.cxx b/src/G4BeamTestSC4SD.cxx
new file mode 100644
index 0000000..f9ef766
--- /dev/null
+++ b/src/G4BeamTestSC4SD.cxx
@@ -0,0 +1,90 @@
+
+#include <sstream>
+
+#include "G4UnitsTable.hh"
+
+#include "G4BeamTestSC4SD.h"
+#include "G4HCofThisEvent.hh"
+#include "G4Step.hh"
+#include "G4ThreeVector.hh"
+#include "G4SDManager.hh"
+#include "G4ios.hh"
+
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+
+G4BeamTestSC4SD::G4BeamTestSC4SD(const G4String& name, const G4String& hitsCollectionName)
+ : G4VSensitiveDetector(name)
+ , fHitsCollection(NULL)
+{
+ collectionName.insert(hitsCollectionName);
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+
+G4BeamTestSC4SD::~G4BeamTestSC4SD()
+{
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+
+void G4BeamTestSC4SD::Initialize(G4HCofThisEvent* hce)
+{
+ // Create hits collection
+
+ fHitsCollection
+ = new G4BeamTestSiHitsCollection(SensitiveDetectorName, collectionName[0]);
+
+ // Add this collection in hce
+
+ G4int hcID
+ = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]);
+ hce->AddHitsCollection( hcID, fHitsCollection );
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+
+G4bool G4BeamTestSC4SD::ProcessHits(G4Step* aStep, G4TouchableHistory*)
+{
+ // G4double stepl = 10;
+ // aStep->SetStepLength(stepl);
+
+ G4String name=aStep->GetTrack()->GetDefinition()->GetParticleName();
+
+
+ if (name == "pi-" || name == "e-") {
+ G4cout << " Particle_name hit SC4 = " << name << G4endl;
+ // total energy
+ G4double etot = aStep->GetTrack()->GetTotalEnergy();
+ // energy deposit
+ G4double edep = aStep->GetTotalEnergyDeposit();
+
+ G4BeamTestSiHit* newHit = new G4BeamTestSiHit();
+ // G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
+
+ newHit->SetTrackID (aStep->GetTrack()->GetTrackID());
+ newHit->SetEdep(edep);
+ newHit->SetPos (aStep->GetPostStepPoint()->GetPosition());
+ newHit-> SetTime (aStep->GetPreStepPoint()->GetProperTime());
+ fHitsCollection->insert( newHit );
+
+ return true;
+ }
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+
+void G4BeamTestSC4SD::EndOfEvent(G4HCofThisEvent*)
+{
+ if ( verboseLevel>1 ) {
+ G4int nofHits = fHitsCollection->entries();
+ G4cout << G4endl
+ << "-------->Hits Collection: in this event they are " << nofHits
+ << " hits in the tracker chambers: " << G4endl;
+ for ( G4int i=0; i<nofHits; i++ ) (*fHitsCollection)[i]->Print();
+
+ }
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+
diff --git a/src/G4BeamTestSiSD.cxx b/src/G4BeamTestSiSD.cxx
index afa6bf8..194a0fc 100644
--- a/src/G4BeamTestSiSD.cxx
+++ b/src/G4BeamTestSiSD.cxx
@@ -15,8 +15,7 @@ std::ofstream testdata("EVENTDATA.txt");
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-G4BeamTestSiSD::G4BeamTestSiSD(const G4String& name,
- const G4String& hitsCollectionName)
+G4BeamTestSiSD::G4BeamTestSiSD(const G4String& name, const G4String& hitsCollectionName)
: G4VSensitiveDetector(name)
, fHitsCollection(NULL)
{
@@ -25,7 +24,7 @@ G4BeamTestSiSD::G4BeamTestSiSD(const G4String& name,
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-G4BeamTestSiSD::~G4BeamTestSiSD()
+G4BeamTestSiSD::~G4BeamTestSiSD()
{ // testdata << "Arrival time" << " " << "Energy " << " " << "Distance" << std::endl;
}
@@ -38,100 +37,92 @@ void G4BeamTestSiSD::Initialize(G4HCofThisEvent* hce)
testdata << "*" << std::endl;
- fHitsCollection
- = new G4BeamTestSiHitsCollection(SensitiveDetectorName, collectionName[0]);
+ fHitsCollection
+ = new G4BeamTestSiHitsCollection(SensitiveDetectorName, collectionName[0]);
// Add this collection in hce
- G4int hcID
+ G4int hcID
= G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]);
- G4cout << "hcID " << hcID << G4endl;
- hce->AddHitsCollection( hcID, fHitsCollection );
+ hce->AddHitsCollection( hcID, fHitsCollection );
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-G4bool G4BeamTestSiSD::ProcessHits(G4Step* aStep,
- G4TouchableHistory*)
-{
+G4bool G4BeamTestSiSD::ProcessHits(G4Step* aStep, G4TouchableHistory*)
+{
// G4double stepl = 10;
// aStep->SetStepLength(stepl);
G4String name=aStep->GetTrack()->GetDefinition()->GetParticleName();
- /* G4cout << " Particle_name = " << name << G4endl; */
+ // G4cout << " Particle_name = " << name << G4endl;
if (name == "opticalphoton" || name == "gamma") {
- // G4cout << " Particle_name = " << name << G4endl;
+ // G4cout << " Particle_name = " << name << G4endl;
-// total energy
- G4double etot = aStep->GetTrack()->GetTotalEnergy();
-// energy deposit
- G4double edep = aStep->GetTotalEnergyDeposit();
+ // total energy
+ G4double etot = aStep->GetTrack()->GetTotalEnergy();
+ // energy deposit
+ G4double edep = aStep->GetTotalEnergyDeposit();
- if (etot < 2.26 * CLHEP::eV) { // Lower threshold of PMT - 550nm
- // if (etot < 2.48 * CLHEP::eV) { // Lower threshold of PMT - 500nm
- // G4cout << "particle " << name << " under threshold with energy " << etot << G4endl;
- return true;
- }
- if (etot > 3.55 * CLHEP::eV) { // Upper threshold of PMT - 350nm
- // if (etot > 3.10 * CLHEP::eV) { // Upper threshold of PMT - 400nm
- // G4cout << "particle " << name << " over threshold with energy " << etot << G4endl;
- return true;
- }
- // G4cout << "inserting particle " << name << " with energy " << etot << " into record" << G4endl;
+ if (etot < 2.26 * CLHEP::eV) { // Lower threshold of PMT - 550nm
+ // if (etot < 2.48 * CLHEP::eV) { // Lower threshold of PMT - 500nm
+ // G4cout << "particle " << name << " under threshold with energy " << etot << G4endl;
+ return true;
+ }
+ if (etot > 3.55 * CLHEP::eV) { // Upper threshold of PMT - 350nm
+ // if (etot > 3.10 * CLHEP::eV) { // Upper threshold of PMT - 400nm
+ // G4cout << "particle " << name << " over threshold with energy " << etot << G4endl;
+ return true;
+ }
+ // G4cout << "inserting particle " << name << " with energy " << etot << " into record" << G4endl;
-// if (edep==0.) return false;
+ G4BeamTestSiHit* newHit = new G4BeamTestSiHit();
+ // G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
- /* G4cout << " Particle_name_after_edep = " << name << G4endl; */
+ newHit->SetTrackID (aStep->GetTrack()->GetTrackID());
+ newHit->SetEdep(edep);
+ newHit->SetPos (aStep->GetPostStepPoint()->GetPosition());
+ newHit-> SetTime (aStep->GetPreStepPoint()->GetProperTime());
+ fHitsCollection->insert( newHit );
- G4BeamTestSiHit* newHit = new G4BeamTestSiHit();
- // G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
- newHit->SetTrackID (aStep->GetTrack()->GetTrackID());
- newHit->SetEdep(edep);
- newHit->SetPos (aStep->GetPostStepPoint()->GetPosition());
- newHit-> SetTime (aStep->GetPreStepPoint()->GetProperTime());
- fHitsCollection->insert( newHit );
-
+ G4double zbef =aStep->GetPreStepPoint()->GetPosition().z();
+ G4double zaft =aStep->GetPostStepPoint()->GetPosition().z();
+ G4double zdelta = zaft-zbef;
+ //G4double zintial = aStep
- G4double zbef =aStep->GetPreStepPoint()->GetPosition().z();
- G4double zaft =aStep->GetPostStepPoint()->GetPosition().z();
- G4double zdelta = zaft-zbef;
- //G4double zintial = aStep
+ G4double zfix = (zdelta + (zbef- 4.7))*1000;
- G4double zfix = (zdelta + (zbef- 4.7))*1000;
+ // G4double zfix = (zdelta + (zbef- 375));
- // G4double zfix = (zdelta + (zbef- 375));
+ // G4cout << " Z VALUE = " << zfix << G4endl;
+ G4double time = aStep->GetPreStepPoint()->GetProperTime();
+ G4double gtime = aStep->GetPreStepPoint()->GetGlobalTime();
+ // G4cout << "Global Time " << gtime << G4endl;
+ G4double timefix = time - 0.288861775183071;
- // G4cout << " Z VALUE = " << zfix << G4endl;
- G4double time = aStep->GetPreStepPoint()->GetProperTime();
- G4double gtime = aStep->GetPreStepPoint()->GetGlobalTime();
- // G4cout << "Global Time " << gtime << G4endl;
- G4double timefix = time - 0.288861775183071;
+ /* G4cout << "TIME == " << time << G4endl; */
-/* G4cout << "TIME == " << time << G4endl; */
-
- G4double vel = aStep->GetPreStepPoint()->GetVelocity();
-
- G4double vel2 = zfix/time;
+ G4double vel = aStep->GetPreStepPoint()->GetVelocity();
- // G4cout << G4BestUnit(vel, "Speed")<< G4endl;
+ G4double vel2 = zfix/time;
-/*
+ // G4cout << G4BestUnit(vel, "Speed")<< G4endl;
- testdata << timefix << " " << edep << " " << zfix << std::endl;
-
- h_bragg9->Fill(zfix,edep);
- enetime->Fill(timefix, edep);
- // newHit->Print();
- sally->Fill(zfix,vel);
+ /*
- G4cout << "DISTANCE = " << zfix << " Velocity = " << vel2 << " Velocityold = " << vel <<G4endl;
+ testdata << timefix << " " << edep << " " << zfix << std::endl;
+ h_bragg9->Fill(zfix,edep);
+ enetime->Fill(timefix, edep);
+ // newHit->Print();
+ sally->Fill(zfix,vel);
-*/
- return true;
+ G4cout << "DISTANCE = " << zfix << " Velocity = " << vel2 << " Velocityold = " << vel <<G4endl;
+ */
+ return true;
}
}
@@ -139,10 +130,10 @@ G4bool G4BeamTestSiSD::ProcessHits(G4Step* aStep,
void G4BeamTestSiSD::EndOfEvent(G4HCofThisEvent*)
{
- if ( verboseLevel>1 ) {
+ if ( verboseLevel>1 ) {
G4int nofHits = fHitsCollection->entries();
G4cout << G4endl
- << "-------->Hits Collection: in this event they are " << nofHits
+ << "-------->Hits Collection: in this event they are " << nofHits
<< " hits in the tracker chambers: " << G4endl;
for ( G4int i=0; i<nofHits; i++ ) (*fHitsCollection)[i]->Print();
diff --git a/src/G4BeamTestTank.cxx b/src/G4BeamTestTank.cxx
index 02f6751..c4bc6b8 100644
--- a/src/G4BeamTestTank.cxx
+++ b/src/G4BeamTestTank.cxx
@@ -28,8 +28,8 @@
G4BeamTestTank::G4BeamTestTank()
{
// Get tank dimensions
- // tankThickness_ = 0.0*CLHEP::cm; // TODO(shivesh) : check thickness
- tankThickness_ = 0.44 * 2.54 *CLHEP::cm; // TODO(shivesh) : check thickness
+ // tankThickness_ = 0.0*CLHEP::cm;
+ tankThickness_ = 0.44 * 2.54 *CLHEP::cm;
tankHeight_ = 76.83 * 2.54 * CLHEP::cm;
innerRadius_ = 32 * 2.54 * CLHEP::cm;
outerRadius_ = innerRadius_ + tankThickness_;
@@ -55,112 +55,103 @@ G4VPhysicalVolume* G4BeamTestTank::InstallTank(G4VPhysicalVolume* mother, const
// User limits (energy cutoffs)
// Do not create photons or electrons below cherenkov threshold
// See also corresponding UserSpecialCuts in Physicslist !!!!
- // TODO(shivesh): Maybe do all of this as stepping action ??????
G4UserLimits* energyLimit = new G4UserLimits();
energyLimit->SetUserMinEkine(2.26 * CLHEP::eV); // Lower threshold of PMT - 550nm
- // std::string tankName=boost::lexical_cast<std::string>(tankKey_);
- std::string tankName = "BTT";
-
- // Define plastic frame TODO(shivesh): plastic or polyethelene?
G4Material* plastic = G4Material::GetMaterial("Plastic");
- G4Tubs* solidTank = new G4Tubs(("solid_tank_" + tankName).c_str(),
- 0.0 * CLHEP::m, outerRadius_, 0.5 * tankHeight_,
- 0.0 * CLHEP::deg, 360.0 * CLHEP::deg);
- tankLog_ = new G4LogicalVolume(solidTank, plastic,
- ("log_tank_" + tankName).c_str(), 0, 0, 0);
+ G4Tubs* solidTank = new G4Tubs(
+ "solid_tank_",
+ 0.0 * CLHEP::m, outerRadius_, 0.5 * tankHeight_,
+ 0.0 * CLHEP::deg, 360.0 * CLHEP::deg
+ );
+ tankLog_ = new G4LogicalVolume(
+ solidTank, plastic, "log_tank_", 0, 0, 0
+ );
// Define water volume
G4Material* water = G4Material::GetMaterial("Water");
- G4Tubs* solidWater = new G4Tubs(("solid_water_" + tankName).c_str(),
- 0.0 * CLHEP::m, innerRadius_, 0.5 * fillHeight_,
- 0.0 * CLHEP::deg, 360.0 * CLHEP::deg);
+ G4Tubs* solidWater = new G4Tubs(
+ "solid_water_",
+ 0.0 * CLHEP::m, innerRadius_, 0.5 * fillHeight_,
+ 0.0 * CLHEP::deg, 360.0 * CLHEP::deg
+ );
G4LogicalVolume* logWater =
- new G4LogicalVolume(solidWater, water, ("log_water_" + tankName).c_str(), 0, 0, 0);
- G4ThreeVector physWaterPosition(0, 0, -0.5*tankHeight_ + tankThickness_ + 0.5*fillHeight_);
+ new G4LogicalVolume(solidWater, water, "log_water_", 0, 0, 0);
+ G4ThreeVector physWaterPosition(
+ 0, 0, -0.5*tankHeight_ + tankThickness_ + 0.5*fillHeight_
+ );
G4VPhysicalVolume* physWater ATTRIBUTE_UNUSED =
- new G4PVPlacement(0, physWaterPosition, logWater,
- ("water_" + tankName).c_str(), tankLog_, false, 0);
+ new G4PVPlacement(
+ 0, physWaterPosition, logWater, "water_", tankLog_, false, 0
+ );
// Define air volume
G4Material* air = G4Material::GetMaterial("Air");
- G4cout << "airHeight_ = " << airHeight_ << G4endl;
- G4cout << "fillHeight_ = " << fillHeight_ << G4endl;
- G4Tubs* solidAir = new G4Tubs(("solid_air_" + tankName).c_str(),
- 0.0 * CLHEP::m, innerRadius_, 0.5 * airHeight_,
- 0.0 * CLHEP::deg, 360.0 * CLHEP::deg);
+ G4Tubs* solidAir = new G4Tubs(
+ "solid_air_",
+ 0.0 * CLHEP::m, innerRadius_, 0.5 * airHeight_,
+ 0.0 * CLHEP::deg, 360.0 * CLHEP::deg
+ );
G4LogicalVolume* logAir =
- new G4LogicalVolume(solidAir, air, ("log_air_" + tankName).c_str(), 0, 0, 0);
- G4ThreeVector physAirPosition(0, 0, -0.5 * tankHeight_ + 0.5 * CLHEP::cm + fillHeight_ +
- 0.5 * airHeight_);
- G4cout << "physAirPosition = " << physAirPosition << G4endl;
+ new G4LogicalVolume(solidAir, air, "log_air_", 0, 0, 0);
+ G4ThreeVector physAirPosition(0, 0,
+ -0.5 * tankHeight_ + 0.5 * CLHEP::cm + fillHeight_ + 0.5 * airHeight_
+ );
G4VPhysicalVolume* physAir_UNUSED =
- new G4PVPlacement(0, physAirPosition, logAir,
- ("air_" + tankName).c_str(), tankLog_, false, 0);
+ new G4PVPlacement(
+ 0, physAirPosition, logAir, "air_", tankLog_, false, 0
+ );
// Define glass sphere & effective DOM material splitted in upper and lower part
G4Material* glass = G4Material::GetMaterial("Glass");
G4Material* effectiveDOM = G4Material::GetMaterial("effectiveDOM");
- // std::map<OMKey, G4ThreeVector> domPosIce;
- // std::map<OMKey, G4ThreeVector>::const_iterator om_iter;
- // for(om_iter=relDomPositions_.begin(); om_iter!=relDomPositions_.end(); ++om_iter)
- // {
- /* const OMKey& omKey = om_iter->first; */
-
G4ThreeVector upperDOMpos(0, 0, -0.5 * airHeight_);
G4ThreeVector lowerDOMpos(0, 0, 0.5 * fillHeight_);
- G4cout << "upperDOMpos = " << upperDOMpos << G4endl;
- G4cout << "lowerDOMpos = " << lowerDOMpos << G4endl;
-
- // domPosIce[omKey] = lowerDOMpos;
-
- // std::string omName=boost::lexical_cast<std::string>(omKey);
- std::string omName="BTOM";
- G4Sphere *upperglasssphere = new G4Sphere (("solid_dom_up_" + omName).c_str(),
+ G4Sphere *upperglasssphere = new G4Sphere ("solid_dom_up_",
0.0 * CLHEP::m, glassOuterRadius_,
0.0 * CLHEP::deg, 360.0 * CLHEP::deg,
0.0 * CLHEP::deg, 90.0 * CLHEP::deg);
- G4Sphere *lowerglasssphere = new G4Sphere (("solid_dom_lo_" + omName).c_str(),
+ G4Sphere *lowerglasssphere = new G4Sphere ("solid_dom_lo_",
0.0 * CLHEP::m, glassOuterRadius_,
0.0 * CLHEP::deg, 360.0 * CLHEP::deg,
90.0 * CLHEP::deg, 180.0 * CLHEP::deg);
G4double domInnerRadius = glassOuterRadius_ - glassThickness_;
- G4Sphere *upperdomsphere = new G4Sphere (("solid_inside_dom_up_" + omName).c_str(),
+ G4Sphere *upperdomsphere = new G4Sphere ("solid_inside_dom_up_",
0.0 * CLHEP::m, domInnerRadius,
0.0 * CLHEP::deg, 360.0 * CLHEP::deg,
0.0 * CLHEP::deg, 90.0 * CLHEP::deg);
- G4Sphere *lowerdomsphere = new G4Sphere (("solid_inside_dom_lo_" + omName).c_str(),
+ G4Sphere *lowerdomsphere = new G4Sphere ("solid_inside_dom_lo_",
0.0 * CLHEP::m, domInnerRadius,
0.0 * CLHEP::deg, 360.0 * CLHEP::deg,
90.0 * CLHEP::deg, 180.0 * CLHEP::deg);
G4LogicalVolume* logUpperGlass =
new G4LogicalVolume(upperglasssphere, glass,
- ("log_dom_up_" + omName).c_str(), 0, 0, 0);
+ "log_dom_up_", 0, 0, 0);
G4LogicalVolume* logLowerGlass =
new G4LogicalVolume(lowerglasssphere, glass,
- ("log_dom_lo_" + omName).c_str(), 0, 0, 0);
+ "log_dom_lo_", 0, 0, 0);
G4LogicalVolume* logUpperDOM =
new G4LogicalVolume(upperdomsphere, effectiveDOM,
- ("log_inside_dom_up_" + omName).c_str(), 0, 0, 0);
+ "log_inside_dom_up_", 0, 0, 0);
G4LogicalVolume* logLowerDOM =
new G4LogicalVolume(lowerdomsphere, effectiveDOM,
- ("log_inside_dom_lo_" + omName).c_str(), 0, 0, 0);
+ "log_inside_dom_lo_", 0, 0, 0);
G4VPhysicalVolume* physUpperGlass ATTRIBUTE_UNUSED =
new G4PVPlacement(0, upperDOMpos, logUpperGlass,
- ("dom_up_" + omName).c_str(), logAir, false, 0);
+ "dom_up_", logAir, false, 0);
G4VPhysicalVolume* physLowerGlass ATTRIBUTE_UNUSED =
new G4PVPlacement(0, lowerDOMpos, logLowerGlass,
- ("dom_lo_" + omName).c_str(), logWater, false, 0);
+ "dom_lo_", logWater, false, 0);
G4VPhysicalVolume* physUpperDOM ATTRIBUTE_UNUSED =
new G4PVPlacement(0, G4ThreeVector(0,0,0), logUpperDOM,
- ("inside_dom_up_" + omName).c_str(), logUpperGlass, false, 0);
+ "inside_dom_up_", logUpperGlass, false, 0);
G4VPhysicalVolume* physLowerDOM ATTRIBUTE_UNUSED =
new G4PVPlacement(0, G4ThreeVector(0,0,0), logLowerDOM,
- ("inside_dom_lo_" + omName).c_str(), logLowerGlass, false, 0);
+ "inside_dom_lo_", logLowerGlass, false, 0);
// apply energy limits
logUpperGlass->SetUserLimits(energyLimit);
@@ -171,24 +162,21 @@ G4VPhysicalVolume* G4BeamTestTank::InstallTank(G4VPhysicalVolume* mother, const
// Define sensitive detector
G4SDManager* sdManager = G4SDManager::GetSDMpointer();
- iceSD_ = new G4BeamTestSiSD(("ice_SD_" + tankName).c_str(), "HitsCollection");
+ iceSD_ = new G4BeamTestSiSD("ice_SD_", "HitsCollection");
sdManager->AddNewDetector(iceSD_);
// logLowerDOM->SetSensitiveDetector(iceSD_);
logLowerGlass->SetSensitiveDetector(iceSD_);
// Instantiation of a set of visualization attributes with red colour
G4VisAttributes * tankVisAtt = new G4VisAttributes(G4Colour(1,0,0));
- // Set the forced wireframe style
- //snowVisAtt->SetForceWireFrame(true);
// Assignment of the visualization attributes to the logical volume
tankLog_->SetVisAttributes(tankVisAtt);
- G4ThreeVector tankPos = position_ - origin - mother->GetTranslation();
- G4cout << "tankPos = " << tankPos << G4endl;
-
- G4VPhysicalVolume* tankPhys = new G4PVPlacement(0, tankPos, tankLog_,
- ("tank_" + tankName).c_str(),
- mother->GetLogicalVolume(), false, 0);
+ // G4ThreeVector tankPos = position_ - origin - mother->GetTranslation();
+ G4ThreeVector tankPos = origin;
+ G4VPhysicalVolume* tankPhys = new G4PVPlacement(
+ 0, tankPos, tankLog_, "tank_", mother->GetLogicalVolume(), false, 0
+ );
// apply energy limits
tankLog_->SetUserLimits(energyLimit);
@@ -197,84 +185,3 @@ G4VPhysicalVolume* G4BeamTestTank::InstallTank(G4VPhysicalVolume* mother, const
return tankPhys;
}
-
-
-/* double G4BeamTestTank::GetNumCherenkov(#<{(| const OMKey& omKey |)}>#) */
-/* { */
-/* return std::max(iceSD_->GetNumCherenkov(#<{(| omKey |)}>#), 0.); */
-/* } */
-/* */
-/* */
-/* double G4BeamTestTank::GetNumCherenkovWeight(#<{(| const OMKey& omKey |)}>#) */
-/* { */
-/* return std::max(iceSD_->GetNumCherenkovWeight(#<{(| omKey |)}>#), 0.); */
-/* } */
-/* */
-/* */
-/* double G4BeamTestTank::GetEDep_G4(#<{(| const OMKey& omKey |)}>#) */
-/* { */
-/* return std::max(iceSD_->GetEDep(#<{(| omKey |)}>#), 0.); */
-/* } */
-/* */
-/* */
-/* double G4BeamTestTank::GetTime_G4(#<{(| const OMKey& omKey |)}>#) */
-/* { */
-/* return iceSD_->GetTime(#<{(| omKey |)}>#); */
-/* } */
-/* */
-
-/* double G4BeamTestTank::GetEDep_I3(const OMKey& omKey) */
-/* { */
-/* return std::max(iceSD_->GetEDep(omKey), 0.) / CLHEP::keV * I3Units::keV; */
-/* } */
-
-
-/* double G4BeamTestTank::GetTime_I3(const OMKey& omKey) */
-/* { */
-/* return (iceSD_->GetTime(omKey) / CLHEP::s) * I3Units::s; */
-/* } */
-
-
-/* I3Position G4BeamTestTank::GetPos_I3() */
-/* { */
-/* I3Position pos((position_.x() / CLHEP::m) * I3Units::m, */
-/* (position_.y() / CLHEP::m) * I3Units::m, */
-/* (position_.z() / CLHEP::m) * I3Units::m); */
-/* return pos; */
-/* } */
-
-
-/* double G4BeamTestTank::GetX_I3() */
-/* { */
-/* return (position_.x() / CLHEP::m) * I3Units::m; */
-/* } */
-/* */
-/* */
-/* double G4BeamTestTank::GetY_I3() */
-/* { */
-/* return (position_.y() / CLHEP::m) * I3Units::m; */
-/* } */
-/* */
-/* */
-/* double G4BeamTestTank::GetZ_I3() */
-/* { */
-/* return (position_.z() / CLHEP::m) * I3Units::m; */
-/* } */
-/* */
-
-/* double G4BeamTestTank::GetTankHeight_I3() */
-/* { */
-/* return (tankHeight_ / CLHEP::m) */
-/* } */
-/* */
-/* */
-/* double G4BeamTestTank::GetTankRadius_I3() */
-/* { */
-/* return (outerRadius_ / CLHEP::m) */
-/* } */
-/* */
-/* */
-/* double G4BeamTestTank::GetFillHeight_I3() */
-/* { */
-/* return (fillHeight_ / CLHEP::m) */
-/* } */
diff --git a/src/G4BeamTestUserSteppingAction.cxx b/src/G4BeamTestUserSteppingAction.cxx
index f81ab62..816dc4a 100644
--- a/src/G4BeamTestUserSteppingAction.cxx
+++ b/src/G4BeamTestUserSteppingAction.cxx
@@ -40,7 +40,7 @@ void G4BeamTestUserSteppingAction::UserSteppingAction(const G4Step* step)
//check if particle energy is below threshold; if true, kill the particle
G4double energy = track->GetTotalEnergy();
if(energy < threshold){
- G4cout << "SteppingAction: killing particle " << particle << " with energy " << energy << " < " << threshold << G4endl;
+ if (energy > 0) G4cout << "SteppingAction: killing particle " << particle << " with energy " << energy << " < " << threshold << G4endl;
track->SetTrackStatus(fStopAndKill);
}
}
diff --git a/src/G4BeamTestUserTrackingAction.cxx b/src/G4BeamTestUserTrackingAction.cxx
index e41ccc3..09a7e59 100644
--- a/src/G4BeamTestUserTrackingAction.cxx
+++ b/src/G4BeamTestUserTrackingAction.cxx
@@ -11,36 +11,35 @@ void G4BeamTestUserTrackingAction::PreUserTrackingAction(const G4Track*){}
void G4BeamTestUserTrackingAction::PostUserTrackingAction(const G4Track* track)
{
- const G4LogicalVolume *volume = track->GetLogicalVolumeAtVertex();
- G4UserLimits *limit = volume->GetUserLimits();
- if(!limit) G4cout << "----> G4LogicalVolume: " << volume->GetName() << " has no defined G4UserLimit" << G4endl;
- G4double threshold = limit->GetUserMinEkine(*track);
- G4double max_threshold = 3.54;
- G4TrackVector* secondaries = fpTrackingManager->GimmeSecondaries();
- if(secondaries)
- {
- size_t nSeco = secondaries->size();
- if(nSeco>0)
- {
- for(size_t i=0;i<nSeco;i++)
- {
- //check if secondary particle is a gamma
- G4String particle = (*secondaries)[i]->GetDefinition()->GetParticleName();
- if(particle == "gamma" || particle == "opticalphoton")
- {
- //check if particle energy is below threshold; if true, kill the particle
- G4double energy = (*secondaries)[i]->GetTotalEnergy();
- // if(energy < threshold){
- // G4cout << "TrackingAction: killing particle " << particle << " with energy " << energy << " < " << threshold << G4endl;
- // (*secondaries)[i]->SetTrackStatus(fStopAndKill);
- // }
+ // const G4LogicalVolume *volume = track->GetLogicalVolumeAtVertex();
+ // G4UserLimits *limit = volume->GetUserLimits();
+ // if(!limit) G4cout << "----> G4LogicalVolume: " << volume->GetName() << " has no defined G4UserLimit" << G4endl;
+ // G4double threshold = limit->GetUserMinEkine(*track);
+ // G4double max_threshold = 3.54;
+ // G4TrackVector* secondaries = fpTrackingManager->GimmeSecondaries();
+ // if(secondaries)
+ // {
+ // size_t nSeco = secondaries->size();
+ // if(nSeco>0)
+ // {
+ // for(size_t i=0;i<nSeco;i++)
+ // {
+ // //check if secondary particle is a gamma
+ // G4String particle = (*secondaries)[i]->GetDefinition()->GetParticleName();
+ // if(particle == "gamma" || particle == "opticalphoton")
+ // {
+ // //check if particle energy is below threshold; if true, kill the particle
+ // G4double energy = (*secondaries)[i]->GetTotalEnergy();
+ // if(energy < threshold){
+ // G4cout << "TrackingAction: killing particle " << particle << " with energy " << energy << " < " << threshold << G4endl;
+ // (*secondaries)[i]->SetTrackStatus(fStopAndKill);
+ // }
// if (energy > max_threshold * CLHEP::eV){
// G4cout << "TrackingAction: killing particle " << particle << " with energy " << energy << " > " << max_threshold << G4endl;
- // (*secondaries)[i]->SetTrackStatus(fStopAndKill);
-
+ // (*secondaries)[i]->SetTrackStatus(fStopAndKill);
// }
- }
- }
- }
- }
+ // }
+ // }
+ // }
+ // }
}
diff --git a/src/G4TankIceSD.cxx b/src/G4TankIceSD.cxx
deleted file mode 100644
index 2974837..0000000
--- a/src/G4TankIceSD.cxx
+++ /dev/null
@@ -1,166 +0,0 @@
-#include "G4TankIceSD.h"
-#include "G4BeamTestTank.h"
-
-#include <G4Step.hh>
-#include <G4HCofThisEvent.hh>
-#include <G4TouchableHistory.hh>
-#include <G4ios.hh>
-#include <iostream>
-#include <ios>
-#include <fstream>
-
-#include <G4VProcess.hh>
-#include <G4OpticalPhoton.hh>
-#include "G4Poisson.hh"
-
-/* #include <icetray/I3Units.h> */
-
-#include "G4Version.hh"
-
-#if G4VERSION_NUMBER >= 950
-// The G4MaterialPropertyVector is gone since 4.9.5.
-// It has been typedef'd to a G4UnorderedPhysicsVector
-// with a different interface. Try to support both
-// versions with an ifdef.
-#define MATERIAL_PROPERTY_VECTOR_IS_PHYSICS_VECTOR
-#endif
-
-
-G4TankIceSD::G4TankIceSD(G4String name/* , const std::map<OMKey, G4ThreeVector>& domPositions */)
- : G4VSensitiveDetector(name)/* , domPositions_(domPositions) */
-{}
-
-
-G4TankIceSD::~G4TankIceSD() {}
-
-
-void G4TankIceSD::Initialize(G4HCofThisEvent* HCE)
-{
- // std::map<OMKey, G4ThreeVector>::const_iterator om_iter;
- // for(om_iter=domPositions_.begin(); om_iter!=domPositions_.end(); ++om_iter)
- // {
- sumEdep_ = 0;
- cogTime_ = 0;
- cherenkovCounter_ = 0;
- cherenkovCounterWeight_ = 0;
- // }
-}
-
-
-G4bool G4TankIceSD::ProcessHits(G4Step* aStep, G4TouchableHistory*)
-{
- // Get energy deposit/time --> CHECK WHAT THESE QUANTITIES ACTUALLY MEAN !!!!!!!!!!!!!!!!!!!!!!!!
- G4double edep = aStep->GetTotalEnergyDeposit();
- G4double time = aStep->GetPreStepPoint()->GetGlobalTime();
- G4double cherenkovNumber = GetCerenkovNumber(aStep);
-
- if(edep<=0 && cherenkovNumber<=0) return false;
-
- // Get Position relative to ice center
- G4ThreeVector preStepPoint = aStep->GetPreStepPoint()->GetPosition();
- G4ThreeVector postStepPoint = aStep->GetPostStepPoint()->GetPosition();
- G4TouchableHandle theTouchable = aStep->GetPreStepPoint()->GetTouchableHandle();
- G4ThreeVector worldPosition = (preStepPoint + postStepPoint) / 2.0;
- G4ThreeVector localPosition = theTouchable->GetHistory()->
- GetTopTransform().TransformPoint(worldPosition);
-
- // std::map<OMKey, G4ThreeVector>::const_iterator om_iter;
- // for(om_iter=domPositions_.begin(); om_iter!=domPositions_.end(); ++om_iter)
- // {
- G4double radius = sqrt(pow(domPositions_.x() - localPosition.x(), 2) +
- pow(domPositions_.y() - localPosition.y(), 2));
- G4double height = (domPositions_.z() - localPosition.z());
-
- sumEdep_ += edep;
- cogTime_ += edep*time;
- cherenkovCounterWeight_ += GetProbability(radius, height) * cherenkovNumber;
- cherenkovCounter_ += cherenkovNumber;
- // }
- return true;
-}
-
-
-void G4TankIceSD::EndOfEvent(G4HCofThisEvent*)
-{
- // std::map<OMKey, G4ThreeVector>::const_iterator om_iter;
- // for(om_iter=domPositions_.begin(); om_iter!=domPositions_.end(); ++om_iter)
- // {
- if(sumEdep_>0)
- {
- cogTime_ /= sumEdep_;
- }
- // }
-}
-
-
-G4double G4TankIceSD::GetCerenkovNumber(G4Step* aStep)
-{
- // same function as in "source/processes/electromagnetic/xrays/src/G4Cerenkov.cc" but without
- // cerenkov angle and only for an energy independent refraction index
- G4Track* aTrack = aStep->GetTrack();
- const G4DynamicParticle* aParticle = aTrack->GetDynamicParticle();
- const G4Material* aMaterial = aTrack->GetMaterial();
-
- G4StepPoint* pPreStepPoint = aStep->GetPreStepPoint();
- G4StepPoint* pPostStepPoint = aStep->GetPostStepPoint();
-
- G4MaterialPropertiesTable* aMaterialPropertiesTable = aMaterial->GetMaterialPropertiesTable();
- if (!aMaterialPropertiesTable) return 0.0;
-
- G4MaterialPropertyVector* Rindex = aMaterialPropertiesTable->GetProperty("RINDEX");
- if (!Rindex) return 0.0;
-
- const G4double Rfact = 369.81 / (CLHEP::eV * CLHEP::cm); // TODO(shivesh): check this
- const G4double charge = aParticle->GetDefinition()->GetPDGCharge();
- const G4double beta = (pPreStepPoint->GetBeta() + pPostStepPoint->GetBeta()) / 2.;
- if (beta <= 0.0) return 0.0;
- G4double BetaInverse = 1. / beta;
-
- // Min and Max photon energies
-#ifdef MATERIAL_PROPERTY_VECTOR_IS_PHYSICS_VECTOR
- G4double Pmin = Rindex->GetMinLowEdgeEnergy();
- G4double Pmax = Rindex->GetMaxLowEdgeEnergy();
-#else
- G4double Pmin = Rindex->GetMinPhotonEnergy();
- G4double Pmax = Rindex->GetMaxPhotonEnergy();
-#endif
-
- // Min and Max Refraction Indices
-#ifdef MATERIAL_PROPERTY_VECTOR_IS_PHYSICS_VECTOR
- G4double nMin = Rindex->GetMinValue();
- G4double nMax = Rindex->GetMaxValue();
-#else
- G4double nMin = Rindex->GetMinProperty();
- G4double nMax = Rindex->GetMaxProperty();
-#endif
- if (nMin!=nMax)
- {
- /* log_error("Support for energy dependent refraction index not yet implemented!"); */
- G4cout << "Support for energy dependent refraction index not yet implemented!";
- return 0.0;
- }
-
- // If n(Pmax) < 1/Beta -- no photons generated
- if (nMax < BetaInverse) return 0.0;
-
- G4double MeanNumberOfPhotons = Rfact * charge/CLHEP::eplus * charge/CLHEP::eplus * (Pmax - Pmin) *
- (1. - BetaInverse * BetaInverse / nMin / nMin);
-
- G4double step_length = aStep->GetStepLength();
-
- return MeanNumberOfPhotons * step_length;
-}
-
-
-// TODO(shivesh): what is this
-G4double G4TankIceSD::GetProbability(G4double radius, G4double height)
-{
- height = 0.90 - height / CLHEP::m;
- radius = radius / CLHEP::m;
- G4double p0 = 0.0340 + 0.0120 * height;
- G4double p1 = 0.0400 + 0.000622 * exp(-1.13 + 8.03 * height);
- G4double p2 = 0.917 + 2.34 * exp(-1.82 + 3.12 * height);
- G4double p3 = -0.00325 - 0.00199 * height - 0.0121 * height*height;
-
- return (p0 + p3 * radius) + (p1 - p0) * exp(-p2*p2 * radius*radius);
-}
generated by cgit v1.2.3 (git 2.25.1) at 2025-12-16 13:55:32 +0000