[Tdis] [EXTERNAL] Inconsistencies in the TDIS proposal: Some things I found

Tue Aug 11 09:19:16 EDT 2020

Hi Wally,

Thank you for looking into this.

Your suggestion to write my own implementation is a very good one, and 
exactly what I did. This is where I found the discrepancies to the 
proposal (except the ones in the proposal itself).  I think there are 
two questions left for that:

1) Is my formula for F2p correct, or the one in the code I attached. 
I.e, is F2 =  x * Sum quarkcharge_i^2 
quarkdistributionfunction_i(x,Q^2)  (what I think) or is there an 
additional factor 2 (what the old code says).

2) Are there some other cuts I didn't realize

I think all other points I raised are simple mistakes in the original 
proposal or otherwise understood artifacts.

  Other than that, your code + my code  produce the same results. Q2 can 
probably only be answered by the person who originally came up with the 
table in the proposal. But if you could write down what you think the 
formula for F2p is, it would be a good cross check that I didn't just 
misunderstand the textbooks.

Best,

Jan

On 8/11/2020 12:30 AM, Wally Melnitchouk wrote:
> Hi Jan,
>
> Thanks for your emails, and apologies for the delay.
>
> The programs you appended look familiar -- I think I wrote the core of 
> them back in the 1990s -- and Tim developed them some half dozen years 
> ago (and possibly other people before and since) when we were working 
> on some related applications.
>
> The best suggestion I would give anyone interested in implementing the 
> meson exchange contributions to the DIS structure functions would be 
> to write some script yourself from scratch. The formulas are actually 
> very simple, if you look at the analytic expressions. The original 
> ancient Fortran could probably now be written in a few lines of python 
> code and done much more efficiently, and writing it yourself would 
> give you a much clearer idea of what is being computed.
>
> I am happy to discuss the formulas with you or anyone else if you have 
> any questions -- maybe more efficiently over bluejeans/zoom.
>
> Best regards,
> Wally
>
>
> ------------------------------------------------------------------------
> *From:* Jan C. Bernauer <jan.bernauer at stonybrook.edu>
> *Sent:* Monday, August 10, 2020 2:03 PM
> *To:* tdis at jlab.org <tdis at jlab.org>; Hobbs, Timothy 
> <TJHobbs at mail.smu.edu>; Wally Melnitchouk <wmelnitc at jlab.org>
> *Subject:* Re: [Tdis] [EXTERNAL] Inconsistencies in the TDIS proposal: 
> Some things I found
>
> Hi all,
>
> An update from the trenches of the generator writing:
>
> 1) G4SBS seems to agree on the definition of F2. Would be nice if 
> somebody who knows what they are doing (Tim, Wally?) can double check 
> though. So if there is no other effect, I'm pretty sure we see ~5 
> times more than we thought.
>
> 2) G4SBS also agrees with my estimate of the total DIS rate. It's 
> actually slightly higher, but I think that's because the simulated 
> e-arm acceptance is slightly bigger. I'll update my analysis with 
> that.  I'm pretty sure it's very close to 2 times what the proposal 
> has, so I'd hazard a guess and a factor of 2 comes from the number of 
> molecules to the number of protons.   It's also assuming that we run 
> at 77K --  I think we talked about running at room temperature, that 
> makes it 4 times worse of course.  (2 times worse than what is in). 
> Can somebody (Carlos?) run g4sbs in DIS configuration and see the 
> total rate between x 0.05 and 0.2? I assume a particle is accepted if 
> Harm.SBSGem.Track.P[0]>1  (for very few it's not track 0).
>
> So I think, if we run room temp, we'll get half the DIS rate we 
> assumed, but about 2.5 times the TDIS rate (but in a different 
> distribution, as far as I can tell). I'll make some updated reach 
> plots, but I don't have background numbers. Does anyone have those for 
> all bins?
>
> Best,
>
> Jan
>
>
> On 8/6/2020 11:43 AM, Jan C. Bernauer wrote:
>> Hi Tim, Hi Wally, Hi TDIS,
>>
>> Sorry for the lengthy email.
>>
>> Some updates and questions to my talk: 
>> https://hallaweb.jlab.org/wiki/images/2/2f/TDIS_sim.pdf 
>> <https://hallaweb.jlab.org/wiki/images/2/2f/TDIS_sim.pdf> where I 
>> found what I believe are internal inconsistencies, and differences to 
>> my own implementation of F_2, f_\pi, and a full MC.
>>
>> Wally, Tim, I also send this to you, because I believe you might have 
>> additional insight to the source of some of these figures/values.
>>
>>
>> 0) Does anybody have the exact bins used for the projected result 
>> figures? Especially the one as a function of t. That would be very 
>> helpful!
>>
>> 1) I do see about a factor of 2 more inclusive CS than in the 
>> proposal. I think it's likely that this error is on my end, I have to 
>> verify with g4SBS. In any case, the code version I got from Carlos 
>> (Thanks!), authored from Tim and Wally, do not calculate that as far 
>> as I can tell. It's not helping us in any case, we would just need 
>> less current I'd assume.
>>
>> 2) It would be great to figure out in which configuration the code 
>> was for the generation of table 6 and 7. I assume it was this code I 
>> attach here?  Is that the code version after a factor 2 has already 
>> been found? I believe so, because I can get the F_2^{\pi p} plots out 
>> that are in the new proposal, which are higher than the old proposal.
>>
>> Some things which I found in the code which do not match the experiment:
>>
>> Theta_e is ~12, not 35 degrees. That only affects the calculation of 
>> Q^2, which is slightly affecting the proton PDFs. Can't explain a big 
>> difference.
>>
>> The code was set up for pi^+, I'm looking at pi^0, so I changed the 
>> isospin factor  to 1. I changed some other integration ranges 
>> (ymax=1, xmax=0, km1, km2 ), and also implemented a cosph cut (proton 
>> theta<70)
>>
>> With these changes, I match exactly (on a log scale :) )  my own 
>> implementation  and the plots in the new proposal.
>>
>> 3) BUT the ratio was still off. I traced it down to  F_2^p. The 
>> relevant lines are here:
>>
>>          CALL SETCTQ6(1)  ! CTEQ 'MS-bar' SCHEME.
>>          u_pro = CTQ6PDF (1, REAL(x), SQRT(REAL(Q2)))
>>          ubar_pro = CTQ6PDF (-1, REAL(x), SQRT(REAL(Q2)))
>>          d_pro = CTQ6PDF (2, REAL(x), SQRT(REAL(Q2)))
>>          dbar_pro = CTQ6PDF (-2, REAL(x), SQRT(REAL(Q2)))
>>
>>          F2neu = 2.*x * ((4./9.)*(d_pro + dbar_pro)
>>      &             + (1./9.)*(u_pro + ubar_pro))
>>
>>
>>  a) LHAPDF 6 and this code has a different definition for what PID=1 
>> and 2 are. That stumped me for a while, but I'm pretty sure u_pro is 
>> indeed the proton u PDF.
>>
>> b) For the neutron, d_neu=u_pro and vice versa, so the line actually 
>> reads 2 *x * ( 4/9 ( u_neu +ubar_neu) + 2/9 ( d_neu +dbar_neu)). I 
>> changed it back to be correct for the proton.
>>
>> c) BUT: I do not have the 2 there. Where does that come from? I'm not 
>> super versed in PDFs, but it is my understanding that that shouldn't 
>> be there.  If not, can anybody please explain?
>>
>> ( d) There is also a small difference in the code here and LHAPDF for 
>> x<0.1.  10% or so. That must be in the underlying PDF or Q^2 evolution)
>>
>> In any case, without the 2, my code and this code essentially agree 
>> on F_2^p (with Q^2=1, very close to the plots in the proposal. With 
>> Q^2 changing, slightly different from the proposal, maybe 30% at 
>> most, but both codes the same way. Makes me think that the proposal 
>> line comes from a different program, which might explain the 
>> discrepancy with the 2).
>>
>>
>> 4) Going back to table 6. We already know that it was not updated 
>> from the first proposal, so it likely already has a factor of 2 
>> missing in the F/F ratio. With this additional factor 2, we are 
>> getting very close to what my program has, 4 is close enough to 5 
>> that I would believe the rest is acceptance, slightly different cuts, 
>> etc, or the first factor 2 was actually 2.5 or something.
>>
>> With these changes, both my code as well as the code I got from 
>> Carlos, modified as described above, gives a ratio F/F of 550 for the 
>> first line in the text. This assumes we accept k between 60 and 500 
>> Mev, x between 0.05 and 0.2, at around 12 degrees, with 30 to 70 deg 
>> proton angle.   This is also the number my MC gets, and roughly what 
>> I would get looking at the plots.  Or is there some other cut that 
>> should be applied? Cut on z (y in the code)?
>>
>> TLDR: If all my assumptions are correct, we see indeed 5.5 times more 
>> TDIS events (per DIS event) than we thought!
>>
>> Let me know what you think!
>>
>> Best,
>>
>> Jan
>>
>>
>> Attached: TDIS_orig.f, code I got from Carlos. TDIS.f: Code with my 
>> modifications.
>>
>>
>>
>>
>>
>> _______________________________________________
>> Tdis mailing list
>> Tdis at jlab.org  <mailto:Tdis at jlab.org>
>> https://mailman.jlab.org/mailman/listinfo/tdis
> -- 
> Dr. Jan C. Bernauer
> Assistant Professor
> Department of Physics and Astronomy
> Stony Brook University
> Stony Brook, NY 11794-3800

-- 
Dr. Jan C. Bernauer
Assistant Professor
Department of Physics and Astronomy
Stony Brook University
Stony Brook, NY 11794-3800

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://mailman.jlab.org/pipermail/tdis/attachments/20200811/d03021df/attachment-0001.html>