photoproduction in proton-Lead

Hi,

1) If you use the heavy ion mode (nb_protron !=1) then the parameter ebeam1 is the energy of the beam (not the energy of the proton/neutron within the heavy ion).
2) lpp=2 is really a photon from a proton. Not a photon from an heavy ion.

I'm not an expert in heavy ion, but from the syntax, I would have expect to provide the elastic component of the full ion.
(which should be possible to extend the current formula to cover).

Your second syntax correspond more to the elastic photon from a proton regardless is the heavy-ion brakes or not. (they might have overlap between the two, I guess). Then I would have expect to re-scale this by 82 (the number of proton) and not by 208 (what is that number?)...

In any case to avoid confusion, I will forbid to use lpp=2 with heavy ion mode for the moment.

Cheers,

Olivier

Dear Olivier,

Thanks for the reply.

Can the photon flux be modeled as in this paper?

https://journals.aps.org/prd/abstract/10.1103/PhysRevD.88.054025

so when it is lpp=2

we could get a photon from a heavy ion.

Which will be different to when the photon comes from the proton.

Another question is if it is correct to express a gamma-proton interaction in pPb as:

---------------------------------------------------------
---------------------------------------------------------
define myp = g u d c s b u~ d~ c~ s~ b~
define myt = t t~
define myw = w+ w-

a myp > myt myw
---------------------------------------------------------
---------------------------------------------------------

Or should one do something like

myp myp > t t~ (via photoproduction)

Which I wouldn't know how to express.

Any help would be appreciated.

Cheers

Javier

Hi,

Not sure what the paper does. Here is a patch to implement the component where the lead is kept intact:
=== modified file 'madgraph/various/banner.py'
--- madgraph/various/banner.py 2019-11-14 22:43:14 +0000
+++ madgraph/various/banner.py 2019-11-21 17:20:58 +0000
@@ -3109,12 +3109,12 @@
             self.get_default('lhaid', log_level=20)

         # if heavy ion mode use for one beam, forbis lpp!=1
- if self['lpp1'] !=1:
- if self['nb_proton1'] != 1 or self['nb_neutron1'] !=0:
- raise InvalidRunCard, "Heavy ion mode is only supported for lpp1=1"
- if self['lpp2'] !=1:
- if self['nb_proton2'] != 1 or self['nb_neutron2'] !=0:
- raise InvalidRunCard, "Heavy ion mode is only supported for lpp1=1"
+ if self['lpp1'] not in [1,2]:
+ if self['nb_proton1'] !=1 or self['nb_neutron1'] !=0:
+ raise InvalidRunCard, "Heavy ion mode is only supported for lpp1=1/2"
+ if self['lpp2'] not in [1,2]:
+ if self['nb_proton2'] !=1 or self['nb_neutron2'] !=0:
+ raise InvalidRunCard, "Heavy ion mode is only supported for lpp2=1/2"
=== modified file 'Template/LO/Source/PDF/PhotonFlux.f'
--- Template/LO/Source/PDF/PhotonFlux.f 2018-04-20 08:02:12 +0000
+++ Template/LO/Source/PDF/PhotonFlux.f 2019-11-21 17:40:59 +0000
@@ -38,19 +38,30 @@

       end

- real*8 function epa_proton(x,q2max)
+ real*8 function epa_proton(x,q2max,beamid)
       integer i
+ integer beamid
       real*8 x,phi_f
       real*8 xin
       real*8 alpha,qz
       real*8 f, qmi,qma, q2max
       real*8 PI
+
+ integer nb_proton(2), nb_neutron(2)
+ common/to_heavyion_pdg/ nb_proton, nb_neutron
+ double precision mass_ion(2)
+ common/to_heavyion_mass/mass_ion
+
       data PI/3.14159265358979323846/
-
       data xin/0.938/ ! proton mass in GeV

       alpha = .0072992701
       qz = 0.71
+
+ if (nb_proton(beamid).ne.1.or.nb_neutron(beamid).ne.0)then
+ xin = mass_ion(beamid)
+ alpha = alpha * nb_proton(beamid)
+ endif

 C // x = omega/E = (E-E')/E
       if (x.lt.1) then

Now I do not think that I will publish such patch seems it seems to me that it is associate to some phase-space integration isssue (you need to probe value of x which are very-very small.
But fixing that will need to assign a student/someone on this (and have a paper out of this)

Here is the script that I was using to test the features (even if it fails with zero result).

import model sm-no_b_mass
generate a a > e+ e-
output
launch
set lhc 204*14
set PbPb
set lpp 2

Note that you might need to play with the factorization scale for such process...

Cheers,

Olivier

> On 21 Nov 2019, at 01:02, Javier <email address hidden> wrote:
>
> Question #685493 on MadGraph5_aMC@NLO changed:
> https://answers.launchpad.net/mg5amcnlo/+question/685493
>
> Javier posted a new comment:
> Dear Olivier,
>
> Thanks for the reply.
>
> Can the photon flux be modeled as in this paper?
>
> https://journals.aps.org/prd/abstract/10.1103/PhysRevD.88.054025
>
> so when it is lpp=2
>
> we could get a photon from a heavy ion.
>
> Which will be different to when the photon comes from the proton.
>
> Another question is if it is correct to express a gamma-proton
> interaction in pPb as:
>
> ---------------------------------------------------------
> ---------------------------------------------------------
> define myp = g u d c s b u~ d~ c~ s~ b~
> define myt = t t~
> define myw = w+ w-
>
> a myp > myt myw
> ---------------------------------------------------------
> ---------------------------------------------------------
>
> Or should one do something like
>
> myp myp > t t~ (via photoproduction)
>
> Which I wouldn't know how to express.
>
> Any help would be appreciated.
>
> Cheers
>
> Javier
>
> --
> You received this question notification because you are an answer
> contact for MadGraph5_aMC@NLO.

Hello Olivier,

Thanks for the code,

I tried to implement in MG5_aMC_v2_6_7,

But I don't see the lines with "minus" sign preceding, the ones that are supped to be removed from banner.py in your patch.

Once modified banner.py and PhotonFlux.f, is there a general way to recompile all the changes at once?

With zero result, you mean it wont run successfully?

Thanks

Javier

HI,

Just do not apply any patch to banner.py then.

After the patch, you will need to regenerate an output and that code will the be compiled.

> With zero result, you mean it wont run successfully?

It did not in my case at least.
Two potential reasons:
- the phase-space integration is not optimal (certainly true)
- the lead break too easily

Cheers,

Olivier

> On 26 Nov 2019, at 02:03, Javier <email address hidden> wrote:
>
> Question #685493 on MadGraph5_aMC@NLO changed:
> https://answers.launchpad.net/mg5amcnlo/+question/685493
>
> Javier posted a new comment:
> Hello Olivier,
>
> Thanks for the code,
>
> I tried to implement in MG5_aMC_v2_6_7,
>
> But I don't see the lines with "minus" sign preceding, the ones that are
> supped to be removed from banner.py in your patch.
>
> Once modified banner.py and PhotonFlux.f, is there a general way to
> recompile all the changes at once?
>
> With zero result, you mean it wont run successfully?
>
> Thanks
>
> Javier
>
> --
> You received this question notification because you are an answer
> contact for MadGraph5_aMC@NLO.

I have push the patch in this branch:;
lp:~maddevelopers/mg5amcnlo/elastic_fromheavyion

I have updated the patch to fix some phase-issue and now I have a non zero result for the following script::

generate a a > t t~
output
launch
set PbPb
set lpp 2
set fixed_scale 500
set ebeam 208*1500

Now it is clear that the systematics module is not behaving in a physical way (but this was expected) and I have no clue if such result is correct or not (pointless to say that the scale choice for the breaking of the lead might be set too high/....)

Cheers,

Olivier

Hello Olivier,

Many Thanks,

I already made sure it also runs for me.

Just wanted to ask, for the case

set PbPb

How can we set set the energy of each Pb nuclei to 2560 = 2.560 GeV?

With the current set up it gives aa (312000.0 x 312000.0 GeV).

I tried doing

set ebeam 2560

but it makes it crash.

Many Thanks and will continue testing.

Javier

ebeam is the energy of the lead beam.

so setting the energy of each beam to 2560 means that you have an energy of
2560/208 per nucleon (around 10GeV) so this is quite low energy per nucleon.
Is it really what you want?

Cheers,

Olivier

> On 12 Dec 2019, at 20:42, Javier <email address hidden> wrote:
>
> 2560

Now you should try to remove this test:

      if (x.eq.0d0) then
         pdg2pdf=0d0
         return
      elseif (x.lt.0d0 .or. (x*nb_hadron).gt.1d0) then
         if (nb_hadron.eq.1.or.x.lt.0d0) then
            write (*,*) 'PDF not supported for Bjorken x ', x*nb_hadron
            open(unit=26,file='../../../error',status='unknown')
            write(26,*) 'Error: PDF not supported for Bjorken x ',x*nb_hadron
            stop 1
         else
            pdg2pdf = 0d0
            return
         endif
      endif

might not be valid for photo-production (and check if they are no similar check in the code.
If you have something working do not hesitate to make the push in the above version for other user to be able to use it later.

Cheers,

Olivier

> On 12 Dec 2019, at 21:00, Olivier Mattelaer <email address hidden> wrote:
>
> ebeam is the energy of the lead beam.
>
> so setting the energy of each beam to 2560 means that you have an energy of
> 2560/208 per nucleon (around 10GeV) so this is quite low energy per nucleon.
> Is it really what you want?
>
> Cheers,
>
> Olivier
>
>> On 12 Dec 2019, at 20:42, Javier <<email address hidden> <mailto:<email address hidden>>> wrote:
>>
>> 2560
>

Thanks Olivier,

So for example a lead-proton collision could be:

generate a myp > t t~
output
launch
set Pbp
set lpp1 2
set fixed_scale 500
set ebeam1 2560*208
set ebeam2 6500

This could in principle simulate photon from the lead.

I also made this working copy to run.

I still have to check the latest comment you sent on removing that bit code.

to cross-check for example:

I found the resulting cross-section is sensitive to the value of the "fixed_scale".

For this case of Pbp at ~ 8 TeV, I recover something close to the expected ttbar photoproduction cross-section for gamma-proton, by setting this scale to ~22.

Though if I want the cross-section for Pb-gamma it doesnt increase as expected, it keeps a similar value.

I'll keep checking, but please let me know your comments.

Cheers

Javier

Hello Olivier,

I'm using the following configuration:

import model sm
define myp = g u d c s b u~ d~ c~ s~ b~
define myt = t t~
define myw = w+ w-
set group_subprocesses False
#
generate a myp > myt myw
output zdir_singlet_photo_Pbp_gammap_02
launch
set nevents 1000
set Pbp
set lpp1 2
set fixed_scale 1500
set ebeam1 208*2560
set ebeam2 6500
#set ebeam1 208*1500
#set ebeam1 2560
#set ebeam2 6500

I'm wondering how is this 'fixed_scale' related with the lead breaking and why it affects the XS values and kinematic distributions in general?

Thanks

Javier

This is the cutoff parameter in this paper:
V.M.Budnev et al., Phys.Rep. 15C (1975) 181
(the formalism used when you set lpp1 to 2)

Cheers,

Olivier

> On 22 Apr 2020, at 09:28, Javier <email address hidden> wrote:
>
> Question #685493 on MadGraph5_aMC@NLO changed:
> https://answers.launchpad.net/mg5amcnlo/+question/685493
>
> Javier posted a new comment:
> Hello Olivier,
>
> I'm using the following configuration:
>
> import model sm
> define myp = g u d c s b u~ d~ c~ s~ b~
> define myt = t t~
> define myw = w+ w-
> set group_subprocesses False
> #
> generate a myp > myt myw
> output zdir_singlet_photo_Pbp_gammap_02
> launch
> set nevents 1000
> set Pbp
> set lpp1 2
> set fixed_scale 1500
> set ebeam1 208*2560
> set ebeam2 6500
> #set ebeam1 208*1500
> #set ebeam1 2560
> #set ebeam2 6500
>
>
> I'm wondering how is this 'fixed_scale' related with the lead breaking and why it affects the XS values and kinematic distributions in general?
>
> Thanks
>
> Javier
>
> --
> You received this question notification because you are an answer
> contact for MadGraph5_aMC@NLO.