anchor.diff:      logical :: smooth_nsq
       logical :: using_blocked_flow
       real    :: N_squared
  +    integer :: nth_file=1
   
       !! unfortunately these have to be allocated every call because we could be calling on both the high-res
       !! domain and the low res domain (to "remove" the linear winds)
  
positive.code:                         sweight = calc_weight(spd_values, spos, nexts, curspd)
                          nweight = calc_weight(nsq_values, npos, nextn, curnsq)
  
                          ! perform linear interpolation between LUT values
                          if (k<=ny) then
                              wind_first =      nweight  * (dweight * u_LUT(spos, dpos,npos, i,j,k) + (1-dweight) * u_LUT(spos, nextd,npos, i,j,k))   &
  
negative.code:         endif
  
          if (reverse) print*, "WARNING using fixed nsq for linear wind removal: 3e-6"
          ! $omp parallel firstprivate(nx,nxu,ny,nyv,nz, reverse, vsmooth, winsz, using_blocked_flow), default(none), &
          ! $omp private(i,j,k,step, uk, vi, east, west, north, south, top, bottom, u1d, v1d), &
          ! $omp private(spos, dpos, npos, nexts,nextd, nextn,n, smoothz, u, v, blocked), &
          ! $omp private(wind_first, wind_second, curspd, curdir, curnsq, sweight,dweight, nweight), &
          ! $omp shared(domain, spd_values, dir_values, nsq_values, u_LUT, v_LUT, linear_mask), &
          ! $omp shared(u_perturbation, v_perturbation, linear_update_fraction, linear_contribution, nsq_calibration), &
          ! $omp shared(min_stability, max_stability, n_dir_values, n_spd_values, n_nsq_values, smooth_nsq)
          !
          ! $omp do
          do k=1,ny
              do j=1,nz
                  do i=1,nx
  