Implementation of NucleusGenerator #435
Conversation
Implement of NucleusGenerator to replace the old FermiMover and VertexGenerator. NucleusGenerator is inherited from EventRecordVisitorI. This interface facilitates the use of various implementations of nuclear models. Now I create two implementations: NucleusGenINCL and NucleusGenTraditional. NucleusGenINCL nuclear model is used for both initial state and FSI model to have a consistent physics results. NucleusGenTraditional combines the FermiMover and VertexGenerator.
LiangLiu212
changed the title
| void INCLNucleus::init(){ | ||
| LOG("NuclData", pINFO) << "init()"; | ||
| theConfig_->init(); | ||
| theConfig_->setINCLXXDataFilePath("/root/inclxx/inclxx-v6.33.1-e5857a1/data"); // FIXME:: using config to set path |
There was a problem hiding this comment.
the FIXME comment is important here. We can't have code that relies on any particular installation location.
In some obscure places in GENIE we have some hardcoded paths, I believe, and it irks me to no end.
| // initialize according process Event in INCL | ||
|
|
||
| GHepParticle * nucleus = evrec->TargetNucleus(); | ||
| G4INCL::ParticleSpecies targetSpecies = G4INCL::ParticleSpecies(nucleus->Name()); |
There was a problem hiding this comment.
I'm a bit surprised that GHepParticle returns a "name" in a format that is exactly compatible with G4INCL. How robust is this? I'm always a bit leery of string handling of things like this -- is there a PDG isotope code interface available. Also, are all possibllities returned by GHepParticle (say more obscure isotopes) handled by G4INCL, or is this open for a hard, obscure, difficult to diagnose crash?
There was a problem hiding this comment.
Oh, also is this also robust against there not being a nuclear target? Say nu_e + e scattering. Or is that handled upstream and would never get to this point?
| theConfig_->setTargetA(targetSpecies.theA); | ||
| theConfig_->setTargetZ(targetSpecies.theZ); | ||
| theConfig_->setTargetS(targetSpecies.theS); | ||
| LOG("NuclData", pINFO) << "initialize()"; |
There was a problem hiding this comment.
Long term pINFO is probably too verbose for such milestone marker messages. GENIE already is very talkative and it makes seeing the wheat amongst the chaff very difficult
|
|
||
| if(nucleus_){ | ||
| if(!nucleus_->getStore()->getParticles().empty()){ | ||
| if(nucleus_->getA() == evrec->TargetNucleus()->A() && nucleus_->getZ() == evrec->TargetNucleus()->Z()){ |
There was a problem hiding this comment.
Here at least you're not using both nucleus_ and nucleus together. Perhaps the previously fetched pointer needn't be explicitly saved to a local variable.
| // FIXME the last two parameters need to be configed | ||
| // theConfig_, G4INCL::ParticleTable::getMaximumNuclearRadius(G4INCL::Proton, targetSpecies.theA, targetSpecies.theZ) | ||
| // G4INCL::NType | ||
| nucleus_ = new G4INCL::Nucleus(targetSpecies.theA, targetSpecies.theZ, targetSpecies.theS, |
There was a problem hiding this comment.
Didn't earlier code use nucleus_ to decide stuff? Doesn't this invalidate that? Secondly, is this a potential memory leak as the previous one just has its pointer lost? The use of "new" here implies to me that this is an ownership pointer, not just one of convenience.
| // INCL nucleus model sample all nucleons with r-p correlation | ||
| INCLNucleus *incl_nucleus = INCLNucleus::Instance(); | ||
| G4INCL::Nucleus *incl_nuc = incl_nucleus->getNuclues(); | ||
| TLorentzVector p4tgt; |
There was a problem hiding this comment.
MIght be useful here to explicitly call out, in a comment, the units used by INCL and GENIE to explain the use of 1000.
| p4->SetPz(p3.Pz()/1000.); | ||
| p4->SetE (hit_nucleon_energy/1000.); | ||
| LOG("NucleusGenINCL", pINFO) << "Momentum"; | ||
| p4->Print(); |
There was a problem hiding this comment.
Again consider PrintUtils Vec3AsString() so that can be part of the LOG line and controlled by the message service.
| LOG("HadronicVtx", pFATAL) | ||
| << "No particle with [A = " << A << ", Z = " << Z | ||
| << ", pdgc = " << ipdgc << "] in PDGLibrary!"; | ||
| assert(remn); |
There was a problem hiding this comment.
The earlier comment here about assert being a no-op for non-debug code makes this section not do what I think you want it to. Even if you're hoping that the branch never gets triggered it's best that it does what it is intended to do.
| fVertexGenerator = nullptr; | ||
| fFermiMover = dynamic_cast<const EventRecordVisitorI *> (this->SubAlg("FermiMover")); | ||
| fVertexGenerator = dynamic_cast<const EventRecordVisitorI *> (this->SubAlg("VertexGenerator")); | ||
|
|
There was a problem hiding this comment.
Just dropping a marker here, because I'm not sure of how SubAlg's work, that it should be checked that fFermiMover and fVertexGenerator get Configured() if there are non-default cases.
| fNucleusGen = nullptr; | ||
| //fNuclModel = dynamic_cast<const NucleusGenerator *> (this->SubAlg(nuclkey)); | ||
| // | ||
| LOG("NucleusGenerator", pINFO) << "Hello world!"; |
Implement of NucleusGenerator to replace the old FermiMover and VertexGenerator.
NucleusGenerator is inherited from EventRecordVisitorI. This interface facilitates the use of various implementations of nuclear models.
Now I create two implementations: NucleusGenINCL and NucleusGenTraditional.
NucleusGenINCL nuclear model is used for both initial state and FSI model to have a consistent physics results.
NucleusGenTraditional combines the FermiMover and VertexGenerator.