Tue 21 Aug 15:36:34 BST 2018

13 files changed, 328 insertions, 681 deletions

diff --git a/src/G4BeamTestDetectorConstruction.cxx b/src/G4BeamTestDetectorConstruction.cxx
index 6eeaa54..cf75551 100644
--- a/src/G4BeamTestDetectorConstruction.cxx
+++ b/src/G4BeamTestDetectorConstruction.cxx
@@ -30,71 +30,31 @@ G4VPhysicalVolume* G4BeamTestDetectorConstruction::Construct()
 
   CreateMaterials();
 
-  // Determine bottom z position of snow volume
-  G4double zSnowBottom(tankList_.at(0)->GetPos().z());
-  BOOST_FOREACH(G4BeamTestTank* tank, tankList_)
-  {
-    // z position of bottom of tank
-    G4double z = tank->GetPos().z() - 0.5*tank->GetTankHeight_G4();
-    zSnowBottom = std::min(z, zSnowBottom);
-  }
-
-  // Subtract safety margin
-  zSnowBottom -= 1.0*CLHEP::m;
-
-  // Triangulate snow surface
-  G4Delaunay delaunay;
-  BOOST_FOREACH(G4BeamTestTank* tank, tankList_)
-  {
-    // z position of snow surface
-    G4double z = tank->GetPos().z() + 0.5 * tank->GetTankHeight_G4() + tank->GetSnowHeight_G4();
-
-    delaunay.AddPoint(tank->GetDelaunayPoint1().x(),
-                      tank->GetDelaunayPoint1().y(),
-                      z - zSnowBottom);
-    delaunay.AddPoint(tank->GetDelaunayPoint2().x(),
-                      tank->GetDelaunayPoint2().y(),
-                      z - zSnowBottom);
-    delaunay.AddPoint(tank->GetDelaunayPoint3().x(),
-                      tank->GetDelaunayPoint3().y(),
-                      z - zSnowBottom);
-  }
-
-  // Create tesselated snow volume
-  G4TessellatedSolid* solidSnow = new G4TessellatedSolid("solid_snow");
-  delaunay.BuildSolid(solidSnow, 50.0*CLHEP::m, 100.0*CLHEP::m);
+  /* // World origin in IceCube coordinates */
+  /* origin_.set(delaunay.GetOrigin().x(), delaunay.GetOrigin().y(), zSnowBottom + zHalfLength); */
 
   // Determine World dimensions
-  G4double xHalfLength = 0.5 * (delaunay.GetXmax() - delaunay.GetXmin());
-  G4double yHalfLength = 0.5 * (delaunay.GetYmax() - delaunay.GetYmin());
-  G4double zHalfLength = 0.5 * (delaunay.GetZmax() + 20.0 * CLHEP::m); // 20 m of atmosphere
-
-  // World origin in IceCube coordinates
-  origin_.set(delaunay.GetOrigin().x(), delaunay.GetOrigin().y(), zSnowBottom + zHalfLength);
+  G4double xWorld = 20.0 * CLHEP::m;
+  G4double yWorld = 20.0 * CLHEP::m;
+  G4double zWorld = 20.0 * CLHEP::m;
 
   // Create world volume
-  G4Box* world_box = new G4Box("solid_world", xHalfLength, yHalfLength, zHalfLength);
+  G4Box* world_box = new G4Box("solid_world", xWorld, yWorld, zWorld);
   G4LogicalVolume* worldLog =
-      new G4LogicalVolume(world_box, G4Material::GetMaterial("Air"), "log_world", 0, 0, 0);
+      new G4LogicalVolume(world_box, G4Material::GetMaterial("G4_AIR"), "log_world", 0, 0, 0);
   G4VPhysicalVolume* worldPhys =
       new G4PVPlacement(0, G4ThreeVector(), worldLog, "world", 0, false, 0);
 
-  // Snow layer
-  G4LogicalVolume* snowLog =
-      new G4LogicalVolume(solidSnow, G4Material::GetMaterial("Snow"), "log_snow", 0, 0, 0);
-  G4VPhysicalVolume* snowPhys =
-      new G4PVPlacement(0, G4ThreeVector(0, 0, -zHalfLength), snowLog, "snow", worldLog, 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);
+  /* // 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 tanks
-  BOOST_FOREACH(G4BeamTestTank* tank, tankList_)
-  {
-    tank->InstallTank(snowPhys, origin_);
-  }
+  // Install tank
+  /* BOOST_FOREACH(G4BeamTestTank* tank, tankList_) */
+  /* { */
+  tank_->InstallTank(worldPhys, origin_);
+  // }
 
   // User limits (energy cutoffs)
   // Do not create photons or electrons below cherenkov threshold
@@ -102,7 +62,7 @@ G4VPhysicalVolume* G4BeamTestDetectorConstruction::Construct()
   G4UserLimits* energyLimit = new G4UserLimits();
   energyLimit->SetUserMinEkine(280.0 * CLHEP::keV);  // Cherenkov threshold of electrons in ice TODO(shivesh)
   worldLog->SetUserLimits(energyLimit);
-  snowLog->SetUserLimits(energyLimit);
+  /* snowLog->SetUserLimits(energyLimit); */
   
   return worldPhys;
 }
@@ -128,7 +88,7 @@ void G4BeamTestDetectorConstruction::CreateMaterials()
 void G4BeamTestDetectorConstruction::CreateAir()
 {
   G4NistManager* nistManager = G4NistManager::Instance();
-  nistManager->ConstructNewGasMaterial("Air","G4_AIR"/* , (273.15 - 40.0)*CLHEP::kelvin, 670.0E-3*CLHEP::bar */);
+  nistManager->ConstructNewGasMaterial("Air","G4_AIR" , (273.15 + 20.0)*CLHEP::kelvin, 1.0*CLHEP::atmosphere);
 }
 
 /*****************************************************************/
@@ -136,10 +96,19 @@ void G4BeamTestDetectorConstruction::CreateAir()
 void G4BeamTestDetectorConstruction::CreateWater()
 {
   G4NistManager* nistManager = G4NistManager::Instance();
-  // G4Material* ice = new G4Material("Water", 1.0 * CLHEP::g / CLHEP::cm3, 2, kStateLiquid);
-  // ice->AddElement(nistManager->FindOrBuildElement("H"), 2);
-  // ice->AddElement(nistManager->FindOrBuildElement("O"), 1);
-  nistManager->ConstructNewGasMaterial("Water","G4_WATER");
+  G4Material* water = new G4Material("Water", 1.0 * CLHEP::g / CLHEP::cm3, 2, kStateLiquid);
+  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_refr[water_bins] = {1.33, 1.33};
+  G4MaterialPropertiesTable *mpt_water = new G4MaterialPropertiesTable ();
+  mpt_water->AddProperty ("RINDEX", water_ephot, water_refr, water_bins);
+  water->SetMaterialPropertiesTable(mpt_water);
 }
 
 /*****************************************************************/
diff --git a/src/G4BeamTestEMPhysics.cxx b/src/G4BeamTestEMPhysics.cxx
index d36207a..55ecf99 100644
--- a/src/G4BeamTestEMPhysics.cxx
+++ b/src/G4BeamTestEMPhysics.cxx
@@ -1,16 +1,3 @@
-/*
- * copyright  (C) 2010
- * The Icecube Collaboration
- *
- * $Id: G4BeamTestEMPhysics.cxx 86420 2012-03-20 16:00:37Z kislat $
- *
- * @version $Revision: 86420 $
- * @date $LastChangedDate: 2012-03-20 16:00:37 +0000 (Tue, 20 Mar 2012) $
- * @author Fabian Kislat <fabian.kislat@desy.de>, Peter Nießen (tanktop)
- * Last changed by: $LastChangedBy: kislat $
- */
-
-
 #include <globals.hh>
 #include <G4ParticleDefinition.hh>
 #include <G4ParticleTable.hh>
diff --git a/src/G4BeamTestGeneralPhysics.cxx b/src/G4BeamTestGeneralPhysics.cxx
index 2def9b4..93304a8 100644
--- a/src/G4BeamTestGeneralPhysics.cxx
+++ b/src/G4BeamTestGeneralPhysics.cxx
@@ -1,14 +1,3 @@
-/*
- * copyright  (C) 2010
- * The Icecube Collaboration
- *
- * $Id: G4BeamTestGeneralPhysics.cxx 152849 2017-01-20 21:44:25Z jgonzalez $
- *
- * @version $Revision: 152849 $
- * @date $LastChangedDate: 2017-01-20 21:44:25 +0000 (Fri, 20 Jan 2017) $
- * @author Fabian Kislat <fabian.kislat@desy.de> Last changed by: $LastChangedBy: jgonzalez $
- */
-
 #include "G4BeamTestGeneralPhysics.h"
 
 #include <iomanip>   
diff --git a/src/G4BeamTestHadronPhysics.cxx b/src/G4BeamTestHadronPhysics.cxx.backup
index 5150b67..1e37ed5 100644
--- a/src/G4BeamTestHadronPhysics.cxx
+++ b/src/G4BeamTestHadronPhysics.cxx.backup
@@ -1,16 +1,3 @@
-/*
- * copyright  (C) 2010
- * The Icecube Collaboration
- *
- * $Id: G4BeamTestHadronPhysics.cxx 152814 2017-01-19 21:34:52Z jgonzalez $
- *
- * @version $Revision: 152814 $
- * @date $LastChangedDate: 2017-01-19 21:34:52 +0000 (Thu, 19 Jan 2017) $
- * @author Fabian Kislat <fabian.kislat@desy.de> Last changed by: $LastChangedBy: jgonzalez $
- *
- * Copied from tanktop
- */
-
 #include <iomanip>   
 
 #include <globals.hh>
diff --git a/src/G4BeamTestIonPhysics.cxx b/src/G4BeamTestIonPhysics.cxx.backup
index 3c7e553..a8cc573 100644
--- a/src/G4BeamTestIonPhysics.cxx
+++ b/src/G4BeamTestIonPhysics.cxx.backup
@@ -1,15 +1,3 @@
-/*
- * copyright  (C) 2010
- * The Icecube Collaboration
- *
- * $Id: G4BeamTestIonPhysics.cxx 86420 2012-03-20 16:00:37Z kislat $
- *
- * @version $Revision: 86420 $
- * @date $LastChangedDate: 2012-03-20 16:00:37 +0000 (Tue, 20 Mar 2012) $
- * @author Fabian Kislat <fabian.kislat@desy.de>, Peter Nießen (tanktop)
- * Last changed by: $LastChangedBy: kislat $
- */
-
 #include <iomanip>   
 
 #include <globals.hh>
diff --git a/src/G4BeamTestMuonPhysics.cxx b/src/G4BeamTestMuonPhysics.cxx
index 92264a2..2c5404e 100644
--- a/src/G4BeamTestMuonPhysics.cxx
+++ b/src/G4BeamTestMuonPhysics.cxx
@@ -1,16 +1,3 @@
-/*
- * copyright  (C) 2010
- * The Icecube Collaboration
- *
- * $Id: G4BeamTestMuonPhysics.cxx 86420 2012-03-20 16:00:37Z kislat $
- *
- * @version $Revision: 86420 $
- * @date $LastChangedDate: 2012-03-20 16:00:37 +0000 (Tue, 20 Mar 2012) $
- * @author Fabian Kislat <fabian.kislat@desy.de>, Peter Nießen (tanktop)
- * Last changed by: $LastChangedBy: kislat $
- */
-
-
 #include <globals.hh>
 #include <G4ParticleDefinition.hh>
 #include <G4ParticleTable.hh>
@@ -78,7 +65,7 @@ void G4BeamTestMuonPhysics::ConstructProcess()
   pManager->AddProcess(&muMinusBremsstrahlung_,     -1,  3, 3);
   pManager->AddProcess(&muMinusPairProduction_,     -1,  4, 4);
 
-  pManager->AddRestProcess(&muMinusCaptureAtRest_);
+  pManager->AddRestProcess(&muMinusCapture_);
 
   // Tau Plus Physics
   pManager = G4TauPlus::TauPlus()->GetProcessManager();
diff --git a/src/G4BeamTestPhysicsList.cxx b/src/G4BeamTestPhysicsList.cxx
index 5a5befd..0c27b0b 100644
--- a/src/G4BeamTestPhysicsList.cxx
+++ b/src/G4BeamTestPhysicsList.cxx
@@ -1,15 +1,3 @@
-/**
- * Copyright (C) 2009
- * The IceCube collaboration
- * ID: $Id: G4BeamTestPhysicsList.cxx 152849 2017-01-20 21:44:25Z jgonzalez $
- *
- * @file G4BeamTestPhysicsList.cxx
- * @version $Rev: 152849 $
- * @date $Date: 2017-01-20 21:44:25 +0000 (Fri, 20 Jan 2017) $
- * @author Tilo Waldenmaier, Thomas Melzig, Javier Gonzalez
- */
-
-
 #include <globals.hh>
 #include <G4Version.hh>
 #include "G4BeamTestPhysicsList.h"
diff --git a/src/G4BeamTestRunManager.cxx b/src/G4BeamTestRunManager.cxx
index 753deb8..5ed2e05 100644
--- a/src/G4BeamTestRunManager.cxx
+++ b/src/G4BeamTestRunManager.cxx
@@ -1,20 +1,8 @@
-/**
- * Copyright (C) 2009
- * The IceCube collaboration
- * ID: $Id: G4BeamTestRunManager.cxx 162348 2018-04-25 20:09:46Z nega $
- *
- * @file G4BeamTestRunManager.cxx
- * @version $Rev: 162348 $
- * @date $Date: 2018-04-25 21:09:46 +0100 (Wed, 25 Apr 2018) $
- * @author Tilo Waldenmaier
- */
-
-
 // On Sun, to prevent conflict with ObjectSpace, G4Timer.hh has to be
 // loaded *before* globals.hh...
 #include "G4Timer.hh"
 
-#include <g4-tankresponse/g4classes/G4BeamTestRunManager.h>
+#include "G4BeamTestRunManager.h"
 #include <G4ParticleGun.hh>
 #include <G4Run.hh>
 
diff --git a/src/G4BeamTestTank.cxx b/src/G4BeamTestTank.cxx
index 4f12c7a..2fc5bed 100644
--- a/src/G4BeamTestTank.cxx
+++ b/src/G4BeamTestTank.cxx
@@ -1,25 +1,3 @@
-/**
- * Copyright (C) 2009
- * The IceCube collaboration
- * ID: $Id: G4BeamTestTank.cxx 158930 2017-10-20 01:17:41Z cweaver $
- *
- * @file G4BeamTestTank.cxx
- * @version $Rev: 158930 $
- * @date $Date: 2017-10-20 02:17:41 +0100 (Fri, 20 Oct 2017) $
- * @author Tilo Waldenmaier, Thomas Melzig
- */
-
-
-#include <g4-tankresponse/g4classes/G4BeamTestTank.h>
-#include <g4-tankresponse/g4classes/G4TankIceSD.h>
-#include <icetray/I3Units.h>
-
-#include <icetray/OMKey.h>
-#include <dataclasses/TankKey.h>
-#include <dataclasses/geometry/I3TankGeo.h>
-#include <dataclasses/geometry/I3Geometry.h>
-#include <dataclasses/I3Position.h>
-
 #include <G4LogicalVolume.hh>
 #include <G4PVPlacement.hh>
 #include <G4Material.hh>
@@ -33,7 +11,10 @@
 
 #include <G4UserLimits.hh>
 
-#include <boost/foreach.hpp>
+/* #include <boost/foreach.hpp> */
+
+#include "G4BeamTestTank.h"
+#include "G4TankIceSD.h"
 
 //prevent gcc to make something stupid with pretended unused variables
 #ifdef __GNUC__
@@ -43,121 +24,26 @@
 #endif
 
 
-G4BeamTestTank::G4BeamTestTank(const TankKey& tankKey, const I3Geometry& geometry):
-tankKey_(tankKey), geometry_(geometry)
+G4BeamTestTank::G4BeamTestTank()
 {
-  const I3StationGeoMap& stationMap = geometry.stationgeo;
-  unsigned int tankID = tankKey.tank==TankKey::TankA?0:1;
-  I3StationGeoMap::const_iterator station_iter = stationMap.find(tankKey.string);
-  
-  if(station_iter==stationMap.end())
-  {
-    log_fatal("The requested station %d in not in the geometry!", tankKey.string);
-    return;
-  }
-  
-  if(station_iter->second.size()<tankID)
-  {
-    log_fatal("The number of tanks in station %d is not correct!", tankKey.string);
-    return;
-  }
-  const I3TankGeo& tankGeo = station_iter->second.at(tankID);
-  
   // Get tank dimensions
-  tankThickness_ = 0.5*CLHEP::cm;
-  tankHeight_  = (tankGeo.tankheight / I3Units::m) * CLHEP::m + tankThickness_;
-  innerRadius_ = (tankGeo.tankradius / I3Units::m) * CLHEP::m;
+  tankThickness_ = 0.5*CLHEP::cm; // TODO(shivesh) : check thickness
+  tankHeight_  = 76.83 * 2.54 * CLHEP::cm;
+  innerRadius_ = 32 * 2.54 * CLHEP::cm;
   outerRadius_ = innerRadius_ + tankThickness_;
-  
-  // Get fill and snow heights
-  fillHeight_    = (tankGeo.fillheight / I3Units::m) * CLHEP::m;
-  snowHeight_    = (tankGeo.snowheight / I3Units::m) * CLHEP::m;
-  perliteHeight_ = tankHeight_ - tankThickness_ - fillHeight_;
-  
+
+  // Get fill height
+  fillHeight_ = 55.1 * 2.54 * CLHEP::cm;
+
   // Set DOM dimensions
   glassOuterRadius_ = 6.5 * 2.54 * CLHEP::cm; // 6.5" outer glass sphere radius
   glassThickness_   = 0.5 * 2.54 * CLHEP::cm; // 0.5" glass sphere thickness
-  
-  // Calculate tank position (tank center)
-  // tankGeo.position corresponds to the average position of the two DOMs in a tank
-  position_.set((tankGeo.position.GetX() / I3Units::m) * CLHEP::m,
-                (tankGeo.position.GetY() / I3Units::m) * CLHEP::m,
-                (tankGeo.position.GetZ() / I3Units::m) * CLHEP::m - fillHeight_ + 0.5*tankHeight_);
-  
-  // Get positions of the doms relativ to tank center
-  BOOST_FOREACH(const OMKey& omKey, tankGeo.omKeyList_)
-  {
-    I3OMGeoMap::const_iterator omGeo_iter = geometry_.omgeo.find(omKey);
-    if(omGeo_iter==geometry_.omgeo.end())
-    {
-      log_error_stream(omKey << " is missing in Tank " << tankKey_);
-      continue;
-    }
-    
-    G4ThreeVector relPos(omGeo_iter->second.position.GetX() - tankGeo.position.GetX(),
-                         omGeo_iter->second.position.GetY() - tankGeo.position.GetY(),
-                         omGeo_iter->second.position.GetZ() - tankGeo.position.GetZ());
-    
-    relDomPositions_[omKey] = (relPos / I3Units::m) * CLHEP::m;
-  }
-  
-  //
-  // Calculate Delaunay points
-  //
-  G4ThreeVector triangleDir(NAN, NAN, NAN);
-  switch(station_iter->second.size())
-  {
-    case 1: // Single tank 
-    {
-      // Vector orthogonal to DOM positions
-      triangleDir.set(relDomPositions_.begin()->second.y(),
-                      -relDomPositions_.begin()->second.x(),
-                      0.0);
-      break;
-    }
-    case 2: // Two tanks
-    {
-      const I3TankGeo& neighborGeo = station_iter->second.at(tankID==0?1:0);
-      G4ThreeVector neighborPos(neighborGeo.position.GetX(),
-                                neighborGeo.position.GetY(),
-                                neighborGeo.position.GetZ());
-      
-      // Convert to G4 units
-      neighborPos *= CLHEP::m/I3Units::m;
-      
-      // Same z position same as other tank
-      neighborPos.setZ(position_.z());
-      
-      triangleDir = position_ - neighborPos;
-      break;
-    }
-    default:
-    {
-      log_fatal("Invalid number of tanks (%zu) in station %d!",
-                station_iter->second.size(), tankKey_.string);
-      break;
-    }
-  }
-  
-  // side length
-  double triangleLength = 10.0 * CLHEP::m;
-  
-  triangleDir.setMag(0.5*triangleLength/cos(CLHEP::pi/6.0));
-  delaunayPoint1_ = position_ + triangleDir;
-  
-  // Rotate by 120 deg
-  triangleDir.rotateZ(CLHEP::pi/1.5);
-  delaunayPoint2_ = position_ + triangleDir;
-  
-  // Rotate by 120 deg
-  triangleDir.rotateZ(CLHEP::pi/1.5);
-  delaunayPoint3_ = position_ + triangleDir;
 }
 
 
 G4BeamTestTank::~G4BeamTestTank()
 {
-  
+
 }
 
 
@@ -166,217 +52,219 @@ 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 !!!!
-  // Maybe do all of this as stepping action ??????
+  // TODO(shivesh): Maybe do all of this as stepping action ??????
   G4UserLimits* energyLimit = new G4UserLimits();
-  energyLimit->SetUserMinEkine(280.0 * CLHEP::keV);  // Cherenkov threshold of electrons in ice
-  
-  std::string tankName=boost::lexical_cast<std::string>(tankKey_);
-  
-  // Define plastic frame
+  energyLimit->SetUserMinEkine(264.1 * CLHEP::keV);  // Cherenkov threshold of electrons in water
+
+  // 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);
-  
-  // Define ice volume
-  G4Material* ice = G4Material::GetMaterial("Ice");
-  G4Tubs* solidIce = new G4Tubs(("solid_ice_" + tankName).c_str(),
+
+  // 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);
-  G4LogicalVolume* logIce =
-      new G4LogicalVolume(solidIce, ice, ("log_ice_" + tankName).c_str(), 0, 0, 0);
-  G4ThreeVector physIcePosition(0, 0, -0.5*tankHeight_ + tankThickness_ + 0.5*fillHeight_);
-  G4VPhysicalVolume* physIce ATTRIBUTE_UNUSED =
-      new G4PVPlacement(0, physIcePosition, logIce,
-                        ("ice_" + tankName).c_str(), tankLog_, false, 0);
-  
-  // Define perlite volume
-  G4Material* perlite = G4Material::GetMaterial("Perlite");
-  G4Tubs* solidPerlite = new G4Tubs(("solid_perlite_" + tankName).c_str(),
-                                    0.0 * CLHEP::m, innerRadius_, 0.5 * perliteHeight_,
+  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_);
+  G4VPhysicalVolume* physWater ATTRIBUTE_UNUSED =
+      new G4PVPlacement(0, physWaterPosition, logWater,
+                        ("water_" + tankName).c_str(), tankLog_, false, 0);
+
+  // Define air volume
+  G4Material* air = G4Material::GetMaterial("Air");
+  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);
-  G4LogicalVolume* logPerlite =
-      new G4LogicalVolume(solidPerlite, perlite, ("log_perlite_" + tankName).c_str(), 0, 0, 0);
-  G4ThreeVector physPerlitePosition(0, 0, -0.5 * tankHeight_ + 0.5 * CLHEP::cm + fillHeight_ +
-                                    0.5 * perliteHeight_);
-  G4VPhysicalVolume* physPerlite ATTRIBUTE_UNUSED =
-      new G4PVPlacement(0, physPerlitePosition, logPerlite,
-                        ("perlite_" + tankName).c_str(), tankLog_, false, 0);
-  
+  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_);
+  G4VPhysicalVolume* physAir_UNUSED =
+      new G4PVPlacement(0, physAirPosition, logAir,
+                        ("air_" + tankName).c_str(), 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(om_iter->second.x(), om_iter->second.y(), -0.5 * perliteHeight_);
-    G4ThreeVector lowerDOMpos(om_iter->second.x(), om_iter->second.y(),  0.5 * fillHeight_);
-    
-    domPosIce[omKey] = lowerDOMpos;
-    
-    std::string omName=boost::lexical_cast<std::string>(omKey);
-    
-    G4Sphere *upperglasssphere = new G4Sphere (("solid_dom_up_" + omName).c_str(),
-                                               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(),
-                                               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(),
-                                             0.0 * CLHEP::m, domInnerRadius,
+
+  // 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_);
+
+  // 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(),
+                                             0.0 * CLHEP::m, glassOuterRadius_,
                                              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(),
-                                             0.0  * CLHEP::m, domInnerRadius,
+  G4Sphere *lowerglasssphere = new G4Sphere (("solid_dom_lo_" + omName).c_str(),
+                                             0.0  * CLHEP::m, glassOuterRadius_,
                                              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);
-    G4LogicalVolume* logLowerGlass =
-        new G4LogicalVolume(lowerglasssphere, glass,
-                            ("log_dom_lo_" + omName).c_str(), 0, 0, 0);
-    G4LogicalVolume* logUpperDOM =
-        new G4LogicalVolume(upperdomsphere, effectiveDOM,
-                            ("log_inside_dom_up_" + omName).c_str(), 0, 0, 0);
-    G4LogicalVolume* logLowerDOM =
-        new G4LogicalVolume(lowerdomsphere, effectiveDOM,
-                            ("log_inside_dom_lo_" + omName).c_str(), 0, 0, 0);
-    G4VPhysicalVolume* physUpperGlass ATTRIBUTE_UNUSED =
-        new G4PVPlacement(0, upperDOMpos, logUpperGlass,
-                          ("dom_up_" + omName).c_str(), logPerlite, false, 0);
-    G4VPhysicalVolume* physLowerGlass ATTRIBUTE_UNUSED =
-        new G4PVPlacement(0, lowerDOMpos, logLowerGlass,
-                          ("dom_lo_" + omName).c_str(), logIce, false, 0);
-    G4VPhysicalVolume* physUpperDOM ATTRIBUTE_UNUSED =
-        new G4PVPlacement(0, G4ThreeVector(0,0,0), logUpperDOM,
-                          ("inside_dom_up_" + omName).c_str(), logUpperGlass, false, 0);
-    G4VPhysicalVolume* physLowerDOM ATTRIBUTE_UNUSED =
-        new G4PVPlacement(0, G4ThreeVector(0,0,0), logLowerDOM,
-                          ("inside_dom_lo_" + omName).c_str(), logLowerGlass, false, 0);
-    
-    // apply energy limits
-    logUpperGlass->SetUserLimits(energyLimit);
-    logLowerGlass->SetUserLimits(energyLimit);
-    logUpperDOM->SetUserLimits(energyLimit);
-    logLowerDOM->SetUserLimits(energyLimit);
-  }
-  
-  // Define sensitive detector
+
+  G4double domInnerRadius = glassOuterRadius_ - glassThickness_;
+  G4Sphere *upperdomsphere = new G4Sphere (("solid_inside_dom_up_" + omName).c_str(),
+                                           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(),
+                                           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);
+  G4LogicalVolume* logLowerGlass =
+      new G4LogicalVolume(lowerglasssphere, glass,
+                          ("log_dom_lo_" + omName).c_str(), 0, 0, 0);
+  G4LogicalVolume* logUpperDOM =
+      new G4LogicalVolume(upperdomsphere, effectiveDOM,
+                          ("log_inside_dom_up_" + omName).c_str(), 0, 0, 0);
+  G4LogicalVolume* logLowerDOM =
+      new G4LogicalVolume(lowerdomsphere, effectiveDOM,
+                          ("log_inside_dom_lo_" + omName).c_str(), 0, 0, 0);
+  G4VPhysicalVolume* physUpperGlass ATTRIBUTE_UNUSED =
+      new G4PVPlacement(0, upperDOMpos, logUpperGlass,
+                        ("dom_up_" + omName).c_str(), logAir, false, 0);
+  G4VPhysicalVolume* physLowerGlass ATTRIBUTE_UNUSED =
+      new G4PVPlacement(0, lowerDOMpos, logLowerGlass,
+                        ("dom_lo_" + omName).c_str(), logWater, false, 0);
+  G4VPhysicalVolume* physUpperDOM ATTRIBUTE_UNUSED =
+      new G4PVPlacement(0, G4ThreeVector(0,0,0), logUpperDOM,
+                        ("inside_dom_up_" + omName).c_str(), logUpperGlass, false, 0);
+  G4VPhysicalVolume* physLowerDOM ATTRIBUTE_UNUSED =
+      new G4PVPlacement(0, G4ThreeVector(0,0,0), logLowerDOM,
+                        ("inside_dom_lo_" + omName).c_str(), logLowerGlass, false, 0);
+
+  // apply energy limits
+  logUpperGlass->SetUserLimits(energyLimit);
+  logLowerGlass->SetUserLimits(energyLimit);
+  logUpperDOM->SetUserLimits(energyLimit);
+  logLowerDOM->SetUserLimits(energyLimit);
+  // }
+
+  // Define sensitive detector TODO(shivesh): make the PMT the SD
   G4SDManager* sdManager = G4SDManager::GetSDMpointer();
-  iceSD_ = new G4TankIceSD(("ice_SD_" + tankName).c_str(), domPosIce);
+  iceSD_ = new G4TankIceSD(("ice_SD_" + tankName).c_str());
   sdManager->AddNewDetector(iceSD_);
-  logIce->SetSensitiveDetector(iceSD_);
-  
+  logWater->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 
+  // 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();
-  
+
   G4VPhysicalVolume* tankPhys = new G4PVPlacement(0, tankPos, tankLog_,
                                                   ("tank_" + tankName).c_str(),
                                                   mother->GetLogicalVolume(), false, 0);
-  
+
   // apply energy limits
   tankLog_->SetUserLimits(energyLimit);
-  logPerlite->SetUserLimits(energyLimit);
-  logIce->SetUserLimits(energyLimit);
-  
-  return tankPhys;
-}
-
-
-double G4BeamTestTank::GetNumCherenkov(const OMKey& omKey)
-{
-  return std::max(iceSD_->GetNumCherenkov(omKey), 0.);
-}
-
+  logAir->SetUserLimits(energyLimit);
+  logWater->SetUserLimits(energyLimit);
 
-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;
+  return tankPhys;
 }
 
 
-double G4BeamTestTank::GetY_I3()
+double G4BeamTestTank::GetNumCherenkov(/* const OMKey& omKey */)
 {
-  return (position_.y() / CLHEP::m) * I3Units::m;
+  return std::max(iceSD_->GetNumCherenkov(/* omKey */), 0.);
 }
 
 
-double G4BeamTestTank::GetZ_I3()
+double G4BeamTestTank::GetNumCherenkovWeight(/* const OMKey& omKey */)
 {
-  return (position_.z() / CLHEP::m) * I3Units::m;
+  return std::max(iceSD_->GetNumCherenkovWeight(/* omKey */), 0.);
 }
 
 
-double G4BeamTestTank::GetTankHeight_I3()
+double G4BeamTestTank::GetEDep_G4(/* const OMKey& omKey */)
 {
-  return (tankHeight_ / CLHEP::m) * I3Units::m;
+  return std::max(iceSD_->GetEDep(/* omKey */), 0.);
 }
 
 
-double G4BeamTestTank::GetTankRadius_I3()
+double G4BeamTestTank::GetTime_G4(/* const OMKey& omKey */)
 {
-  return (outerRadius_ / CLHEP::m) * I3Units::m;
+  return iceSD_->GetTime(/* omKey */);
 }
 
 
-double G4BeamTestTank::GetSnowHeight_I3()
-{
-  return (snowHeight_ / CLHEP::m) * I3Units::m;
-}
+/* 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 b16098c..e064628 100644
--- a/src/G4BeamTestUserSteppingAction.cxx
+++ b/src/G4BeamTestUserSteppingAction.cxx
@@ -1,18 +1,4 @@
-/**
- * Copyright (C) 2011
- * The IceCube collaboration
- * ID: $Id$
- *
- * @file G4BeamTestUserSteppingAction.cxx
- * @version $Revision$
- * @date $Date$
- * @author Thomas Melzig
- *
- * $LastChangedBy$
- */
-
-
-#include <g4-tankresponse/g4classes/G4BeamTestUserSteppingAction.h>
+#include "G4BeamTestUserSteppingAction.h"
 
 #include "G4Step.hh"
 #include "G4Track.hh"
diff --git a/src/G4BeamTestUserTrackingAction.cxx b/src/G4BeamTestUserTrackingAction.cxx
index d554152..c6317df 100644
--- a/src/G4BeamTestUserTrackingAction.cxx
+++ b/src/G4BeamTestUserTrackingAction.cxx
@@ -1,18 +1,4 @@
-/**
- * Copyright (C) 2011
- * The IceCube collaboration
- * ID: $Id$
- *
- * @file G4BeamTestUserTrackingAction.cxx
- * @version $Revision$
- * @date $Date$
- * @author Thomas Melzig
- *
- * $LastChangedBy$
- */
-
-
-#include <g4-tankresponse/g4classes/G4BeamTestUserTrackingAction.h>
+#include "G4BeamTestUserTrackingAction.h"
 
 #include "G4Track.hh"
 #include "G4UserLimits.hh"
diff --git a/src/G4Interface.cxx b/src/G4Interface.cxx
index 36ba3d9..9f44d3e 100644
--- a/src/G4Interface.cxx
+++ b/src/G4Interface.cxx
@@ -1,23 +1,12 @@
-/**
- * Copyright (C) 2009
- * The IceCube collaboration
- * ID: $Id: G4Interface.cxx 149388 2016-08-18 21:50:04Z jgonzalez $
- *
- * @file G4Interface.cxx
- * @version $Rev: 149388 $
- * @date $Date: 2016-08-18 22:50:04 +0100 (Thu, 18 Aug 2016) $
- * @author Tilo Waldenmaier
- */
-
-#include <dataclasses/physics/I3Particle.h>
-
-#include <g4-tankresponse/g4classes/G4Interface.h>
-#include <g4-tankresponse/g4classes/G4BeamTestDetectorConstruction.h>
-#include <g4-tankresponse/g4classes/G4BeamTestPhysicsList.h>
-#include <g4-tankresponse/g4classes/G4BeamTestUserTrackingAction.h>
-#include <g4-tankresponse/g4classes/G4BeamTestUserSteppingAction.h>
-
-#include <icetray/I3Logging.h>
+/* #include <dataclasses/physics/I3Particle.h> */
+
+#include "G4Interface.h"
+#include "G4BeamTestDetectorConstruction.h"
+#include "G4BeamTestPhysicsList.h"
+#include "G4BeamTestUserTrackingAction.h"
+#include "G4BeamTestUserSteppingAction.h"
+
+/* #include <icetray/I3Logging.h> */
 
 #ifdef G4VIS_USE
 #include <G4VisExecutive.hh>
@@ -26,7 +15,7 @@
 #include <G4ParticleGun.hh>
 #include <G4ParticleTable.hh>
 #include <G4ParticleDefinition.hh>
-#include <G4UImanager.hh>
+/* #include <G4UImanager.hh> */
 
 
 G4Interface* G4Interface::g4Interface_ = NULL;
@@ -40,11 +29,8 @@ G4Interface::G4Interface(const std::string& visMacro):
   // Visualization manager
 #ifdef G4VIS_USE
   visManager_ = NULL;
-  if(!visMacro_.empty())
-  {
-    visManager_ = new G4VisExecutive();
-    visManager_->Initialize();
-  }
+  visManager_ = new G4VisExecutive();
+  visManager_->Initialize();
 #endif
 }
 
@@ -63,7 +49,8 @@ void G4Interface::InstallTank(G4BeamTestTank* tank)
 {
   if(initialized_)
   {
-    log_fatal("G4Interface aleady initialized. Cannot install tank!");
+    /* log_fatal("G4Interface aleady initialized. Cannot install tank!"); */
+    G4cout << "G4Interface aleady initialized. Cannot install tank!" << G4endl;
     return;
   }
 
@@ -92,117 +79,37 @@ void G4Interface::InitializeEvent()
 }
 
 
-void G4Interface::InjectParticle(const I3Particle& particle)
+void G4Interface::InjectParticle(
+        const std::string& particleName, const G4ThreeVector& particlePosition,
+        const G4ThreeVector& particleDirection, const G4double particleEnergy
+        )
 {
   if(!eventInitialized_)
   {
-    log_fatal("No event initialized. Cannot inject particle!");
+    /* log_fatal("No event initialized. Cannot inject particle!"); */
+    G4cout << "No event initialized. Cannot inject particle!" << G4endl;
     return;
   }
 
   G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable();
   G4ParticleDefinition* particleDef = NULL;
-  switch(particle.GetType())
-  {
-    case I3Particle::Gamma:
-       particleDef = particleTable->FindParticle("gamma");
-       break;
-    case I3Particle::EMinus:
-       particleDef = particleTable->FindParticle("e-");
-       break;
-    case I3Particle::EPlus:
-       particleDef = particleTable->FindParticle("e+");
-       break;
-    case I3Particle::MuMinus:
-       particleDef = particleTable->FindParticle("mu-");
-       break;
-    case I3Particle::MuPlus:
-       particleDef = particleTable->FindParticle("mu+");
-       break;
-    case I3Particle::PPlus:
-       particleDef = particleTable->FindParticle("proton");
-       break;
-    case I3Particle::PMinus:
-       particleDef = particleTable->FindParticle("anti_proton");
-       break;
-    case I3Particle::Neutron:
-       particleDef = particleTable->FindParticle("neutron");
-       break;
-#ifdef I3PARTICLE_SUPPORTS_PDG_ENCODINGS
-    case I3Particle::NeutronBar:
-#else
-    case 25:
-#endif
-       particleDef = particleTable->FindParticle("anti_neutron");
-       break;
-    case I3Particle::PiPlus:
-       particleDef = particleTable->FindParticle("pi+");
-       break;
-    case I3Particle::PiMinus:
-       particleDef = particleTable->FindParticle("pi-");
-       break;
-    case I3Particle::Pi0:
-       particleDef = particleTable->FindParticle("pi0");
-       break;
-    case I3Particle::KPlus:
-       particleDef = particleTable->FindParticle("kaon+");
-       break;
-    case I3Particle::KMinus:
-       particleDef = particleTable->FindParticle("kaon-");
-       break;
-    case I3Particle::K0_Long:
-       particleDef = particleTable->FindParticle("kaon0L");
-       break;
-    case I3Particle::K0_Short:
-       particleDef = particleTable->FindParticle("kaon0S");
-       break;
-    case I3Particle::NuE:
-       particleDef = particleTable->FindParticle("nu_e");
-       break;
-    case I3Particle::NuEBar:
-       particleDef = particleTable->FindParticle("anti_nu_e");
-       break;
-    case I3Particle::NuMu:
-       particleDef = particleTable->FindParticle("nu_mu");
-       break;
-    case I3Particle::NuMuBar:
-       particleDef = particleTable->FindParticle("anti_nu_mu");
-       break;
-    case I3Particle::NuTau:
-       particleDef = particleTable->FindParticle("nu_tau");
-       break;
-    case I3Particle::NuTauBar:
-       particleDef = particleTable->FindParticle("anti_nu_tau");
-       break;
-    default:
-      log_warn("Man, check out that strange particle \"%s\" ?!", particle.GetTypeString().c_str());
-      return;
-  }
-  
-  // Particle position in G4 units
-  G4ThreeVector position((particle.GetX() / I3Units::m) * CLHEP::m,
-                         (particle.GetY() / I3Units::m) * CLHEP::m,
-                         (particle.GetZ() / I3Units::m) * CLHEP::m);
+  particleDef = particleTable->FindParticle(particleName);
 
-  // Transform I3 coorinates to world system
-  position -= detector_->GetWorldOrigin();
-
-  G4ThreeVector direction(particle.GetDir().GetX(),
-                          particle.GetDir().GetY(),
-                          particle.GetDir().GetZ());
+  // Transform coordinates to world system
+  G4ThreeVector position = particlePosition - detector_->GetWorldOrigin();
 
   G4ParticleGun gun(1);
   gun.SetParticleDefinition(particleDef);
-  gun.SetParticleEnergy((particle.GetEnergy() / I3Units::GeV) * CLHEP::GeV);
+  gun.SetParticleEnergy(particleEnergy);
   gun.SetParticlePosition(position);
-  gun.SetParticleMomentumDirection(direction);
+  gun.SetParticleMomentumDirection(particleDirection);
 
-  log_trace("Injecting %s: x=%.2f m, y=%.2f m, z=%.2f m, E=%.3f MeV",
-            particle.GetTypeString().c_str(),
+  G4cout << "Injecting %s: x=%.2f m, y=%.2f m, z=%.2f m, E=%.3f MeV",
+            particleName.c_str(),
             position.x() / CLHEP::m,
             position.y() / CLHEP::m,
             position.z() / CLHEP::m,
-            gun.GetParticleEnergy() / CLHEP::MeV);
+            gun.GetParticleEnergy() / CLHEP::MeV;
 
   runManager_.InjectParticle(&gun);
 }
@@ -227,66 +134,73 @@ void G4Interface::Initialize()
 {
   if(initialized_)
   {
-    log_error("G4Interface has already been initialized. Ignoring this call!");
+    /* log_error("G4Interface has already been initialized. Ignoring this call!"); */
+    G4cout << "G4Interface has already been initialized. Ignoring this call!" << G4endl;
     return;
   }
 
-  log_debug("Init geometry ...");
+  /* log_debug("Init geometry ..."); */
+  G4cout << "Init geometry ..." << G4endl;
   runManager_.SetUserInitialization(detector_);
 
-  log_debug("Init physics list ...");
+  /* log_debug("Init physics list ..."); */
+  G4cout << "Init physics list ..." << G4endl;
   runManager_.SetUserInitialization(new G4BeamTestPhysicsList());
 
-  log_debug("Init UserTrackingAction ...");
+  /* log_debug("Init UserTrackingAction ..."); */
+  G4cout << "Init UserTrackingAction ..." << G4endl;
   runManager_.SetUserAction(new G4BeamTestUserTrackingAction());
 
-  log_debug("Init UserSteppingAction ...");
+  /* log_debug("Init UserSteppingAction ..."); */
+  G4cout << "Init UserSteppingAction ..." << G4endl;
   runManager_.SetUserAction(new G4BeamTestUserSteppingAction());
 
   // Initialize G4 kernel
-  log_debug("Init run manager ...");
+  /* log_debug("Init run manager ..."); */
+  G4cout << "Init run manager ..." << G4endl;
   runManager_.Initialize();
 
-  // Set verbosity
-  int verboseLevel = 0;
-  switch (GetIcetrayLogger()->LogLevelForUnit("G4Interface"))
-  {
-    case I3LOG_FATAL:
-    case I3LOG_ERROR:
-    case I3LOG_WARN:
-    case I3LOG_INFO:
-    case I3LOG_NOTICE:
-    default:
-      verboseLevel = 0;
-      break;
-    case I3LOG_DEBUG:
-      verboseLevel = 1;
-      break;
-    case I3LOG_TRACE:
-      verboseLevel = 2;
-      break;
-  }
-
-  runManager_.SetVerboseLevel(verboseLevel);
-  G4EventManager::GetEventManager()->SetVerboseLevel(verboseLevel);
-  G4EventManager::GetEventManager()->GetStackManager()->SetVerboseLevel(verboseLevel);
-  G4EventManager::GetEventManager()->GetTrackingManager()->SetVerboseLevel(verboseLevel);
-#ifdef G4VIS_USE
-  if(visManager_) visManager_->SetVerboseLevel(verboseLevel);
-#endif
+  /* // Set verbosity */
+  /* int verboseLevel = 0; */
+  /* switch (GetIcetrayLogger()->LogLevelForUnit("G4Interface")) */
+  /* { */
+  /*   case I3LOG_FATAL: */
+  /*   case I3LOG_ERROR: */
+  /*   case I3LOG_WARN: */
+  /*   case I3LOG_INFO: */
+  /*   case I3LOG_NOTICE: */
+  /*   default: */
+  /*     verboseLevel = 0; */
+  /*     break; */
+  /*   case I3LOG_DEBUG: */
+  /*     verboseLevel = 1; */
+  /*     break; */
+  /*   case I3LOG_TRACE: */
+  /*     verboseLevel = 2; */
+  /*     break; */
+  /* } */
+
+  // TODO(shivesh): verbosity
+/*   runManager_.SetVerboseLevel(verboseLevel); */
+/*   G4EventManager::GetEventManager()->SetVerboseLevel(verboseLevel); */
+/*   G4EventManager::GetEventManager()->GetStackManager()->SetVerboseLevel(verboseLevel); */
+/*   G4EventManager::GetEventManager()->GetTrackingManager()->SetVerboseLevel(verboseLevel); */
+/* #ifdef G4VIS_USE */
+/*   if(visManager_) visManager_->SetVerboseLevel(verboseLevel); */
+/* #endif */
 
   // Execute visualization macro (if specified)
-  if(!visMacro_.empty())
-  {
-    G4UImanager* uim = G4UImanager::GetUIpointer();
-
-    // Checking geometry
-    uim->ApplyCommand("/geometry/test/grid_test");
-
-    // Execute visualization macro
-    std::string visCmd = "/control/execute " + visMacro_;
-    uim->ApplyCommand(visCmd.c_str());
-  }
+  /* if(!visMacro_.empty()) */
+  /* { */
+  /*   G4UImanager* uim = G4UImanager::GetUIpointer(); */
+  /*  */
+  /*   // Checking geometry */
+  /*   uim->ApplyCommand("/geometry/test/grid_test"); */
+  /*  */
+  /*   // Execute visualization macro */
+  /*   std::string visCmd = "/control/execute " + visMacro_; */
+  /*   uim->ApplyCommand(visCmd.c_str()); */
+  /* } */
 
   initialized_ = true;
 }
diff --git a/src/G4TankIceSD.cxx b/src/G4TankIceSD.cxx
index 4a1be0b..2974837 100644
--- a/src/G4TankIceSD.cxx
+++ b/src/G4TankIceSD.cxx
@@ -1,17 +1,5 @@
-/**
- * Copyright (C) 2009
- * The IceCube collaboration
- * ID: $Id: G4TankIceSD.cxx 152814 2017-01-19 21:34:52Z jgonzalez $
- *
- * @file G4TankIceSD.cxx
- * @version $Rev: 152814 $
- * @date $Date: 2017-01-19 21:34:52 +0000 (Thu, 19 Jan 2017) $
- * @author Tilo Waldenmaier, Thomas Melzig
- */
-
-
-#include <g4-tankresponse/g4classes/G4TankIceSD.h>
-#include <g4-tankresponse/g4classes/G4BeamTestTank.h>
+#include "G4TankIceSD.h"
+#include "G4BeamTestTank.h"
 
 #include <G4Step.hh>
 #include <G4HCofThisEvent.hh>
@@ -25,7 +13,7 @@
 #include <G4OpticalPhoton.hh>
 #include "G4Poisson.hh"
 
-#include <icetray/I3Units.h>
+/* #include <icetray/I3Units.h> */
 
 #include "G4Version.hh"
 
@@ -38,8 +26,8 @@
 #endif
 
 
-G4TankIceSD::G4TankIceSD(G4String name, const std::map<OMKey, G4ThreeVector>& domPositions)
-  : G4VSensitiveDetector(name), domPositions_(domPositions)
+G4TankIceSD::G4TankIceSD(G4String name/* , const std::map<OMKey, G4ThreeVector>& domPositions */)
+  : G4VSensitiveDetector(name)/* , domPositions_(domPositions) */
 {}
 
 
@@ -48,14 +36,14 @@ 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_[om_iter->first] = 0;
-    cogTime_[om_iter->first] = 0;
-    cherenkovCounter_[om_iter->first] = 0;
-    cherenkovCounterWeight_[om_iter->first] = 0;
-  }
+  // 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;
+  // }
 }
 
 
@@ -76,32 +64,32 @@ G4bool G4TankIceSD::ProcessHits(G4Step* aStep, G4TouchableHistory*)
   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(om_iter->second.x() - localPosition.x(), 2) +
-                           pow(om_iter->second.y() - localPosition.y(), 2));
-    G4double height = (om_iter->second.z() - localPosition.z());
-
-    sumEdep_[om_iter->first] += edep;
-    cogTime_[om_iter->first] += edep*time;
-    cherenkovCounterWeight_[om_iter->first] += GetProbability(radius, height) * cherenkovNumber;
-    cherenkovCounter_[om_iter->first] += cherenkovNumber;
-  }
+  // 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)
+  // std::map<OMKey, G4ThreeVector>::const_iterator om_iter;
+  // for(om_iter=domPositions_.begin(); om_iter!=domPositions_.end(); ++om_iter)
+  // {
+  if(sumEdep_>0)
   {
-    if(sumEdep_[om_iter->first]>0)
-    {
-      cogTime_[om_iter->first] /= sumEdep_[om_iter->first];
-    }
+    cogTime_ /= sumEdep_;
   }
+  // }
 }
 
 
@@ -122,7 +110,7 @@ G4double G4TankIceSD::GetCerenkovNumber(G4Step* aStep)
   G4MaterialPropertyVector* Rindex = aMaterialPropertiesTable->GetProperty("RINDEX");
   if (!Rindex) return 0.0;
 
-  const G4double Rfact = 369.81 / (CLHEP::eV * CLHEP::cm);
+  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;
@@ -147,7 +135,8 @@ G4double G4TankIceSD::GetCerenkovNumber(G4Step* aStep)
 #endif
   if (nMin!=nMax)
   {
-    log_error("Support for energy dependent refraction index not yet implemented!");
+    /* log_error("Support for energy dependent refraction index not yet implemented!"); */
+    G4cout << "Support for energy dependent refraction index not yet implemented!";
     return 0.0;
   }
 
@@ -163,6 +152,7 @@ G4double G4TankIceSD::GetCerenkovNumber(G4Step* aStep)
 }
 
 
+// TODO(shivesh): what is this
 G4double G4TankIceSD::GetProbability(G4double radius, G4double height)
 {
   height = 0.90 - height / CLHEP::m;