anchor.diff:        write(lerr,"(a)")    "ELENS>       Please choose a value beterrn 1 and 20 (included)"
         iErr = .true.
         return
  +    end if
  +    if(st_debug) then
  +      call sixin_echoVal("interp. points",elens_radial_mpoints(melens),"ELENS",iLine)
       end if
        
     case default
  
positive.code:     ! 2) calculate radius
      !    radial position of main beam relative to center of elens beam
      !    rr = sqrt(xx**2+yy**2)
      rr=sqrt(xx*xx+yy*yy)
      ! 3) calculate kick
      !    shape function: spatial charge density depends on type:
      !    0    if r < R1
      !    frr  if R1 < r < R2
      !    1    if r > R2
      if (rr.gt.elens_r1(ielens(ix))) then ! rr <= r1 -> no kick from elens
        if (rr.lt.elens_r2(ielens(ix))) then ! r1 < rr < r2
          select case (elens_type(ielens(ix)))
          case (1)
            ! UNIFORM: eLens with uniform profile
            ! formula: (r^2-r1^2)/(r2^2-r1^2)
            frr=( (rr+elens_r1(ielens(ix)))*(rr-elens_r1(ielens(ix))) )/elens_geo_norm(ielens(ix))
          case (2)
            ! GAUSSIAN: eLens with Gaussian profile
            ! formula: (exp(-r1^2/2sig^2)-exp(-r^2/2sig^2))/(exp(-r1^2/2sig^2)-exp(-r2^2/2sig^2))
            frr=( exp_mb(-0.5*(elens_r1(ielens(ix))/elens_sig(ielens(ix)))**2)    &
                 -exp_mb(-0.5*(rr                  /elens_sig(ielens(ix)))**2) )/ &
                                    elens_geo_norm(ielens(ix))
          case (3)
            ! RADIAL PROFILE: eLens with radial profile as from file
            ! formula: (cumul_J(r)-cumul_J(r1))/(cumul_J(r2)-cumul_J(r1))
            frr=(polinterp( rr, &
                  elens_radial_profile_R(0:elens_radial_profile_nPoints(elens_iRadial(ielens(ix))),elens_iRadial(ielens(ix))), &
                  elens_radial_profile_J(0:elens_radial_profile_nPoints(elens_iRadial(ielens(ix))),elens_iRadial(ielens(ix))), &
                  elens_radial_profile_nPoints(elens_iRadial(ielens(ix)))+1, &
                  elens_radial_mpoints(ielens(ix)), elens_radial_jguess(ielens(ix)) )-elens_radial_fr1(ielens(ix))) &
                  /elens_geo_norm(ielens(ix))
          case default
            write(lerr,"(a,i0,a)") "ELENS> ERROR elens_kick: elens_type=",elens_type(ielens(ix))," not recognized. "
            write(lerr,"(a)")      "ELENS>       Possible values for type are: 1, 2 and 3"
            call prror
          end select
        else ! r1 < r2 <= rr
          frr = one
        endif
        ! 'radial kick'
        frr = (((elens_theta_r2(ielens(ix))*elens_r2(ielens(ix)))/rr)*frr)*moidpsv(jj)
        if(elens_lThetaR2(ielens(ix))) then
          if(elens_I(ielens(ix)) < zero) then
            frr = frr*((rvv(jj)+elens_beta_e(ielens(ix))*beta0)/(one+elens_beta_e(ielens(ix))*beta0))
  
negative.code:       tmpflt=elens_r2(ielens(iElem))
        elens_r2(ielens(iElem)) = elens_r1(ielens(iElem))
        elens_r1(ielens(iElem)) = tmpflt
      end if
      if(elens_r2(ielens(iElem)) < zero) then
        write(lerr,"(a)") "ELENS> ERROR R2<0!"
        iErr = .true.
        return
      else if(elens_r2(ielens(iElem)) < pieni) then
        elens_r2(ielens(iElem))=zero
      end if
      if(elens_r1(ielens(iElem)) < zero) then
        write(lerr,"(a)") "ELENS> ERROR R1<0!"
        iErr = .true.
        return
      else if(elens_r1(ielens(iElem)) < pieni) then
        elens_r1(ielens(iElem))=zero
      end if
      if(elens_lThetaR2(ielens(iElem))) then
        if(elens_len(ielens(iElem)) <= zero) then
  