anchor.diff:                  print*, "  varname: "//trim(options%nsq_calibration_var)
                   call io_read2d(options%nsq_calibration_file,options%nsq_calibration_var,domain%nsq_calibration)
                   nsq_calibration=domain%nsq_calibration
  +                
  +                where(nsq_calibration<1) nsq_calibration = 1+1/((1-1/nsq_calibration)/100)
  +                where(nsq_calibration>1) nsq_calibration = 1+(nsq_calibration-1)/100
               endif
           endif
           
  
positive.code:                     top   = min(j+winsz,nz)
                      south = max(k-winsz,1)
                      north = min(k+winsz,ny)
                      
                      curnsq = sum(log(domain%nsquared(east:west, bottom:top,south:north)))
                      n = (((west-east)+1) * ((top-bottom)+1) * ((north-south)+1))
                      curnsq = curnsq / n
  
                      npos=1
  
negative.code:     type(C_PTR) :: uh_aligned_data, u_h_aligned_data, vh_aligned_data, v_h_aligned_data
      type(C_PTR), allocatable :: uplans(:), vplans(:)
      
      real, allocatable, dimension(:) :: dir_values, nsq_values, spd_values
      ! Look Up Tables for linear perturbation are nspd x n_dir_values x n_nsq_values x nx x nz x ny
      real, allocatable, target, dimension(:,:,:,:,:,:) :: hi_u_LUT, hi_v_LUT, rev_u_LUT, rev_v_LUT
  