Skip to content
Snippets Groups Projects
Commit d231a04f authored by lintao@ihep.ac.cn's avatar lintao@ihep.ac.cn
Browse files

Merge branch 'master' into 'master'

Add option for gtgun to generate charged particles according to costheta distribution

See merge request !170
parents 48bcef3d 23a141a9
No related branches found
No related tags found
No related merge requests found
...@@ -73,17 +73,6 @@ GtGunTool::initialize() { ...@@ -73,17 +73,6 @@ GtGunTool::initialize() {
} }
// others should be empty or specify // others should be empty or specify
if (m_thetamins.value().size()
&& m_thetamins.value().size() != m_particles.value().size()) {
error() << "Mismatched thetamins and particles." << endmsg;
return StatusCode::FAILURE;
}
if (m_thetamaxs.value().size()
&& m_thetamaxs.value().size() != m_particles.value().size()) {
error() << "Mismatched thetamaxs and particles." << endmsg;
return StatusCode::FAILURE;
}
if (m_phimins.value().size() if (m_phimins.value().size()
&& m_phimins.value().size() != m_particles.value().size()) { && m_phimins.value().size() != m_particles.value().size()) {
error() << "Mismatched phimins and particles." << endmsg; error() << "Mismatched phimins and particles." << endmsg;
...@@ -95,6 +84,38 @@ GtGunTool::initialize() { ...@@ -95,6 +84,38 @@ GtGunTool::initialize() {
return StatusCode::FAILURE; return StatusCode::FAILURE;
} }
if (m_usecostheta.value()) {
if (m_costhetamins.value().size() != m_particles.value().size()) {
error() << "Mismatched CosthetaMins and particles." << endmsg;
return StatusCode::FAILURE;
}
if (m_costhetamaxs.value().size() != m_particles.value().size()) {
error() << "Mismatched CosthetaMaxs and particles." << endmsg;
return StatusCode::FAILURE;
}
for (int i=0; i<m_thetamins.value().size(); ++i) {
if ((m_thetamins.value()[i] < -1) || (m_thetamins.value()[i] > 1)) {
error() << "UseCostheta: ThetaMins has values outside the range [-1, 1]." << endmsg;
return StatusCode::FAILURE;
}
if ((m_thetamaxs.value()[i] < -1) || (m_thetamaxs.value()[i] > 1)) {
error() << "UseCostheta: ThetaMaxs has values outside the range [-1, 1]." << endmsg;
return StatusCode::FAILURE;
}
}
} else {
if (m_thetamins.value().size()
&& m_thetamins.value().size() != m_particles.value().size()) {
error() << "Mismatched thetamins and particles." << endmsg;
return StatusCode::FAILURE;
}
if (m_thetamaxs.value().size()
&& m_thetamaxs.value().size() != m_particles.value().size()) {
error() << "Mismatched thetamaxs and particles." << endmsg;
return StatusCode::FAILURE;
}
}
// Time // Time
if (m_times.value().size()==0){ if (m_times.value().size()==0){
for(int i=0; i<m_particles.value().size(); i++) m_times.value().push_back(0); for(int i=0; i<m_particles.value().size(); i++) m_times.value().push_back(0);
...@@ -243,16 +264,26 @@ GtGunTool::mutate(Gen::GenEvent& event) { ...@@ -243,16 +264,26 @@ GtGunTool::mutate(Gen::GenEvent& event) {
return false; return false;
} }
double costheta = 0;
double theta = 0;
if (m_usecostheta.value()) {
costheta = m_costhetamins.value()[i]==m_costhetamaxs.value()[i] ? m_costhetamins.value()[i] : CLHEP::RandFlat::shoot(m_costhetamins.value()[i], m_costhetamaxs.value()[i]);
} else {
theta = m_thetamins.value()[i]==m_thetamaxs.value()[i] ? m_thetamins.value()[i] : CLHEP::RandFlat::shoot(m_thetamins.value()[i], m_thetamaxs.value()[i]);
costheta = cos(theta*acos(-1)/180);
}
double phi = m_phimins.value()[i]==m_phimaxs.value()[i] ? m_phimins.value()[i] : CLHEP::RandFlat::shoot(m_phimins.value()[i], m_phimaxs.value()[i]);
double theta = m_thetamins.value()[i]==m_thetamaxs.value()[i] ? m_thetamins.value()[i] : CLHEP::RandFlat::shoot(m_thetamins.value()[i], m_thetamaxs.value()[i]);
double phi = m_phimins .value()[i]==m_phimaxs .value()[i] ? m_phimins .value()[i] : CLHEP::RandFlat::shoot(m_phimins .value()[i], m_phimaxs .value()[i]);
double costheta = cos(theta*acos(-1)/180);
double phi_ = phi*acos(-1)/180; double phi_ = phi*acos(-1)/180;
double sintheta = sqrt(1.-costheta*costheta); double sintheta = sqrt(1.-costheta*costheta);
double px = p*sintheta*cos(phi_); double px = p*sintheta*cos(phi_);
double py = p*sintheta*sin(phi_); double py = p*sintheta*sin(phi_);
double pz = p*costheta; double pz = p*costheta;
std::cout<<"GenGt p="<<p<<", px="<<px<<",py="<<py<<",pz="<<pz<<",theta="<<theta<<",phi="<<phi<<std::endl; if(m_usecostheta.value()) {
std::cout<<"GenGt p="<<p<<", px="<<px<<",py="<<py<<",pz="<<pz<<",costheta="<<costheta<<",phi="<<phi<<std::endl;
} else {
std::cout<<"GenGt p="<<p<<", px="<<px<<",py="<<py<<",pz="<<pz<<",theta="<<theta<<",phi="<<phi<<std::endl;
}
mcp.setMomentum(edm4hep::Vector3f(px,py,pz)); mcp.setMomentum(edm4hep::Vector3f(px,py,pz));
// mcp.setMomentumAtEndpoint(); // mcp.setMomentumAtEndpoint();
// mcp.setSpin(); // mcp.setSpin();
......
...@@ -57,6 +57,9 @@ private: ...@@ -57,6 +57,9 @@ private:
Gaudi::Property<std::vector<double>> m_phimins{this, "PhiMins"}; Gaudi::Property<std::vector<double>> m_phimins{this, "PhiMins"};
Gaudi::Property<std::vector<double>> m_phimaxs{this, "PhiMaxs"}; Gaudi::Property<std::vector<double>> m_phimaxs{this, "PhiMaxs"};
Gaudi::Property<bool> m_usecostheta{this, "UseCostheta", false};
Gaudi::Property<std::vector<double>> m_costhetamins{this, "CosthetaMins"};
Gaudi::Property<std::vector<double>> m_costhetamaxs{this, "CosthetaMaxs"};
// For time // For time
Gaudi::Property<std::vector<double>> m_times{this, "Times"}; Gaudi::Property<std::vector<double>> m_times{this, "Times"};
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment