Merge lp:~maddevelopers/mg5amcnlo/boost_of_events into lp:~maddevelopers/mg5amcnlo/2.4.1

[Olivier Mattelaer]

This changes the boost performed
from the partonic center of mass frame to the frame requested in the run_card.
It turns out that the previous one was not correctly taken into account the mass of the beam,
which is particularly problematic if one beam is taken at rest (or his energy being of the order of his mass).

I will send you the validation plot that I have done by email.
Rik I put you as reviewer mainly to check if similar problem can occur at NLO.
(actually how much is supported the beam mass @NLO)

Cheers,

Olivier

PS: This is actually a very small review (but this touch a central piece of the code, so we need to be carefull)

[Olivier Mattelaer (olivier-mattelaer)]

Oops j'ai pusher cela dans 2.4.1 par erreur.

Please comment if you see any problem on this.

Cheers,

Olivier

=== modified file 'Template/LO/SubProcesses/unwgt.f'
--- Template/LO/SubProcesses/unwgt.f	2015-02-18 01:33:22 +0000
+++ Template/LO/SubProcesses/unwgt.f	2016-05-05 11:23:32 +0000
@@ -138,6 +138,38 @@
       nw = 0
       maxwgt = 0d0
       end
+C**************************************************************************
+C      HELPING ROUTINE FOR PERFORMING THE Z BOOST OF THE EVENTS
+C**************************************************************************
+      double precision function get_betaz(pin,pout)
+C     compute the boost for the requested transformation
+      implicit none
+      double precision pin(0:3), pout(0:3)
+      double precision denom
+
+      denom = pin(0)*pout(0) - pin(3)*pout(3)
+      if (denom.ne.0d0) then
+         get_betaz = (pin(3) * pout(0) - pout(3) * pin(0)) / denom
+      else
+         get_betaz = 0d0
+      endif
+      return
+      end
+
+      subroutine zboost_with_beta(pin, beta, pout)
+c     apply the boost
+      implicit none
+      double precision pin(0:3), pout(0:3)
+      double precision beta, gamma
+
+      gamma = 1d0/DSQRT(1-beta**2)
+      pout(0) = gamma * pin(0) - gamma * beta * pin(3)
+      pout(1) = pin(1)
+      pout(2) = pin(2)
+      pout(3) = - gamma * beta * pin(0) + gamma * pin(3)
+      return
+      end
+
 
       SUBROUTINE unwgt(px,wgt,numproc)
 C**************************************************************************
@@ -436,7 +468,9 @@
       integer ip, np, ic, nc, jpart(7,-nexternal+3:2*nexternal-3)
       integer ida(2),ito(-nexternal+3:nexternal),ns,nres,ires,icloop
       integer iseed
-      double precision pboost(0:3),pb(0:4,-nexternal+3:2*nexternal-3),eta
+      double precision pboost(0:3),pb(0:4,-nexternal+3:2*nexternal-3)
+      double precision beta, get_betaz
+      double precision ebi(0:3), ebo(0:3)
       double precision ptcltmp(nexternal), pdum(0:3)
 
       integer idup(nexternal,maxproc,maxsproc)
@@ -564,10 +598,21 @@
       if (nincoming.eq.2) then
          if (xbk(1) .gt. 0d0 .and. xbk(1) .lt. 1d0 .and.
      $       xbk(2) .gt. 0d0 .and. xbk(2) .lt. 1d0) then
-            eta = sqrt(xbk(1)*ebeam(1)/
-     $                 (xbk(2)*ebeam(2)))
-            pboost(0)=p(0,1)*(eta + 1d0/eta)
-            pboost(3)=p(0,1)*(eta - 1d0/eta)
+            pboost(0) = -1d0 ! special flag
+            ! construct the beam momenta in each frame and compute the related (z)boost
+            ebi(0) = p(0,1)/xbk(1)
+            ebi(1) = 0
+            ebi(2) = 0
+            ebi(3) = DSQRT(ebi(0)**2-m1**2)
+            ebo(0) = ebeam(1)
+            ebo(1) = 0
+            ebo(2) = 0
+            ebo(3) = DSQRT(ebo(0)**2-m1**2)
+            beta = get_betaz(ebi, ebo)
+c           eta = sqrt(xbk(1)*ebeam(1)/
+c     $                 (xbk(2)*ebeam(2)))
+c            pboost(0)=p(0,1)*(eta + 1d0/eta)
+c            pboost(3)=p(0,1)*(eta - 1d0/eta)
          else if (xbk(1) .gt. 0d0 .and. xbk(1) .lt. 1d0 .and.
      $       xbk(2) .eq. 1d0) then
             pboost(0)=xbk(1)*ebeam(1)+ebeam(2)
@@ -587,12 +632,18 @@
      $           xbk(1),xbk(2)
          endif
       endif
-      do j=1,nexternal
-         call boostx(p(0,j),pboost,pb(0,isym(j,jsym)))
-c      Add mass information in pb(4)
-         pb(4,isym(j,jsym))=pmass(j)
-      enddo
-
+      if (pboost(0).ne.-1d0) then
+         do j=1,nexternal
+            call boostx(p(0,j),pboost,pb(0,isym(j,jsym)))
+            ! Add mass information in pb(4)
+            pb(4,isym(j,jsym))=pmass(j)
+         enddo
+      else
+         do j=1,nexternal
+            call zboost_with_beta(p(0,j),beta,pb(0,isym(j,jsym)))
+            pb(4,isym(j,jsym))=pmass(j)
+         enddo
+      endif
 c
 c     Add info on resonant mothers
 c
 
=== modified file 'UpdateNotes.txt'
--- UpdateNotes.txt	2016-04-06 18:45:52 +0000
+++ UpdateNotes.txt	2016-05-05 11:23:32 +0000
@@ -1,6 +1,9 @@
 Update notes for MadGraph5_aMC@NLO (in reverse time order)
 
-
+2.4.1 (xx/xx/xx)
+        OM: Fix a bug in fix target experiment with PDF on the particle at rest.
+            The cross-section was correct but the z-boost was not performed correctly.
+ 
 2.4.0.beta (06/04/16)
         OM: Allowing the proper NLO reweighting for NLO sample
 	RF: For NLO processes allow for multiple PDF and scales reweighting, directy by inputting lists
 
=== modified file 'VERSION'
--- VERSION	2016-04-06 18:45:52 +0000
+++ VERSION	2016-05-05 11:23:32 +0000
@@ -1,4 +1,4 @@
-version = 2.4.0.beta
+version = 2.4.1.beta
 date = 2016-04-06
 
 
 
=== modified file 'madgraph/interface/common_run_interface.py'
--- madgraph/interface/common_run_interface.py	2016-04-15 11:14:25 +0000
+++ madgraph/interface/common_run_interface.py	2016-05-05 11:23:32 +0000
@@ -2269,19 +2269,21 @@
 
         pdfsetname=set()
         for lhaid in lhaid_list:
-            try:
-                if lhaid in self.lhapdf_pdfsets:
-                    pdfsetname.add(self.lhapdf_pdfsets[lhaid]['filename'])
-                else:
-                    raise MadGraph5Error('lhaid %s not valid input number for the current lhapdf' % lhaid )
-            except KeyError:
-                if self.lhapdf_version.startswith('5'):
-                    raise MadGraph5Error(\
-                        ('invalid lhaid set in th run_card: %d .\nPlease note that some sets' % lhaid) + \
-                         '(eg MSTW 90%CL error sets) \nare not available in aMC@NLO + LHAPDF 5.x.x')
-                else:
-                    logger.debug('%d not found in pdfsets.index' % lhaid)
-
+            if isinstance(lhaid, (int,float)):
+                try:
+                    if lhaid in self.lhapdf_pdfsets:
+                        pdfsetname.add(self.lhapdf_pdfsets[lhaid]['filename'])
+                    else:
+                        raise MadGraph5Error('lhaid %s not valid input number for the current lhapdf' % lhaid )
+                except KeyError:
+                    if self.lhapdf_version.startswith('5'):
+                        raise MadGraph5Error(\
+                            ('invalid lhaid set in th run_card: %d .\nPlease note that some sets' % lhaid) + \
+                             '(eg MSTW 90%CL error sets) \nare not available in aMC@NLO + LHAPDF 5.x.x')
+                    else:
+                        logger.debug('%d not found in pdfsets.index' % lhaid)
+            else:
+                pdfsetname.add(lhaid)
 
         # check if the file exists, otherwise install it:
         # also check that the PDFsets dir exists, otherwise create it.