diff --git a/Examples/options/tut_detsim_SDT_Heed.py b/Examples/options/tut_detsim_SDT_Heed.py
index e3ee257b7aed41ae14648fb3accb82d6e57322f2..d994ae61418c40c14afcd40640223c872ba810c9 100644
--- a/Examples/options/tut_detsim_SDT_Heed.py
+++ b/Examples/options/tut_detsim_SDT_Heed.py
@@ -178,8 +178,8 @@ elif dedxoption == "TrackHeedSimTool":
#dedx_simtool.IonMobility_file ="/junofs/users/wxfang/MyGit/tmp/check_G4FastSim_20210121/CEPCSW/Digitisers/DigiGarfield/IonMobility_He+_He.txt"
dedx_simtool.gas_file ="he_90_isobutane_10.gas"
dedx_simtool.IonMobility_file ="IonMobility_He+_He.txt"
- dedx_simtool.save_mc = True##IF this is False then ...
- dedx_simtool.debug = False
+ dedx_simtool.save_mc = True
+ dedx_simtool.debug = False
dedx_simtool.sim_pulse = True
#dedx_simtool.model='/junofs/users/wxfang/MyGit/tmp/fork_cepcsw_20220418/CEPCSW/Digitisers/SimCurrentONNX/src/model_90He10C4H10_18mm.onnx'
dedx_simtool.model='model_90He10C4H10_18mm.onnx'
diff --git a/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.cpp b/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.cpp
index b6ab38e67ba1b60c2f7bf086256e954f870298e0..0c687f82701832b673c49e2a744566aa1d76b061 100644
--- a/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.cpp
+++ b/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.cpp
@@ -87,8 +87,6 @@ Edm4hepWriterAnaElemTool::BeginOfEventAction(const G4Event* anEvent) {
// reset
m_track2primary.clear();
-
- auto SimPIonCol = m_SimPrimaryIonizationCol.createAndPut();
}
@@ -97,8 +95,6 @@ Edm4hepWriterAnaElemTool::EndOfEventAction(const G4Event* anEvent) {
msg() << "mcCol size (after simulation) : " << mcCol->size() << endmsg;
// save all data
- auto SimPrimaryIonizationCol = m_SimPrimaryIonizationCol.get();
- //msg() << "SimPrimaryIonizationCol size ="<<SimPrimaryIonizationCol->size()<<endmsg;
// create collections.
auto trackercols = m_trackerCol.createAndPut();
auto calorimetercols = m_calorimeterCol.createAndPut();
diff --git a/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.h b/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.h
index 60cc930a4ab3fdf04b57a9591c1b11a28ae29cf4..e4415b17899012530817186c0428581b585bad62 100644
--- a/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.h
+++ b/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.h
@@ -15,7 +15,6 @@
#include "edm4hep/SimTrackerHitCollection.h"
#include "edm4hep/SimCalorimeterHitCollection.h"
#include "edm4hep/CaloHitContributionCollection.h"
-#include "edm4hep/SimPrimaryIonizationClusterCollection.h"
class Edm4hepWriterAnaElemTool: public extends<AlgTool, IAnaElemTool> {
@@ -130,8 +129,6 @@ private:
"DriftChamberHitsCollection",
Gaudi::DataHandle::Writer, this};
- // for ionized electron
- DataHandle<edm4hep::SimPrimaryIonizationClusterCollection> m_SimPrimaryIonizationCol{"SimPrimaryIonizationClusterCollection", Gaudi::DataHandle::Writer, this};
private:
// in order to associate the hit contribution with the primary track,
diff --git a/Simulation/DetSimDedx/src/TrackHeedSimTool.cpp b/Simulation/DetSimDedx/src/TrackHeedSimTool.cpp
index 3bafda261a99e92b6a84182be886a87b55243ac2..f66c582b2a90d19ad86ede717d7e538ca51d9de0 100644
--- a/Simulation/DetSimDedx/src/TrackHeedSimTool.cpp
+++ b/Simulation/DetSimDedx/src/TrackHeedSimTool.cpp
@@ -25,7 +25,10 @@ double TrackHeedSimTool::dndx(double betagamma) {return 0;}
double TrackHeedSimTool::dedx(const G4Step* Step)
{
-
+ if(m_beginEvt){
+ m_SimPrimaryIonizationCol = m_SimPrimaryIonizationColWriter.createAndPut();
+ m_beginEvt = false;
+ }
clock_t t0 = clock();
double de = 0;
float cm_to_mm = 10;
@@ -43,15 +46,12 @@ double TrackHeedSimTool::dedx(const G4Step* Step)
if(g4Track->GetKineticEnergy() <=0) return 0;
if(pdg_code == 11 && (tmp_str_pro=="phot" || tmp_str_pro=="hIoni" || tmp_str_pro=="eIoni" || tmp_str_pro=="muIoni" || tmp_str_pro=="ionIoni" ) ) return m_eps;//skip the electron produced by Ioni, because it is already simulated by TrackHeed
if(m_particle_map.find(pdg_code) == m_particle_map.end() ) return m_eps;
- edm4hep::SimPrimaryIonizationClusterCollection* SimPrimaryIonizationCol = nullptr;
edm4hep::MCParticleCollection* mcCol = nullptr;
try{
- SimPrimaryIonizationCol = const_cast<edm4hep::SimPrimaryIonizationClusterCollection*>(m_SimPrimaryIonizationCol.get());
mcCol = const_cast<edm4hep::MCParticleCollection*>(m_mc_handle.get());
}
catch(...){
- G4cout<<"Error! Can't find collection in event, please check it have been createAndPut() in Begin of event"<<G4endl;
- G4cout<<"SimPrimaryIonizationCol="<<SimPrimaryIonizationCol<<",mcCol="<<mcCol<<G4endl;
+ G4cout<<"Error! Can't find MCParticle collection in event, please check it have been createAndPut() in Begin of event"<<G4endl;
throw "stop here!";
}
@@ -194,7 +194,8 @@ double TrackHeedSimTool::dedx(const G4Step* Step)
clock_t t02 = clock();
while (m_track->GetCluster(xc, yc, zc, tc, nc, ec, extra)) {
//auto chit = SimHitCol->create();
- auto chit = SimPrimaryIonizationCol->create();
+ //auto chit = SimPrimaryIonizationCol->create();
+ auto chit = m_SimPrimaryIonizationCol->create();
chit.setTime(tc);
double cpos[3] = { cm_to_mm*( (xc - init_x)+position_x/CLHEP::cm) , cm_to_mm*((yc - init_y)+position_y/CLHEP::cm), cm_to_mm*(zc + position_z/CLHEP::cm)};
chit.setPosition(edm4hep::Vector3d(cpos));
@@ -395,7 +396,8 @@ StatusCode TrackHeedSimTool::initialize()
m_current_Parent_ID = -1;
m_change_track = false;
m_total_range = 0;
- m_isFirst = false;
+ m_isFirst = true;
+ m_beginEvt = true;
m_tot_edep = 0;
m_tot_length = 0;
m_pa_KE =0;
@@ -555,26 +557,16 @@ float* TrackHeedSimTool::NNPred(std::vector<float>& inputs)
void TrackHeedSimTool::endOfEvent() {
if(m_sim_pulse){
- edm4hep::SimPrimaryIonizationClusterCollection* SimPrimaryIonizationCol = nullptr;
- try{
- SimPrimaryIonizationCol = const_cast<edm4hep::SimPrimaryIonizationClusterCollection*>(m_SimPrimaryIonizationCol.get());
- }
- catch(...){
- G4cout<<"Error! Can't find collection in event, please check it have been createAndPut() in Begin of event"<<G4endl;
- G4cout<<"SimPrimaryIonizationCol="<<SimPrimaryIonizationCol<<G4endl;
- throw "stop here!";
- }
- if(m_debug) G4cout<<"SimPrimaryIonizationCol size="<<SimPrimaryIonizationCol->size()<<G4endl;
+ if(m_debug) G4cout<<"SimPrimaryIonizationCol size="<<m_SimPrimaryIonizationCol->size()<<G4endl;
clock_t t01 = clock();
std::vector<float> inputs;
std::vector<unsigned long> indexs_c;
std::vector<unsigned long> indexs_i;
std::map<unsigned long long, std::vector<std::pair<float, float> > > id_pulse_map;
- for (unsigned long z=0; z<SimPrimaryIonizationCol->size(); z++) {
- for (unsigned long k=0; k<SimPrimaryIonizationCol->at(z).electronCellID_size(); k++) {
- //auto simIon = SimIonizationCol->at(k);
- auto position = SimPrimaryIonizationCol->at(z).getElectronPosition(k);//mm
- auto cellId = SimPrimaryIonizationCol->at(z).getElectronCellID(k);
+ for (unsigned long z=0; z<m_SimPrimaryIonizationCol->size(); z++) {
+ for (unsigned long k=0; k<m_SimPrimaryIonizationCol->at(z).electronCellID_size(); k++) {
+ auto position = m_SimPrimaryIonizationCol->at(z).getElectronPosition(k);//mm
+ auto cellId = m_SimPrimaryIonizationCol->at(z).getElectronCellID(k);
TVector3 Wstart(0,0,0);
TVector3 Wend (0,0,0);
m_segmentation->cellposition(cellId, Wstart, Wend);
@@ -590,9 +582,9 @@ void TrackHeedSimTool::endOfEvent() {
float local_y = 0;
getLocal(wire_x, wire_y, position[0], position[1], local_x, local_y);
//std::cout<<"pos_z="<<pos_z<<",wire_x="<<wire_x<<",wire_y="<<wire_y<<",position[0]="<<position[0]<<",position[1]="<<position[1]<<",local_x="<<local_x<<",local_y="<<local_y<<",dr="<<sqrt(local_x*local_x+local_y*local_y)<<",dr1="<<sqrt( (wire_x-position[0])*(wire_x-position[0])+(wire_y-position[1])*(wire_y-position[1]) )<<std::endl;
- float m_x_scale = 1;
- float m_y_scale = 1;
- local_x = local_x/m_x_scale;//FIXME, default is 18mm x 18mm, the real cell size maybe a bit different. need konw size of cell for each layer, and do normalization
+ //float ext_x_scale = 1;//TODO, default is 18mm x 18mm, the real cell size maybe a bit different. need konw size of cell for each layer, and do normalization using ext_x_scale
+ //float ext_y_scale = 1;
+ local_x = local_x/m_x_scale;
local_y = local_y/m_y_scale;
float noise = CLHEP::RandGauss::shoot(0,1);
inputs.push_back(local_x);
@@ -611,7 +603,7 @@ void TrackHeedSimTool::endOfEvent() {
//tmp_pluse.setValue(tmp_amp);
//tmp_pluse.setType(SimIonizationCol->at(tmp_index).getType());
//tmp_pluse.setSimIonization(SimIonizationCol->at(tmp_index));
- auto ion_time = SimPrimaryIonizationCol->at(indexs_c.at(i)).getElectronTime(indexs_i.at(i));
+ auto ion_time = m_SimPrimaryIonizationCol->at(indexs_c.at(i)).getElectronTime(indexs_i.at(i));
id_pulse_map[indexs_c.at(i)].push_back(std::make_pair(tmp_time+ion_time,tmp_amp) );
}
inputs.clear();
@@ -632,7 +624,7 @@ void TrackHeedSimTool::endOfEvent() {
//tmp_pluse.setType(SimIonizationCol->at(tmp_index).getType());
//tmp_pluse.setSimIonization(SimIonizationCol->at(tmp_index));
//id_pulse_map[SimIonizationCol->at(tmp_index).getCellID()].push_back(std::make_pair(tmp_pluse.getTime(), tmp_pluse.getValue() ) );
- auto ion_time = SimPrimaryIonizationCol->at(indexs_c.at(i)).getElectronTime(indexs_i.at(i));
+ auto ion_time = m_SimPrimaryIonizationCol->at(indexs_c.at(i)).getElectronTime(indexs_i.at(i));
id_pulse_map[indexs_c.at(i)].push_back(std::make_pair(tmp_time+ion_time,tmp_amp) );
}
inputs.clear();
@@ -641,7 +633,7 @@ void TrackHeedSimTool::endOfEvent() {
delete [] res;
}
for(auto iter = id_pulse_map.begin(); iter != id_pulse_map.end(); iter++){
- edm4hep::MutableSimPrimaryIonizationCluster dcIonCls = SimPrimaryIonizationCol->at(iter->first);
+ edm4hep::MutableSimPrimaryIonizationCluster dcIonCls = m_SimPrimaryIonizationCol->at(iter->first);
for(unsigned int i=0; i< iter->second.size(); i++){
auto tmp_time = iter->second.at(i).first ;
auto tmp_amp = iter->second.at(i).second;
@@ -649,13 +641,14 @@ void TrackHeedSimTool::endOfEvent() {
dcIonCls.addToPulseAmplitude(tmp_amp);
}
if(dcIonCls.electronPosition_size() != dcIonCls.pulseTime_size()){
- G4cout<<"Error ion size != pulse size"<<G4endl;
+ G4cout<<"Error: ion size != pulse size"<<G4endl;
throw "stop here!";
}
}
clock_t t02 = clock();
if(m_debug) std::cout<<"time for Pulse Simulation=" << (double)(t02 - t01) / CLOCKS_PER_SEC <<" seconds"<< std::endl;
}
+ reset();
}
StatusCode TrackHeedSimTool::finalize()
diff --git a/Simulation/DetSimDedx/src/TrackHeedSimTool.h b/Simulation/DetSimDedx/src/TrackHeedSimTool.h
index 8696d4f665008c94256338302d4ca145d2879ce4..4568c4f9016842d9aa5d8ab017a88e7cdfc90447 100644
--- a/Simulation/DetSimDedx/src/TrackHeedSimTool.h
+++ b/Simulation/DetSimDedx/src/TrackHeedSimTool.h
@@ -53,11 +53,11 @@ class TrackHeedSimTool: public extends<AlgTool, IDedxSimTool> {
double dndx(double betagamma) override;
void getMom(float ee, float dx, float dy,float dz, float mom[3] );
void reset(){
+ m_beginEvt = true;
m_isFirst = true;
m_previous_track_ID = 0;
m_previous_KE = 0;
m_tot_edep = 0;
- //std::cout<<"m_tot_length="<<m_tot_length<<std::endl;
m_tot_length = 0;
}
void endOfEvent();
@@ -89,7 +89,8 @@ class TrackHeedSimTool: public extends<AlgTool, IDedxSimTool> {
Gaudi::Property<float> m_BField {this, "BField", -3};
Gaudi::Property<float> m_eps { this, "eps" , 1e-6 };//very small value, it is returned dedx for unsimulated step (may needed for SimTrackerHit)
// Output collections
- DataHandle<edm4hep::SimPrimaryIonizationClusterCollection> m_SimPrimaryIonizationCol{"SimPrimaryIonizationClusterCollection", Gaudi::DataHandle::Writer, this};
+ DataHandle<edm4hep::SimPrimaryIonizationClusterCollection> m_SimPrimaryIonizationColWriter{"SimPrimaryIonizationClusterCollection", Gaudi::DataHandle::Writer, this};
+ edm4hep::SimPrimaryIonizationClusterCollection* m_SimPrimaryIonizationCol;
// In order to associate MCParticle with contribution, we need to access MC Particle.
DataHandle<edm4hep::MCParticleCollection> m_mc_handle{"MCParticle", Gaudi::DataHandle::Writer, this};
@@ -102,19 +103,20 @@ class TrackHeedSimTool: public extends<AlgTool, IDedxSimTool> {
Sensor* m_sensor;
std::map<int, std::string> m_particle_map;
- int m_previous_track_ID=0;
- float m_previous_KE=0;
+ int m_previous_track_ID;
+ float m_previous_KE;
int m_current_track_ID;
int m_current_Parent_ID;
int m_pdg_code;
G4StepPoint* m_pre_point;
G4StepPoint* m_post_point;
G4double m_total_range;
- bool m_isFirst=true;
+ bool m_isFirst;
+ bool m_beginEvt;
bool m_change_track;
edm4hep::MCParticle m_mc_paricle;
- float m_tot_edep=0;
- float m_tot_length=0;
+ float m_tot_edep;
+ float m_tot_length;
float m_pa_KE;
G4double m_pre_x ;
@@ -147,12 +149,12 @@ class TrackHeedSimTool: public extends<AlgTool, IDedxSimTool> {
Gaudi::Property<bool> m_sim_pulse { this, "sim_pulse" , true };
Gaudi::Property<std::string> m_model_file{ this, "model", "model_test.onnx"};
Gaudi::Property<int> m_batchsize { this, "batchsize", 100};
- Gaudi::Property<float> m_time_scale { this, "time_scale", 99.0};
- Gaudi::Property<float> m_time_shift { this, "time_shift", 166.4};
- Gaudi::Property<float> m_amp_scale { this, "amp_scale" , 1e-2 };
- Gaudi::Property<float> m_amp_shift { this, "amp_shift" , 0 };
- Gaudi::Property<float> m_x_scale { this, "x_scale" , 5. };// in mm
- Gaudi::Property<float> m_y_scale { this, "y_scale" , 5. };// in mm
+ Gaudi::Property<float> m_time_scale { this, "time_scale", 503.0};
+ Gaudi::Property<float> m_time_shift { this, "time_shift", 814.0};
+ Gaudi::Property<float> m_amp_scale { this, "amp_scale" , 1.15 };
+ Gaudi::Property<float> m_amp_shift { this, "amp_shift" , 0.86 };
+ Gaudi::Property<float> m_x_scale { this, "x_scale" , 9. };// in mm
+ Gaudi::Property<float> m_y_scale { this, "y_scale" , 9. };// in mm
diff --git a/Simulation/DetSimInterface/include/DetSimInterface/IDedxSimTool.h b/Simulation/DetSimInterface/include/DetSimInterface/IDedxSimTool.h
index d4e8e9e4ae8b14b2c748256bcf75a536fa29c88b..4da01d532a69f4cba7078d82fb161fcc34a1634c 100644
--- a/Simulation/DetSimInterface/include/DetSimInterface/IDedxSimTool.h
+++ b/Simulation/DetSimInterface/include/DetSimInterface/IDedxSimTool.h
@@ -28,7 +28,6 @@ public:
virtual double dedx(const G4Step* aStep) = 0;
virtual double dedx(const edm4hep::MCParticle& mc) = 0;
virtual double dndx(double betagamma) = 0;
- virtual void reset() {}
virtual void endOfEvent() {}
};
diff --git a/Simulation/DetSimSD/src/DriftChamberSensitiveDetector.cpp b/Simulation/DetSimSD/src/DriftChamberSensitiveDetector.cpp
index c1d768021e0deb91712baab30235af225ca77821..44bf01d1016ae3578ef192f609a3163b6ea6b483 100644
--- a/Simulation/DetSimSD/src/DriftChamberSensitiveDetector.cpp
+++ b/Simulation/DetSimSD/src/DriftChamberSensitiveDetector.cpp
@@ -72,5 +72,4 @@ DriftChamberSensitiveDetector::ProcessHits(G4Step* step, G4TouchableHistory*) {
void
DriftChamberSensitiveDetector::EndOfEvent(G4HCofThisEvent* HCE) {
m_dedx_simtool->endOfEvent();
- m_dedx_simtool->reset();
}