Hi Gabe -
I am done with cross-section stuff. Think the developers were not quite clear on where each element
needed for the calculation should come from and code was somewhat confused. It is actually quite simple now.
I tried to generate splines via GSL and, looking at a couple of values and comparing with what we had before,
it seems consistent. So everything works again and it is also logically structured and maintainable.
Please make sure you also svn update Conventions (had to add a new phase space enumeration).
Also I was wrong in what I said before - all 4 cross-sections are actually used in the generation code.
Have removed the arrays nevertheless.
There is no cross section calculation in the hadron tensor code any more,
so if you update / clean-up the kinematical selection, everything should be done via
NievesSimoVacasMECPXSec2016. Notice its inputs.
MEC interactions with only neutrino+target codes trigger calculation of full cross-section.
If you set a resonance, this triggers the calculation of the PDD component.
If you set a {np} hit pair it triggers the calculation of the pn component (either in full or PDD mode)
The input interaction is the only source of info for the interaction.
This is opposite to how the generation code was working (taking array
of all 4, decide, then set resonance or hit pair PDG. Now setting these (trial values) comes first.
Make sure you make a clean distinction between trial values (pn,delta) and the ones set in the output
interaction (especially if you abort with trial values set). You can make a copy of the entire Interaction
object to set trial values of delta,pn,T,costheta - then copy these to the Interaction object in the output event.
Also notice that in the acceptance / rejection loop you all you need to do is something like
Kinematics * kinematics = interaction->KinePtr();
kinematics->SetKV(kKVTl, T);
kinematics->SetKV(kKVctl, costh);
double xsec = xsec_model->XSec(fInteraction, kPSTlctl);
If there is any complicated business setting up 4-vectors etc you can ignore / delete.
Notice the second argument in XSec() - ensures that the cross-section used is d2sigma/dT_l dcostheta_l.
(as it should if you throw T_l, costheta_l from a uniform distribution).
If one selects different kinematical variables, then one should adjust the second argument to get the
appropriate cross-section transformation. The code does only d2sigma/dT_l dcostheta_l at the moment,
but at least the error will be flagged.
Well, good luck!
cheers
C
On 23 Mar 2016, at 17:13, Gabriel Nathan Perdue <[log in to unmask]<mailto:[log in to unmask]>> wrote:
Oh, I totally agree we should merge the generators. I'm just saying I wish we had a cleaner way of switching models that didn't require long if-else chains full of strings. But, I propose punting on that for now. When/if we think of what we want to do instead, we can also do it for coherent pions.
I don't fully understand your comments about the four cross sections. They use the values to set some exclusive state tags (was the event due to pionless delta-decay) that at first glance require computing all four values (maybe it is 2v2). But `XSec()` returns just one. I'd hate to call the method all over again with a flag to return a value I already computed and threw away. It isn't obvious to me what the solution is here, but maybe I'm missing something.
I really liked the alternate tune structure you proposed, but given our current workload, it seems hard to implement in time. Does anyone at a university have a bunch of bright undergrad CS majors we could tap for a project like that? I have this feeling we have a number of problems in GENIE that are not really physics problems and it is hard to tap physicists to work on them, but someone who was "just" a really coder could do them.
pax
Gabe
On Mar 23, 2016, at 7:53 AM, [log in to unmask]<mailto:[log in to unmask]> wrote:
Hi Gabe: OK I’ll leave the generator merging/clean-up to you then and I will focus on finishing and testing the
differential and integrated cross-section bits. Regarding the cross-sections: Yes, they have an array of 4 cross-section.
Really, they are 2 + 2, so there two possible computations from what I understand. I already included a `ComputeDelta’
config variable (default off). If ComputeDelta is off it returns element 0. I will tweak the code so that if ComputeDelta is on
it returns element 2. Element 1 is auxiliary for the calculation of element 0, and element 3 is auxiliary for the calculation of
element 2. I will also get rid of their array… Regarding merging the generators: Think we need to bite the bullet and merge
them - I think this really is the only convenient way to switch from one model to another by changing only the MEC
cross-section model. This does reflect a limitation of GENIE: We can’t support multiple tunes (with different model params,
different cross-section, different process generators event) specified and set via a single name. The document I circulated
some time ago outlines some ideas for addressing that, but we also need do schedule the work. Otherwise, changing between
the old and new tune in our fancy new v3.0.0 will be a real joy for our users!
cheers
C.
On 22 Mar 2016, at 22:48, Gabriel Nathan Perdue <[log in to unmask]<mailto:[log in to unmask]><mailto:[log in to unmask]>> wrote:
Sure, folks can follow along if interested.
I had a parade of people in my office today, but I started to merge the generators. There is a big list of TODOs in the code now - we need a hadron utils (or maybe just use the general utils class there) to hold some hadron state generation code that doesn't really belong in the generator class. Also, the big thing is figuring out how to deal with the fact the old MECTensor code expected the XSec method to basically return four(!) cross sections. We can either set the flags inside the new XSec method while all the pieces of the cross section calc are on hand or think of a slightly
different factorization.
Also, I'm handling the merged generator the same way we do Coherent pions, for example. I don't really like this way of doing it - you get long if-else blocks and strings everywhere (I hate strings, the compiler can't catch your spelling mistakes), and we talked about a fix (in the context of the coherent pions), but I can't remember what it was...
At any rate, I'll poke at this more either tomorrow afternoon or morning depending on whether anyone else is working on it. My #1 goal is to get it to run so we can compare some before and after physics distributions and see if things are broken.
pax
Gabe
On Mar 22, 2016, at 4:00 PM, Costas Andreopoulos <[log in to unmask]<mailto:[log in to unmask]><mailto:[log in to unmask]>> wrote:
Also forwarding to the collaboration e-mail list. We've been having discussions about reviews in the EB list out of inertia, but let's switch to the full list of active authors. Everyone is welcome to sit in the reviews, and contribute to (or at least inspect) new developments and the integration of external contributions. There is a wiki page which outlines the scheduled reviews.
cheers
C
--
Dr Constantinos Andreopoulos
Reader (Assoc. Professor) in Particle Physics
Univ. of Liverpool & STFC Rutherford Appleton Lab
+44-(0)7540 847333 (Mobile)
+44-(0)1235 445091 (Office/RAL)
+44-(0)1517 943201 (Office/Liverpool)
http://costas.andreopoulos.eu<http://costas.andreopoulos.eu/>
https://valor.pp.rl.ac.uk<https://valor.pp.rl.ac.uk/>
https://genie.hepforge.org<https://genie.hepforge.org/>
Sent from my iPhone
On 22 Mar 2016, at 19:07, Gabriel Nathan Perdue <[log in to unmask]<mailto:[log in to unmask]><mailto:[log in to unmask]><mailto:[log in to unmask]>> wrote:
Okay, I'll take a look at it.
pax
Gabe
On Mar 22, 2016, at 1:33 PM, [log in to unmask]<mailto:[log in to unmask]><mailto:[log in to unmask]><mailto:[log in to unmask]> wrote:
OK, just submitted the first batch of revisions. I will refrain from detailed comments on this experience, but _OMG_ !!
Clearly we need to do _much_ better in the future so that we steer projects to a better path and avoid this gruelling step.
I do not touch the MECTensor package, which is left as it was.
I copy changes over to the MEC package. Once we are done we can update MEC in the trunk and delete MECTensor.
If there is a digital equivalent of shredding I would be glad to do that with MECTensor.
- The old MECLoadHadTensor is now MECHadronTensor. This class used to do all sorts of things that had little to do
with loading the tensors. It where preppy much everything was ending to. I stripped it to absolutely basics (just keeping
and interrogating the collection of tensors). Very little was left. It is likely this class will need another shake up to generalise
slightly once we start including Martini tensors. But will not worry about this now.
- Added a NievesSimoVacasMECPXSec2016 where I assembled (from several places) all that encapsulated this model
in its 2016 incarnation (using hadron tensors). There is no need to copy over MECTensorsPXSec which will be deleted.
- Added MECUtils for general functions previously stuck in MECHadronTensor
- I joined the InteractionListGenerators into the one that was in the MEC package. Added a “SetDiNucleonCode”
config variable which is true by default. So default behaviour is old one, but can be tweaked to behave as
the “MECTensor” version
- Added a MECXSec as an interface between the GENIE Nieves (and later Martini) model and GSL
- I will fuse the two MECGenerator classes tomorrow so we can go from one MEC model to another by specifying a
change in the cross-section algorithm alone. I also want to change the old MECPXSec to a more meaningful name,
like EmpiricalMECPXSec201? so there is something that indicates the difference between that class and
NievesSimoVacasMECPXSec2016 (and, eventually MartiniMECPXSec201?)
- I have not tested or run any of the changes jut yet, but I am uncertain about a thing or two:
In NievesSimoVacasMECPXSec2016
https://genie.hepforge.org/trac/browser/generator/devel/branches/mec_devel/src/MEC/NievesSimoVacasMECPXSec2016.cxx
lines 54-73.
It is not really clear what is the expected input: The algorithm needs Enu, T(mu - for numu CC) and costhetamu.
But it seems to extract costheta from the f/s lepton 4-momentum.
I suspect this is not correct (or leads to hacks in the generation code).
You want to evaluate the cross-section to generate f/s lepton 4-momentum.
Therefore, you want a value of costheta before filling-up the f/s lepton 4-momentum.
I was expecting a call to Kinematics::GetKV(KineVar_t )
Because it is not clear what is expected, it is not easy to complete the GSL interface. Lines 155-156.
https://genie.hepforge.org/trac/browser/generator/devel/branches/mec_devel/src/MEC/MECXSec.cxx
Also, where is the calculation of kinematics limits for T and costheta (needed in lines 77-78 above)?
Gabe I am not touching MEC today - I need to leave very soon as I am travelling to Liverpool
very early tomorrow morning. It is yours if you want to continue with this e.g. addressing the above and testing
that we can generate return meaningful differential and total xsecs)?
cheers
Costas
--
Dr. Constantinos Andreopoulos
Reader (Associate Professor) in Experimental Particle Physics
University of Liverpool & STFC Rutherford Appleton Laboratory
+44-(0)7540-847333 (Mobile)
+44-(0)1235-445091 (Office/RAL)
+44-(0)1517-943201 (Office/Liverpool)
http://costas.andreopoulos.eu<http://costas.andreopoulos.eu/>
https://valor.pp.rl.ac.uk<https://valor.pp.rl.ac.uk/>
https://genie.hepforge.org<https://genie.hepforge.org/>
Gabriel Perdue
Associate Scientist
Scientific Computing Division
Fermi National Accelerator Laboratory
PO Box 500, MS 234, Batavia, IL 60510, USA
[log in to unmask]<mailto:[log in to unmask]><mailto:[log in to unmask]>
Office: 630-840-6499
Cell: 630-605-8062
Connect with us!
Newsletter<http://news.fnal.gov/fermilab-at-work/> | Facebook<http://www.facebook.com/Fermilab> | Twitter<https://www.twitter.com/fermilab>
Gabriel Perdue
Associate Scientist
Scientific Computing Division
Fermi National Accelerator Laboratory
PO Box 500, MS 234, Batavia, IL 60510, USA
[log in to unmask]<mailto:[log in to unmask]>
Office: 630-840-6499
Cell: 630-605-8062
Connect with us!
Newsletter<http://news.fnal.gov/fermilab-at-work/> | Facebook<http://www.facebook.com/Fermilab> | Twitter<https://www.twitter.com/fermilab>
|
