anchor.diff:    real(RP) :: dataout9(9, 10)  ! results from Velocity Gradient
     real(RP) :: dataout18(18, 10) ! results from Velocity Hessian
   
  +  integer,parameter :: x=0, y=0, z=0, xwidth=16, ywidth=16, zwidth=16
  +  real(RP),dimension(xwidth*ywidth*zwidth*3) :: rawvelocity 
   
     ! Declare the return type of the turblib functions as integer.
     ! This is required for custom error handling (see the README).
     integer :: getvelocity, getforce, getvelocityandpressure, getvelocitygradient
     integer :: getboxfilter, getboxfiltersgs, getboxfiltergradient
     integer :: getvelocitylaplacian, getvelocityhessian
  -  integer :: getpressuregradient, getpressurehessian
  +  integer :: getpressure, getpressuregradient, getpressurehessian
  +  integer :: getrawvelocity, getrawpressure
     integer :: getposition
     ! If working with hdf5 cutout files, uncomment the lines below
     ! to load the file and use the functions locally:
  
positive.code:   end do
  
    write(*,*)
    write(*,'(a)') 'Requesting pressure gradient at 10 points...'
    rc = getpressuregradient(authkey, dataset,  time, FD4Lag4, NoTInt, 10, points, dataout3)
    do i = 1, 10 
  
negative.code:   integer :: i
  
    ! Formatting rules
    character(*), parameter :: format1='(i3,1(a,e13.6))'
    character(*), parameter :: format3='(i3,3(a,e13.6))'
    character(*), parameter :: format4='(i3,4(a,e13.6))'
    character(*), parameter :: format6='(i3,6(a,e13.6))'
  