update1

Detector/DetEcalMatrix/compact/det.xml

@@ -32,7 +32,7 @@
 
   <detectors>
     <detector id="1" name="CaloDetector" type="EcalMatrix" readout="CaloHitsCollection" vis="VisibleGreen" sensitive="true">
-      <position x="0"  y="0"  z="1835*mm+30*cm"/>
+      <position x="200*cm"  y="0"  z="0"/>
       <dimensions dx="30*cm"  dy="30*cm"  dz="30*cm"/>
     </detector>
   </detectors>
@@ -43,7 +43,7 @@
                     grid_size_x="1*cm" 
                     grid_size_y="1*cm" 
                     grid_size_z="1*cm"/>
-      <id>system:8,x:32:-6,y:-6,z:-6</id>
+      <id>system:8,x:-16,y:-16,z:-16</id>
     </readout>
   </readouts>
 

Detector/DetEcalMatrix/src/calorimeter/EcalMatrix.cpp

@@ -16,6 +16,7 @@
 #define MYDEBUG(x) std::cout << __FILE__ << ":" << __LINE__ << ": " << x << std::endl;
 #define MYDEBUGVAL(x) std::cout << __FILE__ << ":" << __LINE__ << ": " << #x << ": " << x << std::endl;
 
+using dd4hep::rec::LayeredCalorimeterData;
 static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector,
                                      xml_h e,
                                      dd4hep::SensitiveDetector sens) {
@@ -40,6 +41,41 @@ static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector,
     dd4hep::Transform3D transform(dd4hep::Rotation3D(),
                                   dd4hep::Position(pos.x(),pos.y(),pos.z()));
     dd4hep::PlacedVolume phv = motherVol.placeVolume(det_vol,transform);
+    //Create caloData object to extend driver with data required for reconstruction
+    LayeredCalorimeterData* caloData = new LayeredCalorimeterData ;
+    caloData->layoutType = LayeredCalorimeterData::BarrelLayout ;
+    caloData->inner_symmetry = 8 ;//nsides
+    caloData->outer_symmetry = 8; //nsides;
+    caloData->inner_phi0 = 0.;
+    caloData->outer_phi0 = 0.; 
+    caloData->gap0 = 0.; //FIXME
+    caloData->gap1 = 0.; //FIXME
+    caloData->gap2 = 0.; //FIXME
+    // extent of the calorimeter in the r-z-plane [ rmin, rmax, zmin, zmax ] in mm.
+    caloData->extent[0] = 10*( pos.x()-dim.dx() );// from cm to mm
+    caloData->extent[1] = 10*( pos.x()+dim.dx() );// from cm to mm
+    caloData->extent[2] = 0 ;
+    caloData->extent[3] = 10*( pos.z()+dim.dz() );// from cm to mm
+    std::cout<<"rmin="<<caloData->extent[0]<<",rmax="<<caloData->extent[1]<<",zmin="<<caloData->extent[2]<<",zmax="<<caloData->extent[3]<<std::endl;
+    dd4hep::Readout readout = sens.readout();
+    dd4hep::Segmentation seg = readout.segmentation();
+    double cellSize_x = ::atof( seg.segmentation()->parameter("grid_size_x")->value().c_str() ) * 10;// from cm to mm
+    double cellSize_y = ::atof( seg.segmentation()->parameter("grid_size_y")->value().c_str() ) * 10;// from cm to mm
+    double cellSize_z = ::atof( seg.segmentation()->parameter("grid_size_z")->value().c_str() ) * 10;// from cm to mm
+    int n_layer = int(2*dim.dx()*10/cellSize_x) ; // here the calorimeter is placed in barrel, so x direaction is layer direction
+    std::cout<<"cellx="<<cellSize_x<<",celly="<<cellSize_y<<",cellz="<<cellSize_z<<",dx="<<dim.dx()*10<<"mm,n_layer="<<n_layer<<std::endl;
+    for(int i=1 ; i <= n_layer; i++)
+    {
+        LayeredCalorimeterData::Layer caloLayer ;
+        caloLayer.distance = caloData->extent[0] + (i-0.5)*cellSize_x; //NEED TO START FROM ORIGIN, to mm
+        caloLayer.sensitive_thickness = cellSize_x ;
+        caloLayer.inner_thickness = cellSize_x ;
+        caloLayer.outer_thickness = cellSize_x ;
+        caloLayer.absorberThickness = cellSize_x;
+        caloLayer.cellSize0 = cellSize_y;
+        caloLayer.cellSize1 = cellSize_z;
+        caloData->layers.push_back(caloLayer); 
+    }
 
     if ( x_det.isSensitive() )   {
         dd4hep::SensitiveDetector sd = sens;
@@ -50,7 +86,7 @@ static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector,
         phv.addPhysVolID("system",x_det.id());
     }
     sdet.setPlacement(phv);
-    
+    sdet.addExtension< LayeredCalorimeterData >(caloData) ; 
     MYDEBUG("create_detector DONE. ");
     return sdet;
 }

Generator/src/GtGunTool.cpp

@@ -15,11 +15,12 @@ GtGunTool::initialize() {
         error() << "Please specify the list of particle names/pdgs" << endmsg;
         return StatusCode::FAILURE;
     }
+    /*
     if (m_energies.value().size() != m_particles.value().size()) {
         error() << "Mismatched energies and particles." << endmsg;
         return StatusCode::FAILURE;
     }
-
+    */
     // others should be empty or specify
     if (m_thetamins.value().size()
         && m_thetamins.value().size() != m_particles.value().size()) {
@@ -78,7 +79,9 @@ GtGunTool::mutate(MyHepMC::GenEvent& event) {
             }
         }
 
-        double energy = m_energies.value()[i];
+        //double energy = m_energies.value()[i];
+        double energy_min = m_energies_min.value()[0];
+        double energy_max = m_energies_max.value()[0];
 
         // create the MC particle
         edm4hep::MCParticle mcp = event.m_mc_vec.create();
@@ -92,16 +95,18 @@ GtGunTool::mutate(MyHepMC::GenEvent& event) {
         // mcp.setEndpoint();
 
         // assume energy is momentum
-        double p = energy;
+        double p = energy_min==energy_max ? energy_max : CLHEP::RandFlat::shoot(energy_min, energy_max);
         
         // direction
         // by default, randomize the direction
-        double costheta = CLHEP::RandFlat::shoot(-1, 1);
-        double phi = 360*CLHEP::RandFlat::shoot()*CLHEP::degree;
+        double theta = m_thetamins.value()[0]==m_thetamaxs.value()[0] ? m_thetamins.value()[0] : CLHEP::RandFlat::shoot(m_thetamins.value()[0], m_thetamaxs.value()[0]);
+        double phi =   m_phimins  .value()[0]==m_phimaxs  .value()[0] ? m_phimins  .value()[0] : CLHEP::RandFlat::shoot(m_phimins  .value()[0], m_phimaxs  .value()[0]);
+        double costheta = cos(theta*acos(-1)/180);
+        phi = phi*acos(-1)/180;
         double sintheta = sqrt(1.-costheta*costheta);
 
         // check if theta min/max is set
-
+        /*
         if (i < m_thetamins.value().size() 
             && i < m_thetamaxs.value().size()) {
             double thetamin = m_thetamins.value()[i];
@@ -129,11 +134,11 @@ GtGunTool::mutate(MyHepMC::GenEvent& event) {
                 << " cos(theta): " << costheta
                 << " phi: " << phi
                 << endmsg;
-        
+        */ 
         double px = p*sintheta*cos(phi);
         double py = p*sintheta*sin(phi);
         double pz = p*costheta;
-
+        std::cout<<"GenGt p="<<p<<", px="<<px<<",py="<<py<<",pz="<<pz<<",costheta="<<costheta<<std::endl;
         mcp.setMomentum(edm4hep::Vector3f(px,py,pz));
         // mcp.setMomentumAtEndpoint();
         // mcp.setSpin();

Generator/src/GtGunTool.h

@@ -35,6 +35,8 @@ private:
     Gaudi::Property<std::vector<std::string>> m_particles{this, "Particles"};
 
     Gaudi::Property<std::vector<double>> m_energies{this, "Energies"};
+    Gaudi::Property<std::vector<double>> m_energies_min{this, "EnergyMin"};
+    Gaudi::Property<std::vector<double>> m_energies_max{this, "EnergyMax"};
 
     Gaudi::Property<std::vector<double>> m_thetamins{this, "ThetaMins"};
     Gaudi::Property<std::vector<double>> m_thetamaxs{this, "ThetaMaxs"};

Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.cpp

@@ -173,7 +173,7 @@ Edm4hepWriterAnaElemTool::EndOfEventAction(const G4Event* anEvent) {
                     ++n_cal_hit;
                     auto edm_calo_hit = calo_col_ptr->create();
                     edm_calo_hit.setCellID(cal_hit->cellID);
-                    edm_calo_hit.setEnergy(cal_hit->energyDeposit);
+                    edm_calo_hit.setEnergy(cal_hit->energyDeposit/CLHEP::GeV);
                     float pos[3] = {cal_hit->position.x()/CLHEP::mm,
                                     cal_hit->position.y()/CLHEP::mm,
                                     cal_hit->position.z()/CLHEP::mm};