ROMS
Loading...
Searching...
No Matches
grid_coords.F File Reference
#include "cppdefs.h"
Include dependency graph for grid_coords.F:

Go to the source code of this file.

Functions/Subroutines

subroutine grid_coords (ng, model)
 

Function/Subroutine Documentation

◆ grid_coords()

subroutine grid_coords ( integer, intent(in) ng,
integer, intent(in) model )

Definition at line 3 of file grid_coords.F.

4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This routine converts initial locations to fractional grid (I,J) !
13! coordinates, if appropriate. !
14! !
15!=======================================================================
16!
17 USE mod_param
18 USE mod_parallel
19# ifdef FLOATS
20 USE mod_floats
21# endif
22 USE mod_grid
23 USE mod_scalars
24 USE roms_interpolate_mod
25!
26# ifdef DISTRIBUTE
27 USE distribute_mod, ONLY : mp_collect
28# endif
29!
30 implicit none
31!
32! Imported variable declarations.
33!
34 integer, intent(in) :: ng, model
35!
36! Local variable declarations.
37!
38 integer :: IstrR, Iend, JstrR, Jend
39 integer :: LBi, UBi, LBj, UBj
40 integer :: i, j, k, l, mc
41
42 real(r8), parameter :: spv = 0.0_r8
43
44# ifdef FLOATS
45 real(r8) :: Xstr, Xend, Ystr, Yend, zfloat
46 logical, dimension(Nfloats(ng)) :: my_thread
47 real(r8), dimension(Nfloats(ng)) :: Iflt, Jflt
48# ifdef SOLVE3D
49 real(r8), dimension(Nfloats(ng)) :: Kflt
50# endif
51# endif
52
53# ifdef STATIONS
54 real(r8), dimension(Nstation(ng)) :: Slon, Slat
55 real(r8), dimension(Nstation(ng)) :: Ista, Jsta
56# endif
57!
58!-----------------------------------------------------------------------
59! Determine searching model grid box and arrays bounds.
60!-----------------------------------------------------------------------
61!
62# ifdef DISTRIBUTE
63 istrr=bounds(ng)%IstrR(myrank)
64 iend =bounds(ng)%Iend (myrank)
65 jstrr=bounds(ng)%JstrR(myrank)
66 jend =bounds(ng)%Jend (myrank)
67# else
68 istrr=0
69 iend =lm(ng)
70 jstrr=0
71 jend =mm(ng)
72# endif
73!
74 lbi=lbound(grid(ng)%h,dim=1)
75 ubi=ubound(grid(ng)%h,dim=1)
76 lbj=lbound(grid(ng)%h,dim=2)
77 ubj=ubound(grid(ng)%h,dim=2)
78
79# ifdef FLOATS
80!
81 xstr=real(bounds(ng)%Istr(myrank),r8)-0.5_r8
82 xend=real(bounds(ng)%Iend(myrank),r8)+0.5_r8
83 ystr=real(bounds(ng)%Jstr(myrank),r8)-0.5_r8
84 yend=real(bounds(ng)%Jend(myrank),r8)+0.5_r8
85!
86!-----------------------------------------------------------------------
87! If applicable, convert initial floats locations (Flon,Flat) to
88! fractional grid coordinates.
89!-----------------------------------------------------------------------
90!
91 IF (spherical) THEN
92 IF (lfloats(ng)) THEN
93 mc=drifter(ng)%Findex(0)
94 IF (drifter(ng)%Findex(0).gt.0) THEN
95 CALL hindices (ng, lbi, ubi, lbj, ubj, &
96 & istrr, iend+1, jstrr, jend+1, &
97 & grid(ng)%angler, &
98 & grid(ng)%lonr, &
99 & grid(ng)%latr, &
100 & 1, mc, 1, 1, &
101 & 1, mc, 1, 1, &
102 & drifter(ng)%Flon, &
103 & drifter(ng)%Flat, &
104 & iflt, jflt, spv, .false.)
105# ifdef DISTRIBUTE
106 CALL mp_collect (ng, model, mc, spv, iflt)
107 CALL mp_collect (ng, model, mc, spv, jflt)
108# endif
109 DO i=1,mc
110 l=drifter(ng)%Findex(i)
111 drifter(ng)%Tinfo(ixgrd,l)=iflt(i)
112 drifter(ng)%Tinfo(iygrd,l)=jflt(i)
113 END DO
114 END IF
115 END IF
116 END IF
117
118# ifdef SOLVE3D
119!
120! Determine which node bounds the initial float location.
121!
122# ifdef DISTRIBUTE
123 IF (lfloats(ng)) THEN
124 DO l=1,nfloats(ng)
125 IF ((xstr.le.drifter(ng)%Tinfo(ixgrd,l)).and. &
126 & (drifter(ng)%Tinfo(ixgrd,l).lt.xend).and. &
127 & (ystr.le.drifter(ng)%Tinfo(iygrd,l)).and. &
128 & (drifter(ng)%Tinfo(iygrd,l).lt.yend)) THEN
129 my_thread(l)=.true.
130 ELSE
131 my_thread(l)=.false.
132 END IF
133 END DO
134 END IF
135# else
136 DO l=1,nfloats(ng)
137 my_thread(l)=.true.
138 END DO
139# endif
140# endif
141!
142!-----------------------------------------------------------------------
143! Set float initial vertical level position, if inside application
144! grid. If the initial float depth (in meters) is not found, release
145! float at the surface model level.
146!-----------------------------------------------------------------------
147!
148 DO l=1,nfloats(ng)
149 IF (lfloats(ng)) THEN
150# ifdef SOLVE3D
151 drifter(ng)%Fz0(l)=spv
152 IF (my_thread(l).and. &
153 & ((drifter(ng)%Tinfo(ixgrd,l).ge.0.5_r8).and. &
154 & (drifter(ng)%Tinfo(iygrd,l).ge.0.5_r8).and. &
155 & (drifter(ng)%Tinfo(ixgrd,l).le. &
156 & real(lm(ng),r8)+0.5_r8).and. &
157 & (drifter(ng)%Tinfo(iygrd,l).le. &
158 & real(mm(ng),r8)+0.5_r8))) THEN
159 zfloat=drifter(ng)%Tinfo(izgrd,l)
160 drifter(ng)%Fz0(l)=zfloat ! Save original value
161 kflt(l)=zfloat
162 IF (zfloat.le.0.0_r8) THEN
163 i=int(drifter(ng)%Tinfo(ixgrd,l)) ! Fractional positions
164 j=int(drifter(ng)%Tinfo(iygrd,l)) ! are still in this cell
165 IF (zfloat.lt.grid(ng)%z_w(i,j,0)) THEN
166 zfloat=grid(ng)%z_w(i,j,0)+5.0_r8
167 drifter(ng)%Fz0(l)=zfloat
168 END IF
169 drifter(ng)%Tinfo(izgrd,l)=real(n(ng),r8)
170 DO k=n(ng),1,-1
171 IF ((grid(ng)%z_w(i,j,k)-zfloat)* &
172 & (zfloat-grid(ng)%z_w(i,j,k-1)).ge.0.0_r8) THEN
173 kflt(l)=real(k-1,r8)+ &
174 & (zfloat-grid(ng)%z_w(i,j,k-1))/ &
175 & grid(ng)%Hz(i,j,k)
176 END IF
177 END DO
178 END IF
179 ELSE
180 kflt(l)=spv
181 END IF
182# else
183 drifter(ng)%Tinfo(izgrd,l)=0.0_r8
184# endif
185 END IF
186 END DO
187# ifdef SOLVE3D
188 IF (lfloats(ng)) THEN
189# ifdef DISTRIBUTE
190 CALL mp_collect (ng, model, nfloats(ng), spv, drifter(ng)%Fz0)
191 CALL mp_collect (ng, model, nfloats(ng), spv, kflt)
192# endif
193 DO l=1,nfloats(ng)
194 drifter(ng)%Tinfo(izgrd,l)=kflt(l)
195 END DO
196 END IF
197# endif
198# endif
199# ifdef STATIONS
200!
201!-----------------------------------------------------------------------
202! If applicable, convert station locations (SposX,SposY) to fractional
203! grid coordinates.
204!-----------------------------------------------------------------------
205!
206 IF (spherical) THEN
207 mc=0
208 DO l=1,nstation(ng)
209 IF (scalars(ng)%Sflag(l).gt.0) THEN
210 mc=mc+1
211 slon(mc)=scalars(ng)%SposX(l)
212 slat(mc)=scalars(ng)%SposY(l)
213 END IF
214 END DO
215 IF (mc.gt.0) THEN
216 CALL hindices (ng, lbi, ubi, lbj, ubj, &
217 & istrr, iend+1, jstrr, jend+1, &
218 & grid(ng)%angler, &
219 & grid(ng)%lonr, &
220 & grid(ng)%latr, &
221 & 1, mc, 1, 1, &
222 & 1, mc, 1, 1, &
223 & slon, slat, &
224 & ista, jsta, &
225 & spv, .false.)
226# ifdef DISTRIBUTE
227 CALL mp_collect (ng, model, mc, spv, ista)
228 CALL mp_collect (ng, model, mc, spv, jsta)
229# endif
230 mc=0
231 DO l=1,nstation(ng)
232 IF (scalars(ng)%Sflag(l).gt.0) THEN
233 mc=mc+1
234 scalars(ng)%SposX(l)=ista(mc)
235 scalars(ng)%SposY(l)=jsta(mc)
236 END IF
237 END DO
238 END IF
239 END IF
240# endif
241 RETURN
subroutine hindices(ng, lbi, ubi, lbj, ubj, is, ie, js, je, angler, xgrd, ygrd, lbm, ubm, lbn, ubn, ms, me, ns, ne, xpos, ypos, ipos, jpos, ijspv, rectangular)
integer, parameter iygrd
Definition mod_floats.F:82
integer, parameter ixgrd
Definition mod_floats.F:81
type(t_drifter), dimension(:), allocatable drifter
Definition mod_floats.F:67
integer, parameter izgrd
Definition mod_floats.F:83
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer, dimension(:), allocatable nfloats
Definition mod_param.F:543
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
integer, dimension(:), allocatable nstation
Definition mod_param.F:557
integer, dimension(:), allocatable lm
Definition mod_param.F:455
integer, dimension(:), allocatable mm
Definition mod_param.F:456
logical spherical
logical, dimension(:), allocatable lfloats
type(t_scalars), dimension(:), allocatable scalars
Definition mod_scalars.F:65

References mod_param::bounds, mod_floats::drifter, mod_grid::grid, hindices(), mod_floats::ixgrd, mod_floats::iygrd, mod_floats::izgrd, mod_scalars::lfloats, mod_param::lm, mod_param::mm, mod_parallel::myrank, mod_param::n, mod_param::nfloats, mod_param::nstation, mod_scalars::scalars, and mod_scalars::spherical.

Referenced by ad_initial(), roms_kernel_mod::adm_initial(), initial(), roms_kernel_mod::nlm_initial(), rp_initial(), tl_initial(), and roms_kernel_mod::tlm_initial().

Here is the call graph for this function:
Here is the caller graph for this function: