From 285158a18b31d00fdc56370bc824213206706909 Mon Sep 17 00:00:00 2001
From: "Pedro A. Jimenez" <jimenez@ucar.edu>
Date: Tue, 27 Jan 2026 15:48:25 -0700
Subject: [PATCH 1/8] Add horizontal interpolation option (nearest neighbor or
 bilinear)

---
 io/namelist_mod.F90  |  5 +++
 io/wrf_mod.F90       | 45 ++++++++++++-------------
 share/interp_mod.F90 | 80 +++++++++++++++++++++++++++++++++++++++++++-
 state/state_mod.F90  | 30 ++++++++---------
 4 files changed, 120 insertions(+), 40 deletions(-)

diff --git a/io/namelist_mod.F90 b/io/namelist_mod.F90
index b91d2b1..cfd628b 100644
--- a/io/namelist_mod.F90
+++ b/io/namelist_mod.F90
@@ -50,6 +50,7 @@ module namelist_mod
 
       real :: fire_wind_height = 6.096        ! "height of uah,vah wind in fire spread formula" "m"
       integer :: wind_vinterp_opt = 0         ! "wind (adjustment factor) interpolation option"
+      integer :: hinterp_opt = 1              ! "Horizontal interpolation from atm to fire (offline option): 1) ngp, 2)bi-linear"
       logical :: fire_lsm_zcoupling = .false. ! "flag to activate reference velocity at a different height from fire_wind_height"
       real :: fire_lsm_zcoupling_ref = 50.0   ! "reference height from wich u at fire_wind_hegiht is calculated using a logarithmic profile" "m"
 
@@ -191,6 +192,7 @@ subroutine Broadcast_nml (this)
       call Broadcast_integer (this%ros_opt)
       call Broadcast_integer (this%emis_opt)
       call Broadcast_integer (this%wind_vinterp_opt)
+      call Broadcast_integer (this%hinterp_opt)
 
       call Broadcast_integer (this%fire_num_ignitions)
 
@@ -410,6 +412,7 @@ subroutine Init_fire_block (this, file_name)
       real :: fmoist_dt = 600                 ! "moisture model time step" "s"
       real :: fire_wind_height = 6.096        ! "height of uah,vah wind in fire spread formula" "m"
       integer :: wind_vinterp_opt = 0         ! "wind (adjustment factor) interpolation option"
+      integer :: hinterp_opt = 1              ! "Horizontal interpolation from atm to fire (offline option): 1) ngp, 2) bi-linear"
       logical :: fire_is_real_perim = .false. ! .false. = point/line ignition, .true. = observed perimeter"
       real :: frac_fburnt_to_smoke = 0.02     ! "parts per unit of burned fuel becoming smoke " "g_smoke/kg_air"
       real :: fuelmc_g = 0.08                 ! Fuel moisture content ground (Dead FMC)
@@ -463,6 +466,7 @@ subroutine Init_fire_block (this, file_name)
             ! objects
           fuel_opt, ros_opt, fmc_opt, emis_opt, &
           wind_vinterp_opt, &
+          hinterp_opt, &
             ! Ignitions
           fire_num_ignitions, &
             ! Ignition 1
@@ -527,6 +531,7 @@ subroutine Init_fire_block (this, file_name)
       this%fmc_opt = fmc_opt
       this%emis_opt = emis_opt
       this%wind_vinterp_opt = wind_vinterp_opt
+      this%hinterp_opt = hinterp_opt
 
       this%fire_num_ignitions = fire_num_ignitions
 
diff --git a/io/wrf_mod.F90 b/io/wrf_mod.F90
index 757db4e..311a3c4 100644
--- a/io/wrf_mod.F90
+++ b/io/wrf_mod.F90
@@ -6,7 +6,8 @@ module wrf_mod
     use netcdf_mod, only : Get_netcdf_var, Get_netcdf_att, Get_netcdf_dim, Is_netcdf_file_present
     use proj_lc_mod, only : proj_lc_t
     use stderrout_mod, only : Print_message, Stop_simulation
-    use interp_mod, only : Interp_profile
+    use interp_mod, only : Interp_profile, Interp_horizontal_nearest, Interp_horizontal_bilinear, &
+        HINTERP_NEAREST, HINTERP_BILINEAR
 
     implicit none
 
@@ -47,7 +48,7 @@ module wrf_mod
       procedure, public :: Get_v3d => Get_meridional_wind_3d
       procedure, public :: Get_v3d_stag => Get_meridional_wind_stag_3d
       procedure, public :: Get_z0 => Get_z0
-      procedure, public :: Interp_var2grid_nearest => Interp_var2grid_nearest
+      procedure, public :: Interp_var2grid => Interp_var2grid
       procedure, public :: Print_domain => Print_domain
       procedure, public :: Update_atm_state => Update_atm_state
     end type wrf_t
@@ -644,7 +645,7 @@ end subroutine Write_latlon_check
 
     end subroutine Get_latcloncs
 
-    subroutine Interp_var2grid_nearest (this, lats_in, lons_in, var_name, vals_out)
+    subroutine Interp_var2grid (this, lats_in, lons_in, var_name, hinterp_opt, vals_out)
 
       use, intrinsic :: iso_fortran_env, only : ERROR_UNIT
 
@@ -653,18 +654,16 @@ subroutine Interp_var2grid_nearest (this, lats_in, lons_in, var_name, vals_out)
       class (wrf_t), intent(in) :: this
       real, dimension(:, :), intent(in) :: lats_in, lons_in
       character (len = *), intent (in) :: var_name
+      integer, intent (in) :: hinterp_opt
       real, dimension(:, :), allocatable, intent(out) :: vals_out
 
       real, parameter :: DEFAULT_INIT = 0.0
-      real, dimension(:, :), allocatable :: ds, var_wrf
-      real :: d, i_real, j_real
+      real, dimension(:, :), allocatable :: var_wrf
       type (proj_lc_t) :: proj
-      integer :: nx, ny, nx_wrf, ny_wrf, i, j, i_wrf, j_wrf
+      integer :: nx, ny, i, j
 
 
         ! Init
-      nx_wrf = this%ide - 1
-      ny_wrf = this%jde - 1
       nx = size (lats_in, dim = 1)
       ny = size (lats_in, dim = 2)
 
@@ -672,11 +671,9 @@ subroutine Interp_var2grid_nearest (this, lats_in, lons_in, var_name, vals_out)
       allocate (vals_out(nx, ny))
       vals_out = DEFAULT_INIT
 
-      allocate (ds(nx, ny))
-      ds = huge (d)
-
       proj = this%Get_projection ()
 
+        ! Get WRF data
       select case (var_name)
         case ('t2')
           var_wrf = this%t2_stag(this%ids:this%ide - 1, this%jds:this%jde - 1)
@@ -705,20 +702,20 @@ subroutine Interp_var2grid_nearest (this, lats_in, lons_in, var_name, vals_out)
       end select
 
         ! Algorithm
-      do j = 1, ny
-        do i = 1, nx
-          call proj%Calc_ij (lats_in(i, j), lons_in(i, j), i_real, j_real)
-          i_wrf = min (max (1, nint (i_real)), nx_wrf)
-          j_wrf = min (max (1, nint (j_real)), ny_wrf)
-          d = (this%lats(i_wrf, j_wrf) - lats_in(i, j)) ** 2 + (this%lons(i_wrf, j_wrf) - lons_in(i, j)) ** 2
-          if (d < ds(i, j)) then
-            vals_out(i, j) = var_wrf(i_wrf, j_wrf)
-            ds(i, j) = d
-          end if
-        end do
-      end do
+      Hinterp: select case (hinterp_opt)
+
+        case (HINTERP_NEAREST)
+          call Interp_horizontal_nearest (var_wrf, proj, lats_in, lons_in, vals_out)
+
+        case (HINTERP_BILINEAR)
+          call Interp_horizontal_bilinear (var_wrf, proj, lats_in, lons_in, vals_out)
+
+        case default
+        call Stop_simulation ('The horizontal interpolation option selected does not exist.')
+
+      end select Hinterp
 
-    end subroutine Interp_var2grid_nearest
+    end subroutine Interp_var2grid
 
     subroutine Interp_wrf2dvar_to_cfbm (wrfatm2dvar, ims, ime, jms, jme, ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe, &
         lats_in, lons_in, proj, vals_out)
diff --git a/share/interp_mod.F90 b/share/interp_mod.F90
index 57932b3..0e8b183 100644
--- a/share/interp_mod.F90
+++ b/share/interp_mod.F90
@@ -1,13 +1,91 @@
   module interp_mod
 
+    use proj_lc_mod, only : proj_lc_t
+
     implicit none
 
     private
 
-    public :: Interp_profile
+    integer, parameter :: HINTERP_NEAREST = 1, HINTERP_BILINEAR = 2
+
+    public :: Interp_profile, Interp_horizontal_nearest, Interp_horizontal_bilinear, HINTERP_NEAREST, HINTERP_BILINEAR
 
   contains
 
+    subroutine Interp_horizontal_nearest (data_in, proj_data_in, lats_out, lons_out, data_out)
+
+    ! Purpose: Nearest neighbor interpolation/extrapolation
+
+      implicit none
+
+      real, dimension(:, :), intent (in) :: data_in
+      type (proj_lc_t), intent (in) :: proj_data_in
+      real, dimension (:, :), intent (in) :: lats_out, lons_out
+!      real, dimension (size (lats_out, dim = 1), size (lats_out, dim = 2)), intent (out) :: data_out
+      real, dimension (:, :), intent (in out) :: data_out
+
+      integer :: i, j, i_in, j_in, nx, ny, nx_in, ny_in
+      real :: i_real, j_real
+
+
+      nx = size (lats_out, dim = 1)
+      ny = size (lats_out, dim = 2)
+      nx_in = size (data_in, dim = 1)
+      ny_in = size (data_in, dim = 2)
+
+      do j = 1, ny
+        do i = 1, nx
+          call proj_data_in%Calc_ij (lats_out(i, j), lons_out(i, j), i_real, j_real)
+          i_in = min (max (1, nint (i_real)), nx_in)
+          j_in = min (max (1, nint (j_real)), ny_in)
+          data_out(i, j) = data_in(i_in, j_in)
+        end do
+      end do
+
+    end subroutine Interp_horizontal_nearest
+
+    subroutine Interp_horizontal_bilinear (data_in, proj_data_in, lats_out, lons_out, data_out)
+
+    ! Purpose: bi-linear interpolation + nearest neighbor extrapolation
+
+      implicit none
+
+      real, dimension(:, :), intent (in) :: data_in
+      type (proj_lc_t), intent (in) :: proj_data_in
+      real, dimension (:, :), intent (in) :: lats_out, lons_out
+!      real, dimension (size (lats_out, dim = 1), size (lats_out, dim = 2)), intent (out) :: data_out
+      real, dimension (:, :), intent (in out) :: data_out
+
+      integer :: i, j, i0, j0, i1, j1, nx, ny, nx_in, ny_in
+      real :: i_real, j_real, di, dj
+
+
+      nx = size (lats_out, dim = 1)
+      ny = size (lats_out, dim = 2)
+      nx_in = size (data_in, dim = 1)
+      ny_in = size (data_in, dim = 2)
+
+      do j = 1, ny
+        do i = 1, nx
+          call proj_data_in%Calc_ij (lats_out(i, j), lons_out(i, j), i_real, j_real)
+
+          i0 = max (1, min (nx_in - 1, int (floor (i_real))))
+          j0 = max (1, min (ny_in - 1, int (floor (j_real))))
+          i1 = i0 + 1
+          j1 = j0 + 1
+
+          di = max (0.0, min (1.0, i_real - real (i0)))
+          dj = max (0.0, min (1.0, j_real - real (j0)))
+
+          data_out(i, j) = (1.0 - di) * (1.0 - dj) * data_in(i0, j0) + &
+              di * (1.0 - dj) * data_in(i1, j0) + &
+              (1.0 - di) * dj * data_in(i0, j1) + &
+              di * dj * data_in(i1, j1)
+        end do
+      end do
+
+    end subroutine Interp_horizontal_bilinear
+
     subroutine Interp_profile (fire_lsm_zcoupling, fire_lsm_zcoupling_ref, fire_wind_height, kfds, kfde, &
         uin, vin, z_at_w, z0f, uout, vout)
 
diff --git a/state/state_mod.F90 b/state/state_mod.F90
index b928fc0..109b143 100644
--- a/state/state_mod.F90
+++ b/state/state_mod.F90
@@ -241,7 +241,7 @@ subroutine Handle_wrfdata_update (this, wrf, config_flags)
         call wrf%Update_atm_state (this%datetime_now)
 
         if (DEBUG_LOCAL) call Print_message ('  Interpolating WRF vars...')
-        call this%interpolate_vars_atm_to_fire(wrf, config_flags)
+        call this%Interpolate_vars_atm_to_fire(wrf, config_flags)
 
         call this%datetime_next_atm_update%Add_seconds (config_flags%interval_atm)
 
@@ -752,8 +752,8 @@ subroutine Interpolate_vars_atm_to_fire (this, wrf, config_flags)
       if (allocated (this%lats) .and. allocated (this%lons)) then
 
         If_start: if (this%datetime_now == this%datetime_start) then
-          call wrf%interp_var2grid_nearest (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
-              this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'fz0', var2d)
+          call wrf%Interp_var2grid (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
+              this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'fz0', config_flags%hinterp_opt, var2d)
               this%fz0(this%ifps:this%ifpe, this%jfps:this%jfpe) = var2d
         endif If_start
 
@@ -764,28 +764,28 @@ subroutine Interpolate_vars_atm_to_fire (this, wrf, config_flags)
           end do
         end do
 
-        call wrf%interp_var2grid_nearest (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
-            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'uf', var2d)
+        call wrf%Interp_var2grid (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
+            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'uf', config_flags%hinterp_opt, var2d)
             this%uf(this%ifps:this%ifpe, this%jfps:this%jfpe) = var2d
 
-        call wrf%interp_var2grid_nearest (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
-            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'vf', var2d)
+        call wrf%Interp_var2grid (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
+            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'vf', config_flags%hinterp_opt, var2d)
             this%vf(this%ifps:this%ifpe, this%jfps:this%jfpe) = var2d
 
-        call wrf%interp_var2grid_nearest (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
-            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 't2', var2d)
+        call wrf%Interp_var2grid (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
+            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 't2', config_flags%hinterp_opt, var2d)
             this%fire_t2(this%ifps:this%ifpe, this%jfps:this%jfpe) = var2d
 
-        call wrf%interp_var2grid_nearest (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
-            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'q2', var2d)
+        call wrf%Interp_var2grid (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
+            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'q2', config_flags%hinterp_opt, var2d)
             this%fire_q2(this%ifps:this%ifpe, this%jfps:this%jfpe) = var2d
 
-        call wrf%interp_var2grid_nearest (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
-            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'psfc', var2d)
+        call wrf%Interp_var2grid (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
+            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'psfc', config_flags%hinterp_opt, var2d)
             this%fire_psfc(this%ifps:this%ifpe, this%jfps:this%jfpe) = var2d
 
-        call wrf%interp_var2grid_nearest (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
-            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'rain', var2d)
+        call wrf%Interp_var2grid (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
+            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'rain', config_flags%hinterp_opt, var2d)
             this%fire_rain(this%ifps:this%ifpe, this%jfps:this%jfpe) = var2d
 
         deallocate (var2d)

From c897e2d64082f98400abb149c62474d0e721fe65 Mon Sep 17 00:00:00 2001
From: "Pedro A. Jimenez" <jimenez@ucar.edu>
Date: Tue, 27 Jan 2026 16:25:49 -0700
Subject: [PATCH 2/8] More elegant set up of horizontal interpolations

---
 io/wrf_mod.F90       | 28 +++++++++++-----------------
 share/interp_mod.F90 | 34 ++++++++++++++++------------------
 state/state_mod.F90  | 40 +++++++++++++++-------------------------
 3 files changed, 42 insertions(+), 60 deletions(-)

diff --git a/io/wrf_mod.F90 b/io/wrf_mod.F90
index 311a3c4..2b00c42 100644
--- a/io/wrf_mod.F90
+++ b/io/wrf_mod.F90
@@ -645,32 +645,25 @@ end subroutine Write_latlon_check
 
     end subroutine Get_latcloncs
 
-    subroutine Interp_var2grid (this, lats_in, lons_in, var_name, hinterp_opt, vals_out)
+    subroutine Interp_var2grid (this, lats_in, lons_in, ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe, &
+        var_name, hinterp_opt, vals_out)
 
       use, intrinsic :: iso_fortran_env, only : ERROR_UNIT
 
       implicit none
 
       class (wrf_t), intent(in) :: this
-      real, dimension(:, :), intent(in) :: lats_in, lons_in
+      integer, intent (in) :: ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe
+      real, dimension(ifms:ifme, jfms:jfme), intent(in) :: lats_in, lons_in
       character (len = *), intent (in) :: var_name
       integer, intent (in) :: hinterp_opt
-      real, dimension(:, :), allocatable, intent(out) :: vals_out
+      real, dimension(ifms:ifme, jfms:jfme), intent(in out) :: vals_out
 
-      real, parameter :: DEFAULT_INIT = 0.0
       real, dimension(:, :), allocatable :: var_wrf
       type (proj_lc_t) :: proj
-      integer :: nx, ny, i, j
 
 
         ! Init
-      nx = size (lats_in, dim = 1)
-      ny = size (lats_in, dim = 2)
-
-      if (allocated (vals_out)) deallocate (vals_out)
-      allocate (vals_out(nx, ny))
-      vals_out = DEFAULT_INIT
-
       proj = this%Get_projection ()
 
         ! Get WRF data
@@ -701,17 +694,18 @@ subroutine Interp_var2grid (this, lats_in, lons_in, var_name, hinterp_opt, vals_
           stop
       end select
 
-        ! Algorithm
+        ! Algorithm(s)
       Hinterp: select case (hinterp_opt)
-
         case (HINTERP_NEAREST)
-          call Interp_horizontal_nearest (var_wrf, proj, lats_in, lons_in, vals_out)
+          call Interp_horizontal_nearest (var_wrf, proj, ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe, &
+              lats_in, lons_in, vals_out)
 
         case (HINTERP_BILINEAR)
-          call Interp_horizontal_bilinear (var_wrf, proj, lats_in, lons_in, vals_out)
+          call Interp_horizontal_bilinear (var_wrf, proj, ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe, &
+              lats_in, lons_in, vals_out)
 
         case default
-        call Stop_simulation ('The horizontal interpolation option selected does not exist.')
+          call Stop_simulation ('The horizontal interpolation option selected does not exist.')
 
       end select Hinterp
 
diff --git a/share/interp_mod.F90 b/share/interp_mod.F90
index 0e8b183..35ec4da 100644
--- a/share/interp_mod.F90
+++ b/share/interp_mod.F90
@@ -12,7 +12,8 @@ module interp_mod
 
   contains
 
-    subroutine Interp_horizontal_nearest (data_in, proj_data_in, lats_out, lons_out, data_out)
+    subroutine Interp_horizontal_nearest (data_in, proj_data_in, ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe, &
+        lats_out, lons_out, data_out)
 
     ! Purpose: Nearest neighbor interpolation/extrapolation
 
@@ -20,21 +21,19 @@ subroutine Interp_horizontal_nearest (data_in, proj_data_in, lats_out, lons_out,
 
       real, dimension(:, :), intent (in) :: data_in
       type (proj_lc_t), intent (in) :: proj_data_in
-      real, dimension (:, :), intent (in) :: lats_out, lons_out
-!      real, dimension (size (lats_out, dim = 1), size (lats_out, dim = 2)), intent (out) :: data_out
-      real, dimension (:, :), intent (in out) :: data_out
+      integer, intent (in) :: ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe
+      real, dimension (ifms:ifme, jfms:jfme), intent (in) :: lats_out, lons_out
+      real, dimension (ifms:ifme, jfms:jfme), intent (in out) :: data_out
 
-      integer :: i, j, i_in, j_in, nx, ny, nx_in, ny_in
+      integer :: i, j, i_in, j_in, nx_in, ny_in
       real :: i_real, j_real
 
 
-      nx = size (lats_out, dim = 1)
-      ny = size (lats_out, dim = 2)
       nx_in = size (data_in, dim = 1)
       ny_in = size (data_in, dim = 2)
 
-      do j = 1, ny
-        do i = 1, nx
+      do j = jfps, jfpe
+        do i = ifps, ifpe
           call proj_data_in%Calc_ij (lats_out(i, j), lons_out(i, j), i_real, j_real)
           i_in = min (max (1, nint (i_real)), nx_in)
           j_in = min (max (1, nint (j_real)), ny_in)
@@ -44,7 +43,8 @@ subroutine Interp_horizontal_nearest (data_in, proj_data_in, lats_out, lons_out,
 
     end subroutine Interp_horizontal_nearest
 
-    subroutine Interp_horizontal_bilinear (data_in, proj_data_in, lats_out, lons_out, data_out)
+    subroutine Interp_horizontal_bilinear (data_in, proj_data_in, ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe, &
+        lats_out, lons_out, data_out)
 
     ! Purpose: bi-linear interpolation + nearest neighbor extrapolation
 
@@ -52,21 +52,19 @@ subroutine Interp_horizontal_bilinear (data_in, proj_data_in, lats_out, lons_out
 
       real, dimension(:, :), intent (in) :: data_in
       type (proj_lc_t), intent (in) :: proj_data_in
-      real, dimension (:, :), intent (in) :: lats_out, lons_out
-!      real, dimension (size (lats_out, dim = 1), size (lats_out, dim = 2)), intent (out) :: data_out
-      real, dimension (:, :), intent (in out) :: data_out
+      integer, intent (in) :: ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe
+      real, dimension (ifms:ifme, jfms:jfme), intent (in) :: lats_out, lons_out
+      real, dimension (ifms:ifme, jfms:jfme), intent (in out) :: data_out
 
-      integer :: i, j, i0, j0, i1, j1, nx, ny, nx_in, ny_in
+      integer :: i, j, i0, j0, i1, j1, nx_in, ny_in
       real :: i_real, j_real, di, dj
 
 
-      nx = size (lats_out, dim = 1)
-      ny = size (lats_out, dim = 2)
       nx_in = size (data_in, dim = 1)
       ny_in = size (data_in, dim = 2)
 
-      do j = 1, ny
-        do i = 1, nx
+      do j = jfps, jfpe
+        do i = ifps, ifpe
           call proj_data_in%Calc_ij (lats_out(i, j), lons_out(i, j), i_real, j_real)
 
           i0 = max (1, min (nx_in - 1, int (floor (i_real))))
diff --git a/state/state_mod.F90 b/state/state_mod.F90
index 109b143..dc18cd7 100644
--- a/state/state_mod.F90
+++ b/state/state_mod.F90
@@ -744,18 +744,15 @@ subroutine Interpolate_vars_atm_to_fire (this, wrf, config_flags)
       type (wrf_t), intent(inout) :: wrf
       type (namelist_t), intent (in) :: config_flags
 
-      real, dimension(:, :), allocatable :: var2d
       integer :: i, j
 
 
         ! We need the fire grid lat/lon
       if (allocated (this%lats) .and. allocated (this%lons)) then
 
-        If_start: if (this%datetime_now == this%datetime_start) then
-          call wrf%Interp_var2grid (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
-              this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'fz0', config_flags%hinterp_opt, var2d)
-              this%fz0(this%ifps:this%ifpe, this%jfps:this%jfpe) = var2d
-        endif If_start
+        if (this%datetime_now == this%datetime_start) call wrf%Interp_var2grid (this%lats, this%lons, &
+            this%ifms, this%ifme, this%jfms, this%jfme, this%ifps, this%ifpe, this%jfps, this%jfpe, 'fz0', &
+            config_flags%hinterp_opt, this%fz0)
 
         do j = 1, wrf%jde
           do i = 1, wrf%ide
@@ -764,31 +761,24 @@ subroutine Interpolate_vars_atm_to_fire (this, wrf, config_flags)
           end do
         end do
 
-        call wrf%Interp_var2grid (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
-            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'uf', config_flags%hinterp_opt, var2d)
-            this%uf(this%ifps:this%ifpe, this%jfps:this%jfpe) = var2d
+        call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
+            this%ifps, this%ifpe, this%jfps, this%jfpe, 'uf', config_flags%hinterp_opt, this%uf)
 
-        call wrf%Interp_var2grid (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
-            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'vf', config_flags%hinterp_opt, var2d)
-            this%vf(this%ifps:this%ifpe, this%jfps:this%jfpe) = var2d
+        call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
+            this%ifps, this%ifpe, this%jfps, this%jfpe, 'vf', config_flags%hinterp_opt, this%vf)
 
-        call wrf%Interp_var2grid (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
-            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 't2', config_flags%hinterp_opt, var2d)
-            this%fire_t2(this%ifps:this%ifpe, this%jfps:this%jfpe) = var2d
+        call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
+            this%ifps, this%ifpe, this%jfps, this%jfpe, 't2', config_flags%hinterp_opt, this%fire_t2)
 
-        call wrf%Interp_var2grid (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
-            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'q2', config_flags%hinterp_opt, var2d)
-            this%fire_q2(this%ifps:this%ifpe, this%jfps:this%jfpe) = var2d
+        call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
+            this%ifps, this%ifpe, this%jfps, this%jfpe, 'q2', config_flags%hinterp_opt, this%fire_q2)
 
-        call wrf%Interp_var2grid (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
-            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'psfc', config_flags%hinterp_opt, var2d)
-            this%fire_psfc(this%ifps:this%ifpe, this%jfps:this%jfpe) = var2d
+        call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
+            this%ifps, this%ifpe, this%jfps, this%jfpe, 'psfc', config_flags%hinterp_opt, this%fire_psfc)
 
-        call wrf%Interp_var2grid (this%lats(this%ifps:this%ifpe, this%jfps:this%jfpe), &
-            this%lons(this%ifps:this%ifpe, this%jfps:this%jfpe), 'rain', config_flags%hinterp_opt, var2d)
-            this%fire_rain(this%ifps:this%ifpe, this%jfps:this%jfpe) = var2d
+        call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
+            this%ifps, this%ifpe, this%jfps, this%jfpe, 'rain', config_flags%hinterp_opt, this%fire_rain)
 
-        deallocate (var2d)
       end if
 
     end subroutine Interpolate_vars_atm_to_fire

From 0181fbfd3681c89f6278118867f67e4d570a2436 Mon Sep 17 00:00:00 2001
From: "Pedro A. Jimenez" <jimenez@ucar.edu>
Date: Tue, 27 Jan 2026 17:12:44 -0700
Subject: [PATCH 3/8] Replace stop statement with stop_simulation sub

---
 io/wrf_mod.F90 | 6 ++----
 1 file changed, 2 insertions(+), 4 deletions(-)

diff --git a/io/wrf_mod.F90 b/io/wrf_mod.F90
index 2b00c42..9869f10 100644
--- a/io/wrf_mod.F90
+++ b/io/wrf_mod.F90
@@ -648,8 +648,6 @@ end subroutine Get_latcloncs
     subroutine Interp_var2grid (this, lats_in, lons_in, ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe, &
         var_name, hinterp_opt, vals_out)
 
-      use, intrinsic :: iso_fortran_env, only : ERROR_UNIT
-
       implicit none
 
       class (wrf_t), intent(in) :: this
@@ -690,8 +688,8 @@ subroutine Interp_var2grid (this, lats_in, lons_in, ifms, ifme, jfms, jfme, ifps
           var_wrf = this%va(this%ids:this%ide - 1, this%jds:this%jde - 1)
 
         case default
-          write (ERROR_UNIT, *) 'Unknown variable name to interpolate'
-          stop
+          call Stop_simulation ('Unknown variable name to interpolate')
+
       end select
 
         ! Algorithm(s)

From 9a8e6ac5717e87c6dae75ab0c77f3296cec2761b Mon Sep 17 00:00:00 2001
From: "Pedro A. Jimenez" <jimenez@ucar.edu>
Date: Tue, 27 Jan 2026 17:20:42 -0700
Subject: [PATCH 4/8] Better logic in offline interpolations

---
 state/state_mod.F90 | 46 ++++++++++++++++++++++-----------------------
 1 file changed, 22 insertions(+), 24 deletions(-)

diff --git a/state/state_mod.F90 b/state/state_mod.F90
index dc18cd7..3265e40 100644
--- a/state/state_mod.F90
+++ b/state/state_mod.F90
@@ -747,39 +747,37 @@ subroutine Interpolate_vars_atm_to_fire (this, wrf, config_flags)
       integer :: i, j
 
 
-        ! We need the fire grid lat/lon
-      if (allocated (this%lats) .and. allocated (this%lons)) then
+      if (.not. allocated (this%lats) .or. .not. allocated (this%lons)) &
+          call Stop_simulation ('Init lats/lons before calling hinterp atm variables')
 
-        if (this%datetime_now == this%datetime_start) call wrf%Interp_var2grid (this%lats, this%lons, &
-            this%ifms, this%ifme, this%jfms, this%jfme, this%ifps, this%ifpe, this%jfps, this%jfpe, 'fz0', &
-            config_flags%hinterp_opt, this%fz0)
+      if (this%datetime_now == this%datetime_start) call wrf%Interp_var2grid (this%lats, this%lons, &
+          this%ifms, this%ifme, this%jfms, this%jfme, this%ifps, this%ifpe, this%jfps, this%jfpe, 'fz0', &
+          config_flags%hinterp_opt, this%fz0)
 
-        do j = 1, wrf%jde
-          do i = 1, wrf%ide
-            call this%Interpolate_profile (config_flags, config_flags%fire_wind_height, this%kfds, this%kfde, &
-                wrf%u3d_stag(i,:,j),wrf%v3d_stag(i,:,j), wrf%phl_stag(i,:,j), wrf%ua(i,j),wrf%va(i,j),wrf%z0_stag(i,j))
-          end do
+      do j = 1, wrf%jde
+        do i = 1, wrf%ide
+          call this%Interpolate_profile (config_flags, config_flags%fire_wind_height, this%kfds, this%kfde, &
+              wrf%u3d_stag(i,:,j),wrf%v3d_stag(i,:,j), wrf%phl_stag(i,:,j), wrf%ua(i,j),wrf%va(i,j),wrf%z0_stag(i,j))
         end do
+      end do
 
-        call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
-            this%ifps, this%ifpe, this%jfps, this%jfpe, 'uf', config_flags%hinterp_opt, this%uf)
-
-        call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
-            this%ifps, this%ifpe, this%jfps, this%jfpe, 'vf', config_flags%hinterp_opt, this%vf)
+      call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
+          this%ifps, this%ifpe, this%jfps, this%jfpe, 'uf', config_flags%hinterp_opt, this%uf)
 
-        call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
-            this%ifps, this%ifpe, this%jfps, this%jfpe, 't2', config_flags%hinterp_opt, this%fire_t2)
+      call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
+          this%ifps, this%ifpe, this%jfps, this%jfpe, 'vf', config_flags%hinterp_opt, this%vf)
 
-        call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
-            this%ifps, this%ifpe, this%jfps, this%jfpe, 'q2', config_flags%hinterp_opt, this%fire_q2)
+      call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
+          this%ifps, this%ifpe, this%jfps, this%jfpe, 't2', config_flags%hinterp_opt, this%fire_t2)
 
-        call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
-            this%ifps, this%ifpe, this%jfps, this%jfpe, 'psfc', config_flags%hinterp_opt, this%fire_psfc)
+      call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
+          this%ifps, this%ifpe, this%jfps, this%jfpe, 'q2', config_flags%hinterp_opt, this%fire_q2)
 
-        call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
-            this%ifps, this%ifpe, this%jfps, this%jfpe, 'rain', config_flags%hinterp_opt, this%fire_rain)
+      call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
+          this%ifps, this%ifpe, this%jfps, this%jfpe, 'psfc', config_flags%hinterp_opt, this%fire_psfc)
 
-      end if
+      call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
+          this%ifps, this%ifpe, this%jfps, this%jfpe, 'rain', config_flags%hinterp_opt, this%fire_rain)
 
     end subroutine Interpolate_vars_atm_to_fire
 

From c924d40e58ffdf82157426525f11aa47c9490129 Mon Sep 17 00:00:00 2001
From: "Pedro A. Jimenez" <jimenez@ucar.edu>
Date: Tue, 27 Jan 2026 18:15:52 -0700
Subject: [PATCH 5/8] Parallel OpenMP horizontal interp atm (offline)

---
 io/wrf_mod.F90      | 34 +++++++++++++++++++++++++++-------
 state/state_mod.F90 | 22 ++++++++++++++--------
 2 files changed, 41 insertions(+), 15 deletions(-)

diff --git a/io/wrf_mod.F90 b/io/wrf_mod.F90
index 9869f10..f6dbe6b 100644
--- a/io/wrf_mod.F90
+++ b/io/wrf_mod.F90
@@ -645,13 +645,14 @@ end subroutine Write_latlon_check
 
     end subroutine Get_latcloncs
 
-    subroutine Interp_var2grid (this, lats_in, lons_in, ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe, &
-        var_name, hinterp_opt, vals_out)
+    subroutine Interp_var2grid (this, lats_in, lons_in, ifms, ifme, jfms, jfme, &
+        num_tiles, i_start, i_end, j_start, j_end, var_name, hinterp_opt, vals_out)
 
       implicit none
 
       class (wrf_t), intent(in) :: this
-      integer, intent (in) :: ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe
+      integer, intent (in) :: ifms, ifme, jfms, jfme, num_tiles
+      integer, dimension (num_tiles), intent (in) :: i_start, i_end, j_start, j_end
       real, dimension(ifms:ifme, jfms:jfme), intent(in) :: lats_in, lons_in
       character (len = *), intent (in) :: var_name
       integer, intent (in) :: hinterp_opt
@@ -659,6 +660,7 @@ subroutine Interp_var2grid (this, lats_in, lons_in, ifms, ifme, jfms, jfme, ifps
 
       real, dimension(:, :), allocatable :: var_wrf
       type (proj_lc_t) :: proj
+      integer :: ij, ifts, ifte, jfts, jfte
 
 
         ! Init
@@ -695,12 +697,30 @@ subroutine Interp_var2grid (this, lats_in, lons_in, ifms, ifme, jfms, jfme, ifps
         ! Algorithm(s)
       Hinterp: select case (hinterp_opt)
         case (HINTERP_NEAREST)
-          call Interp_horizontal_nearest (var_wrf, proj, ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe, &
-              lats_in, lons_in, vals_out)
+          !$OMP PARALLEL DO   &
+          !$OMP PRIVATE (ij, ifts, ifte, jfts, jfte)
+          do ij = 1, num_tiles
+            ifts = i_start(ij)
+            ifte = i_end(ij)
+            jfts = j_start(ij)
+            jfte = j_end(ij)
+            call Interp_horizontal_nearest (var_wrf, proj, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
+                lats_in, lons_in, vals_out)
+          end do
+          !$OMP END PARALLEL DO
 
         case (HINTERP_BILINEAR)
-          call Interp_horizontal_bilinear (var_wrf, proj, ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe, &
-              lats_in, lons_in, vals_out)
+          !$OMP PARALLEL DO   &
+          !$OMP PRIVATE (ij, ifts, ifte, jfts, jfte)
+          do ij = 1, num_tiles
+            ifts = i_start(ij)
+            ifte = i_end(ij)
+            jfts = j_start(ij)
+            jfte = j_end(ij)
+            call Interp_horizontal_bilinear (var_wrf, proj, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
+                lats_in, lons_in, vals_out)
+          end do
+          !$OMP END PARALLEL DO
 
         case default
           call Stop_simulation ('The horizontal interpolation option selected does not exist.')
diff --git a/state/state_mod.F90 b/state/state_mod.F90
index 3265e40..d73323e 100644
--- a/state/state_mod.F90
+++ b/state/state_mod.F90
@@ -751,8 +751,8 @@ subroutine Interpolate_vars_atm_to_fire (this, wrf, config_flags)
           call Stop_simulation ('Init lats/lons before calling hinterp atm variables')
 
       if (this%datetime_now == this%datetime_start) call wrf%Interp_var2grid (this%lats, this%lons, &
-          this%ifms, this%ifme, this%jfms, this%jfme, this%ifps, this%ifpe, this%jfps, this%jfpe, 'fz0', &
-          config_flags%hinterp_opt, this%fz0)
+          this%ifms, this%ifme, this%jfms, this%jfme, config_flags%num_tiles, this%i_start, this%i_end, &
+          this%j_start, this%j_end, 'fz0', config_flags%hinterp_opt, this%fz0)
 
       do j = 1, wrf%jde
         do i = 1, wrf%ide
@@ -762,22 +762,28 @@ subroutine Interpolate_vars_atm_to_fire (this, wrf, config_flags)
       end do
 
       call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
-          this%ifps, this%ifpe, this%jfps, this%jfpe, 'uf', config_flags%hinterp_opt, this%uf)
+          config_flags%num_tiles, this%i_start, this%i_end, this%j_start, this%j_end, &
+          'uf', config_flags%hinterp_opt, this%uf)
 
       call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
-          this%ifps, this%ifpe, this%jfps, this%jfpe, 'vf', config_flags%hinterp_opt, this%vf)
+          config_flags%num_tiles, this%i_start, this%i_end, this%j_start, this%j_end, &
+          'vf', config_flags%hinterp_opt, this%vf)
 
       call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
-          this%ifps, this%ifpe, this%jfps, this%jfpe, 't2', config_flags%hinterp_opt, this%fire_t2)
+          config_flags%num_tiles, this%i_start, this%i_end, this%j_start, this%j_end, &
+          't2', config_flags%hinterp_opt, this%fire_t2)
 
       call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
-          this%ifps, this%ifpe, this%jfps, this%jfpe, 'q2', config_flags%hinterp_opt, this%fire_q2)
+          config_flags%num_tiles, this%i_start, this%i_end, this%j_start, this%j_end, &
+          'q2', config_flags%hinterp_opt, this%fire_q2)
 
       call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
-          this%ifps, this%ifpe, this%jfps, this%jfpe, 'psfc', config_flags%hinterp_opt, this%fire_psfc)
+          config_flags%num_tiles, this%i_start, this%i_end, this%j_start, this%j_end, &
+          'psfc', config_flags%hinterp_opt, this%fire_psfc)
 
       call wrf%Interp_var2grid (this%lats, this%lons, this%ifms, this%ifme, this%jfms, this%jfme, &
-          this%ifps, this%ifpe, this%jfps, this%jfpe, 'rain', config_flags%hinterp_opt, this%fire_rain)
+          config_flags%num_tiles, this%i_start, this%i_end, this%j_start, this%j_end, &
+          'rain', config_flags%hinterp_opt, this%fire_rain)
 
     end subroutine Interpolate_vars_atm_to_fire
 

From 386df196ad1ba820be908451ebc563c2ef84582f Mon Sep 17 00:00:00 2001
From: "Pedro A. Jimenez" <jimenez@ucar.edu>
Date: Tue, 27 Jan 2026 19:19:23 -0700
Subject: [PATCH 6/8] Changing wrong patches indices to tile indices (no
 impact)

---
 share/interp_mod.F90 | 16 ++++++++--------
 1 file changed, 8 insertions(+), 8 deletions(-)

diff --git a/share/interp_mod.F90 b/share/interp_mod.F90
index 35ec4da..c3ef981 100644
--- a/share/interp_mod.F90
+++ b/share/interp_mod.F90
@@ -12,7 +12,7 @@ module interp_mod
 
   contains
 
-    subroutine Interp_horizontal_nearest (data_in, proj_data_in, ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe, &
+    subroutine Interp_horizontal_nearest (data_in, proj_data_in, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
         lats_out, lons_out, data_out)
 
     ! Purpose: Nearest neighbor interpolation/extrapolation
@@ -21,7 +21,7 @@ subroutine Interp_horizontal_nearest (data_in, proj_data_in, ifms, ifme, jfms, j
 
       real, dimension(:, :), intent (in) :: data_in
       type (proj_lc_t), intent (in) :: proj_data_in
-      integer, intent (in) :: ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe
+      integer, intent (in) :: ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte
       real, dimension (ifms:ifme, jfms:jfme), intent (in) :: lats_out, lons_out
       real, dimension (ifms:ifme, jfms:jfme), intent (in out) :: data_out
 
@@ -32,8 +32,8 @@ subroutine Interp_horizontal_nearest (data_in, proj_data_in, ifms, ifme, jfms, j
       nx_in = size (data_in, dim = 1)
       ny_in = size (data_in, dim = 2)
 
-      do j = jfps, jfpe
-        do i = ifps, ifpe
+      do j = jfts, jfte
+        do i = ifts, ifte
           call proj_data_in%Calc_ij (lats_out(i, j), lons_out(i, j), i_real, j_real)
           i_in = min (max (1, nint (i_real)), nx_in)
           j_in = min (max (1, nint (j_real)), ny_in)
@@ -43,7 +43,7 @@ subroutine Interp_horizontal_nearest (data_in, proj_data_in, ifms, ifme, jfms, j
 
     end subroutine Interp_horizontal_nearest
 
-    subroutine Interp_horizontal_bilinear (data_in, proj_data_in, ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe, &
+    subroutine Interp_horizontal_bilinear (data_in, proj_data_in, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
         lats_out, lons_out, data_out)
 
     ! Purpose: bi-linear interpolation + nearest neighbor extrapolation
@@ -52,7 +52,7 @@ subroutine Interp_horizontal_bilinear (data_in, proj_data_in, ifms, ifme, jfms,
 
       real, dimension(:, :), intent (in) :: data_in
       type (proj_lc_t), intent (in) :: proj_data_in
-      integer, intent (in) :: ifms, ifme, jfms, jfme, ifps, ifpe, jfps, jfpe
+      integer, intent (in) :: ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte
       real, dimension (ifms:ifme, jfms:jfme), intent (in) :: lats_out, lons_out
       real, dimension (ifms:ifme, jfms:jfme), intent (in out) :: data_out
 
@@ -63,8 +63,8 @@ subroutine Interp_horizontal_bilinear (data_in, proj_data_in, ifms, ifme, jfms,
       nx_in = size (data_in, dim = 1)
       ny_in = size (data_in, dim = 2)
 
-      do j = jfps, jfpe
-        do i = ifps, ifpe
+      do j = jfts, jfte
+        do i = ifts, ifte
           call proj_data_in%Calc_ij (lats_out(i, j), lons_out(i, j), i_real, j_real)
 
           i0 = max (1, min (nx_in - 1, int (floor (i_real))))

From 0dc061b0169a43a9cc7b420d11aec6a4d51cc788 Mon Sep 17 00:00:00 2001
From: "Pedro A. Jimenez" <jimenez@ucar.edu>
Date: Wed, 28 Jan 2026 06:18:10 -0700
Subject: [PATCH 7/8] passes atm indices to interpolation subroutines

---
 io/wrf_mod.F90       | 10 +++++++---
 share/interp_mod.F90 | 30 ++++++++++++------------------
 2 files changed, 19 insertions(+), 21 deletions(-)

diff --git a/io/wrf_mod.F90 b/io/wrf_mod.F90
index f6dbe6b..8a08ef6 100644
--- a/io/wrf_mod.F90
+++ b/io/wrf_mod.F90
@@ -660,7 +660,7 @@ subroutine Interp_var2grid (this, lats_in, lons_in, ifms, ifme, jfms, jfme, &
 
       real, dimension(:, :), allocatable :: var_wrf
       type (proj_lc_t) :: proj
-      integer :: ij, ifts, ifte, jfts, jfte
+      integer :: ij, ifts, ifte, jfts, jfte, ims, ime, jms, jme
 
 
         ! Init
@@ -694,6 +694,10 @@ subroutine Interp_var2grid (this, lats_in, lons_in, ifms, ifme, jfms, jfme, &
 
       end select
 
+      ims = 1
+      ime = size (var_wrf, dim = 1)
+      jms = 1
+      jme = size (var_wrf, dim = 2)
         ! Algorithm(s)
       Hinterp: select case (hinterp_opt)
         case (HINTERP_NEAREST)
@@ -704,7 +708,7 @@ subroutine Interp_var2grid (this, lats_in, lons_in, ifms, ifme, jfms, jfme, &
             ifte = i_end(ij)
             jfts = j_start(ij)
             jfte = j_end(ij)
-            call Interp_horizontal_nearest (var_wrf, proj, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
+            call Interp_horizontal_nearest (var_wrf, proj, ims, ime, jms, jme, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
                 lats_in, lons_in, vals_out)
           end do
           !$OMP END PARALLEL DO
@@ -717,7 +721,7 @@ subroutine Interp_var2grid (this, lats_in, lons_in, ifms, ifme, jfms, jfme, &
             ifte = i_end(ij)
             jfts = j_start(ij)
             jfte = j_end(ij)
-            call Interp_horizontal_bilinear (var_wrf, proj, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
+            call Interp_horizontal_bilinear (var_wrf, proj, ims, ime, jms, jme, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
                 lats_in, lons_in, vals_out)
           end do
           !$OMP END PARALLEL DO
diff --git a/share/interp_mod.F90 b/share/interp_mod.F90
index c3ef981..f653d50 100644
--- a/share/interp_mod.F90
+++ b/share/interp_mod.F90
@@ -12,63 +12,57 @@ module interp_mod
 
   contains
 
-    subroutine Interp_horizontal_nearest (data_in, proj_data_in, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
+    subroutine Interp_horizontal_nearest (data_in, proj_data_in, ims, ime, jms, jme, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
         lats_out, lons_out, data_out)
 
     ! Purpose: Nearest neighbor interpolation/extrapolation
 
       implicit none
 
-      real, dimension(:, :), intent (in) :: data_in
       type (proj_lc_t), intent (in) :: proj_data_in
-      integer, intent (in) :: ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte
+      integer, intent (in) :: ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, ims, ime, jms, jme
+      real, dimension(ims:ime, jms:jme), intent (in) :: data_in
       real, dimension (ifms:ifme, jfms:jfme), intent (in) :: lats_out, lons_out
       real, dimension (ifms:ifme, jfms:jfme), intent (in out) :: data_out
 
-      integer :: i, j, i_in, j_in, nx_in, ny_in
+      integer :: i, j, i_in, j_in
       real :: i_real, j_real
 
 
-      nx_in = size (data_in, dim = 1)
-      ny_in = size (data_in, dim = 2)
-
       do j = jfts, jfte
         do i = ifts, ifte
           call proj_data_in%Calc_ij (lats_out(i, j), lons_out(i, j), i_real, j_real)
-          i_in = min (max (1, nint (i_real)), nx_in)
-          j_in = min (max (1, nint (j_real)), ny_in)
+          i_in = min (max (ims, nint (i_real)), ime)
+          j_in = min (max (jms, nint (j_real)), ime)
           data_out(i, j) = data_in(i_in, j_in)
         end do
       end do
 
     end subroutine Interp_horizontal_nearest
 
-    subroutine Interp_horizontal_bilinear (data_in, proj_data_in, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
+    subroutine Interp_horizontal_bilinear (data_in, proj_data_in, ims, ime, jms, jme, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
         lats_out, lons_out, data_out)
 
     ! Purpose: bi-linear interpolation + nearest neighbor extrapolation
 
       implicit none
 
-      real, dimension(:, :), intent (in) :: data_in
       type (proj_lc_t), intent (in) :: proj_data_in
-      integer, intent (in) :: ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte
+      integer, intent (in) :: ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, ims, ime, jms, jme
+      real, dimension(ims:ime, jms:jme), intent (in) :: data_in
       real, dimension (ifms:ifme, jfms:jfme), intent (in) :: lats_out, lons_out
       real, dimension (ifms:ifme, jfms:jfme), intent (in out) :: data_out
 
-      integer :: i, j, i0, j0, i1, j1, nx_in, ny_in
+      integer :: i, j, i0, j0, i1, j1
       real :: i_real, j_real, di, dj
 
 
-      nx_in = size (data_in, dim = 1)
-      ny_in = size (data_in, dim = 2)
-
       do j = jfts, jfte
         do i = ifts, ifte
           call proj_data_in%Calc_ij (lats_out(i, j), lons_out(i, j), i_real, j_real)
 
-          i0 = max (1, min (nx_in - 1, int (floor (i_real))))
-          j0 = max (1, min (ny_in - 1, int (floor (j_real))))
+          i0 = max (ims, min (ime - 1, int (floor (i_real))))
+          j0 = max (jms, min (jme - 1, int (floor (j_real))))
           i1 = i0 + 1
           j1 = j0 + 1
 

From a3d22ca87deaba8c9f68b5ff85d8a7b83a5a857b Mon Sep 17 00:00:00 2001
From: "Pedro A. Jimenez" <jimenez@ucar.edu>
Date: Wed, 28 Jan 2026 07:12:50 -0700
Subject: [PATCH 8/8] Adding a sharable hinterp driver

---
 CMakeLists.txt      |  1 +
 io/coupling_mod.F90 | 65 +++++++++++++++++++++++++++++++++++++++++++++
 io/wrf_mod.F90      | 48 +++++++--------------------------
 3 files changed, 76 insertions(+), 38 deletions(-)
 create mode 100644 io/coupling_mod.F90

diff --git a/CMakeLists.txt b/CMakeLists.txt
index 8b17730..818b257 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -74,6 +74,7 @@ list(APPEND _share_files share/proj_lc_mod.F90
 
 # IO
 list(APPEND _io_files io/namelist_mod.F90
+                      io/coupling_mod.F90
                       io/geogrid_mod.F90
                       io/wrf_mod.F90
 		      io/netcdf_mod.F90
diff --git a/io/coupling_mod.F90 b/io/coupling_mod.F90
new file mode 100644
index 0000000..bd912f4
--- /dev/null
+++ b/io/coupling_mod.F90
@@ -0,0 +1,65 @@
+  module coupling_mod
+
+    use proj_lc_mod, only : proj_lc_t
+    use interp_mod, only : Interp_horizontal_nearest, Interp_horizontal_bilinear, HINTERP_NEAREST, HINTERP_BILINEAR
+    use stderrout_mod, only : Stop_simulation
+
+
+    implicit none
+
+    private
+
+    public :: Interp_horizontal
+
+  contains
+
+    subroutine Interp_horizontal (data_in, proj_data_in, ims, ime, jms, jme, ifms, ifme, jfms, jfme, &
+        num_tiles, i_start, i_end, j_start, j_end, hinterp_opt, lats_out, lons_out, data_out)
+
+      implicit none
+
+      type (proj_lc_t), intent (in) :: proj_data_in
+      integer, intent (in) :: ifms, ifme, jfms, jfme, ims, ime, jms, jme, num_tiles, hinterp_opt
+      integer, dimension (num_tiles), intent (in) :: i_start, i_end, j_start, j_end
+      real, dimension(ims:ime, jms:jme), intent (in) :: data_in
+      real, dimension (ifms:ifme, jfms:jfme), intent (in) :: lats_out, lons_out
+      real, dimension (ifms:ifme, jfms:jfme), intent (in out) :: data_out
+
+      integer :: ij, ifts, ifte, jfts, jfte
+
+
+      Hinterp: select case (hinterp_opt)
+        case (HINTERP_NEAREST)
+          !$OMP PARALLEL DO   &
+          !$OMP PRIVATE (ij, ifts, ifte, jfts, jfte)
+          do ij = 1, num_tiles
+            ifts = i_start(ij)
+            ifte = i_end(ij)
+            jfts = j_start(ij)
+            jfte = j_end(ij)
+            call Interp_horizontal_nearest (data_in, proj_data_in, ims, ime, jms, jme, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
+                lats_out, lons_out, data_out)
+          end do
+          !$OMP END PARALLEL DO
+
+        case (HINTERP_BILINEAR)
+          !$OMP PARALLEL DO   &
+          !$OMP PRIVATE (ij, ifts, ifte, jfts, jfte)
+          do ij = 1, num_tiles
+            ifts = i_start(ij)
+            ifte = i_end(ij)
+            jfts = j_start(ij)
+            jfte = j_end(ij)
+            call Interp_horizontal_bilinear (data_in, proj_data_in, ims, ime, jms, jme, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
+                lats_out, lons_out, data_out)
+          end do
+          !$OMP END PARALLEL DO
+
+        case default
+          call Stop_simulation ('The horizontal interpolation option selected does not exist.')
+
+      end select Hinterp
+
+    end subroutine Interp_horizontal
+
+  end module coupling_mod
diff --git a/io/wrf_mod.F90 b/io/wrf_mod.F90
index 8a08ef6..05217d5 100644
--- a/io/wrf_mod.F90
+++ b/io/wrf_mod.F90
@@ -6,8 +6,8 @@ module wrf_mod
     use netcdf_mod, only : Get_netcdf_var, Get_netcdf_att, Get_netcdf_dim, Is_netcdf_file_present
     use proj_lc_mod, only : proj_lc_t
     use stderrout_mod, only : Print_message, Stop_simulation
-    use interp_mod, only : Interp_profile, Interp_horizontal_nearest, Interp_horizontal_bilinear, &
-        HINTERP_NEAREST, HINTERP_BILINEAR
+    use interp_mod, only : Interp_profile
+    use coupling_mod, only : Interp_horizontal
 
     implicit none
 
@@ -645,22 +645,22 @@ end subroutine Write_latlon_check
 
     end subroutine Get_latcloncs
 
-    subroutine Interp_var2grid (this, lats_in, lons_in, ifms, ifme, jfms, jfme, &
-        num_tiles, i_start, i_end, j_start, j_end, var_name, hinterp_opt, vals_out)
+    subroutine Interp_var2grid (this, lats_out, lons_out, ifms, ifme, jfms, jfme, &
+        num_tiles, i_start, i_end, j_start, j_end, var_name, hinterp_opt, data_out)
 
       implicit none
 
       class (wrf_t), intent(in) :: this
       integer, intent (in) :: ifms, ifme, jfms, jfme, num_tiles
       integer, dimension (num_tiles), intent (in) :: i_start, i_end, j_start, j_end
-      real, dimension(ifms:ifme, jfms:jfme), intent(in) :: lats_in, lons_in
+      real, dimension(ifms:ifme, jfms:jfme), intent(in) :: lats_out, lons_out
       character (len = *), intent (in) :: var_name
       integer, intent (in) :: hinterp_opt
-      real, dimension(ifms:ifme, jfms:jfme), intent(in out) :: vals_out
+      real, dimension(ifms:ifme, jfms:jfme), intent(in out) :: data_out
 
       real, dimension(:, :), allocatable :: var_wrf
       type (proj_lc_t) :: proj
-      integer :: ij, ifts, ifte, jfts, jfte, ims, ime, jms, jme
+      integer :: ims, ime, jms, jme
 
 
         ! Init
@@ -694,42 +694,14 @@ subroutine Interp_var2grid (this, lats_in, lons_in, ifms, ifme, jfms, jfme, &
 
       end select
 
+        ! Interpolate
       ims = 1
       ime = size (var_wrf, dim = 1)
       jms = 1
       jme = size (var_wrf, dim = 2)
-        ! Algorithm(s)
-      Hinterp: select case (hinterp_opt)
-        case (HINTERP_NEAREST)
-          !$OMP PARALLEL DO   &
-          !$OMP PRIVATE (ij, ifts, ifte, jfts, jfte)
-          do ij = 1, num_tiles
-            ifts = i_start(ij)
-            ifte = i_end(ij)
-            jfts = j_start(ij)
-            jfte = j_end(ij)
-            call Interp_horizontal_nearest (var_wrf, proj, ims, ime, jms, jme, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
-                lats_in, lons_in, vals_out)
-          end do
-          !$OMP END PARALLEL DO
-
-        case (HINTERP_BILINEAR)
-          !$OMP PARALLEL DO   &
-          !$OMP PRIVATE (ij, ifts, ifte, jfts, jfte)
-          do ij = 1, num_tiles
-            ifts = i_start(ij)
-            ifte = i_end(ij)
-            jfts = j_start(ij)
-            jfte = j_end(ij)
-            call Interp_horizontal_bilinear (var_wrf, proj, ims, ime, jms, jme, ifms, ifme, jfms, jfme, ifts, ifte, jfts, jfte, &
-                lats_in, lons_in, vals_out)
-          end do
-          !$OMP END PARALLEL DO
-
-        case default
-          call Stop_simulation ('The horizontal interpolation option selected does not exist.')
 
-      end select Hinterp
+      call Interp_horizontal (var_wrf, proj, ims, ime, jms, jme, ifms, ifme, jfms, jfme, &
+          num_tiles, i_start, i_end, j_start, j_end, hinterp_opt, lats_out, lons_out, data_out)
 
     end subroutine Interp_var2grid