anchor.diff:  
   c       write(*,*)'in hyhard tau',tau
   
  -	etaLj=gaus(rm,sig) !fm
  +	etaLj=qpgaus(rm,sig) !fm
   
   c       write(*,*)'in hyhard after gauss',tau
   	z=tau*sinh(etaLj)
  
positive.code: *     generate 'njet' number of hard parton-parton scatterings with PYTHIA 
        IMPLICIT DOUBLE PRECISION(A-H, O-Z)
        double precision npar0,nbco0,npart,nbcol,npart0
        double precision gaus
        
        INTEGER PYK,PYCHGE,PYCOMP
        CHARACTER beam*2,targ*2
        external pydata 
        external pyp,pyr,pyk,pyquen,shad1,gaus
        
        common /pyjets/ n,npad,k(4000,5),p(4000,5),v(4000,5)
        common /hyjets/ nhj,nhp,khj(150000,5),phj(150000,5),vhj(150000,5)
  
negative.code: 
  
  * function to generate gauss distribution
        double precision function qpgaus(x0,sig)
        IMPLICIT DOUBLE PRECISION(A-H, O-Z)
   41   u1=pyr(0) 
        u2=pyr(0)  
  