anchor.diff:    real(kind=fPrec), save :: xmax_pencil, ymax_pencil, xmax_nom, ymax_nom, nom_aperture, pencil_aperture
   
     real(kind=fPrec), save :: x_pencil0, y_pencil0, sum, sqsum
  -  real(kind=fPrec), save :: average, sigma, sigsecut, nspxd, xndisp, zpj
  -
  -  real(kind=fPrec), save :: max_tmp, a_tmp1, a_tmp2, ldrift, c_nex2, c_ney2, Nap1pos,Nap2pos,Nap1neg,Nap2neg
  -  real(kind=fPrec), save :: tiltOffsPos1,tiltOffsPos2,tiltOffsNeg1,tiltOffsNeg2
  -  real(kind=fPrec), save :: beamsize1, beamsize2,betax1,betax2,betay1,betay2, alphax1, alphax2,alphay1,alphay2,minAmpl
  +  real(kind=fPrec), save :: nspxd, xndisp, zpj
   
   #ifdef HDF5
     ! Variables to save hdf5 dataset indices
  
positive.code:   use geant4
  #endif
  
    implicit none
  
    integer :: j
  
  #ifdef HDF5
    ! For tracks2
  
negative.code:   if (n_impact.gt.0) then
      average = sum/n_impact
  
      if(sqsum/n_impact >= average**2) then
        sigma = sqrt(sqsum/n_impact - average**2)
      else
        sigma = zero
  