/[MITgcm]/MITgcm/pkg/profiles/profiles_interp.F
ViewVC logotype

Annotation of /MITgcm/pkg/profiles/profiles_interp.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.6 - (hide annotations) (download)
Fri Jul 14 22:12:24 2006 UTC (17 years, 10 months ago) by gforget
Branch: MAIN
CVS Tags: checkpoint58u_post, checkpoint58w_post, checkpoint58r_post, checkpoint58n_post, checkpoint58x_post, checkpoint58t_post, checkpoint58q_post, checkpoint59a, checkpoint59c, checkpoint59b, checkpoint59, checkpoint58o_post, checkpoint58y_post, checkpoint58v_post, checkpoint58s_post, checkpoint58p_post, checkpoint58m_post
Changes since 1.5: +31 -0 lines
adding ptracer and ssh

1 gforget 1.5 C $Header: /u/gcmpack/MITgcm/pkg/profiles/profiles_interp.F,v 1.4 2006/05/06 15:14:01 heimbach Exp $
2 heimbach 1.3 C $Name: $
3    
4     #include "PROFILES_OPTIONS.h"
5 heimbach 1.1
6     C o==========================================================o
7     C | subroutine profiles_interp |
8     C | o 3D interpolation of model counterparts |
9     C | for netcdf profiles data |
10     C | started: Gael Forget 15-March-2006 |
11     C o==========================================================o
12    
13     SUBROUTINE profiles_interp(
14     O traj_cur_out,
15     I lon_cur,
16     I lat_cur,
17     I type_cur,
18     I file_cur,
19     I mytime,
20 gforget 1.5 I bi,
21     I bj,
22 heimbach 1.1 I myThid
23     & )
24    
25     implicit none
26    
27     C ==================== Global Variables ===========================
28     #include "EEPARAMS.h"
29     #include "SIZE.h"
30     #include "GRID.h"
31     #include "DYNVARS.h"
32     #include "PARAMS.h"
33 heimbach 1.3 #ifdef ALLOW_CAL
34 heimbach 1.1 #include "cal.h"
35 heimbach 1.3 #endif
36 heimbach 1.4 #ifdef ALLOW_PROFILES
37 heimbach 1.1 # include "profiles.h"
38     #else
39     integer NLEVELMAX
40     parameter (NLEVELMAX=1)
41     #endif
42 gforget 1.6 #ifdef ALLOW_PTRACERS
43     #include "PTRACERS_SIZE.h"
44     #include "PTRACERS.h"
45     #endif
46 heimbach 1.1 C ==================== Routine Variables ==========================
47     _RL mytime
48     integer mythid
49     integer type_cur,file_cur
50     _RL traj_cur_out(NLEVELMAX)
51     _RL lon_cur,lat_cur
52    
53 heimbach 1.4 #ifdef ALLOW_PROFILES
54 heimbach 1.1
55     C ==================== Local Variables ==========================
56     _RL tab_coeffs1(2,2),tab_coeffs2(2,2),tab_coeffs3(2,2)
57     _RL ponderations(2,2),pondsSUM,distance1,distance2
58     integer i,j,k,kk,kcur,iG,jG,bi,bj
59     _RL traj_cur(nR),mask_cur(nR)
60     integer prof_i,prof_j
61 gforget 1.2 _RL lon_tmp1,lon_tmp2,lon_1,lon_2,lat_1,lat_2,tmp_coeff
62 heimbach 1.1 c-- == end of interface ==
63    
64 gforget 1.2 prof_i=-10
65     prof_j=-10
66     lon_1=-10
67     lon_2=-10
68     lat_1=-10
69     lat_2=-10
70 heimbach 1.1
71     DO j=1,sNy+1
72     DO i=1,sNx+1
73 gforget 1.2
74 heimbach 1.1 cgf value of j, south of the data point:
75     if (type_cur.NE.4) then
76 gforget 1.2 if ((yC(i,j,bi,bj).LE.lat_cur).AND.
77     & (yC(i,j+1,bi,bj).GT.lat_cur)) then
78     prof_j=j
79     lat_1=yC(i,j,bi,bj)
80     lat_2=yC(i,j+1,bi,bj)
81     else
82     prof_j=prof_j
83     lat_1=lat_1
84     lat_2=lat_2
85     endif
86 heimbach 1.1 else
87 gforget 1.2 if ((yG(i,j,bi,bj).LE.lat_cur).AND.
88     & (yG(i,j+1,bi,bj).GT.lat_cur)) then
89     prof_j=j
90     lat_1=yG(i,j,bi,bj)
91     lat_2=yG(i,j+1,bi,bj)
92     else
93     prof_j=prof_j
94     lat_1=lat_1
95     lat_2=lat_2
96     endif
97 heimbach 1.1 endif
98 gforget 1.2
99 heimbach 1.1 cgf value of i, west of the data point:
100     if (type_cur.NE.3) then
101 gforget 1.2 if (xC(i+1,j,bi,bj).LT.xC(1,j,bi,bj)) then
102     lon_tmp2=xC(i+1,j,bi,bj)+360
103 heimbach 1.1 else
104 gforget 1.2 lon_tmp2=xC(i+1,j,bi,bj)
105 heimbach 1.1 endif
106 gforget 1.2 if (xC(i,j,bi,bj).LT.xC(1,j,bi,bj)) then
107     lon_tmp1=xC(i,j,bi,bj)+360
108 heimbach 1.1 else
109 gforget 1.2 lon_tmp1=xC(i,j,bi,bj)
110 heimbach 1.1 endif
111     else
112 gforget 1.2 if (xG(i+1,j,bi,bj).LT.xG(1,j,bi,bj)) then
113     lon_tmp2=xG(i+1,j,bi,bj)+360
114 heimbach 1.1 else
115 gforget 1.2 lon_tmp2=xG(i+1,j,bi,bj)
116 heimbach 1.1 endif
117 gforget 1.2 if (xG(i,j,bi,bj).LT.xG(1,j,bi,bj)) then
118     lon_tmp1=xG(i,j,bi,bj)+360
119 heimbach 1.1 else
120 gforget 1.2 lon_tmp1=xG(i,j,bi,bj)
121 heimbach 1.1 endif
122 gforget 1.2 endif
123     if ((lon_tmp1.LE.lon_cur).AND.
124     &(lon_tmp2.GT.lon_cur)) then
125     prof_i=i
126     lon_1=lon_tmp1
127     lon_2=lon_tmp2
128     else
129     prof_i=prof_i
130     lon_1=lon_1
131     lon_2=lon_2
132 heimbach 1.1 endif
133    
134 gforget 1.2 ENDDO
135     ENDDO
136    
137    
138     if ((prof_i.NE.-10).AND.(prof_j.NE.-10)) then
139     cgf) spatial interpolation
140    
141     distance1=(lat_cur-lat_1)/(lat_2-lat_1)
142     distance2=(lon_cur-lon_1)/(lon_2-lon_1)
143     tab_coeffs2(1,1)=(1-distance2)*(1-distance1)
144     tab_coeffs2(1,2)=distance2*(1-distance1)
145     tab_coeffs2(2,1)=(1-distance2)*distance1
146     tab_coeffs2(2,2)=distance2*distance1
147 heimbach 1.1
148     do k=1,nr
149     if (type_cur.EQ.1) then
150     tab_coeffs1(1,1)=theta(prof_i,prof_j,k,bi,bj) !SO
151     tab_coeffs1(1,2)=theta(prof_i+1,prof_j,k,bi,bj) !SE
152     tab_coeffs1(2,1)=theta(prof_i,prof_j+1,k,bi,bj) !NO
153     tab_coeffs1(2,2)=theta(prof_i+1,prof_j+1,k,bi,bj) !NZ
154     tab_coeffs3(1,1)=maskC(prof_i,prof_j,k,bi,bj) !SO
155     tab_coeffs3(1,2)=maskC(prof_i+1,prof_j,k,bi,bj) !SE
156     tab_coeffs3(2,1)=maskC(prof_i,prof_j+1,k,bi,bj) !NO
157     tab_coeffs3(2,2)=maskC(prof_i+1,prof_j+1,k,bi,bj) !NZ
158     elseif (type_cur.EQ.2) then
159     tab_coeffs1(1,1)=salt(prof_i,prof_j,k,bi,bj) !SO
160     tab_coeffs1(1,2)=salt(prof_i+1,prof_j,k,bi,bj) !SE
161     tab_coeffs1(2,1)=salt(prof_i,prof_j+1,k,bi,bj) !NO
162     tab_coeffs1(2,2)=salt(prof_i+1,prof_j+1,k,bi,bj) !NZ
163     tab_coeffs3(1,1)=maskC(prof_i,prof_j,k,bi,bj) !SO
164     tab_coeffs3(1,2)=maskC(prof_i+1,prof_j,k,bi,bj) !SE
165     tab_coeffs3(2,1)=maskC(prof_i,prof_j+1,k,bi,bj) !NO
166     tab_coeffs3(2,2)=maskC(prof_i+1,prof_j+1,k,bi,bj) !NZ
167     elseif (type_cur.EQ.3) then
168     tab_coeffs1(1,1)=uVel(prof_i,prof_j,k,bi,bj) !SO
169     tab_coeffs1(1,2)=uVel(prof_i+1,prof_j,k,bi,bj) !SE
170     tab_coeffs1(2,1)=uVel(prof_i,prof_j+1,k,bi,bj) !NO
171     tab_coeffs1(2,2)=uVel(prof_i+1,prof_j+1,k,bi,bj) !NZ
172     tab_coeffs3(1,1)=maskW(prof_i,prof_j,k,bi,bj) !SO
173     tab_coeffs3(1,2)=maskW(prof_i+1,prof_j,k,bi,bj) !SE
174     tab_coeffs3(2,1)=maskW(prof_i,prof_j+1,k,bi,bj) !NO
175     tab_coeffs3(2,2)=maskW(prof_i+1,prof_j+1,k,bi,bj) !NZ
176     elseif (type_cur.EQ.4) then
177     tab_coeffs1(1,1)=vVel(prof_i,prof_j,k,bi,bj) !SO
178     tab_coeffs1(1,2)=vVel(prof_i+1,prof_j,k,bi,bj) !SE
179     tab_coeffs1(2,1)=vVel(prof_i,prof_j+1,k,bi,bj) !NO
180     tab_coeffs1(2,2)=vVel(prof_i+1,prof_j+1,k,bi,bj) !NZ
181     tab_coeffs3(1,1)=maskS(prof_i,prof_j,k,bi,bj) !SO
182     tab_coeffs3(1,2)=maskS(prof_i+1,prof_j,k,bi,bj) !SE
183     tab_coeffs3(2,1)=maskS(prof_i,prof_j+1,k,bi,bj) !NO
184     tab_coeffs3(2,2)=maskS(prof_i+1,prof_j+1,k,bi,bj) !NZ
185 gforget 1.6 elseif (type_cur.EQ.5) then
186     #ifdef ALLOW_PTRACERS
187     cgf if this gets used, an additional common block could be defined, containing
188     cgf the pTracer number (now 1, hard-coded), that would be read from the .nc input file
189     tab_coeffs1(1,1)=pTracer(prof_i,prof_j,k,bi,bj,1) !SO
190     tab_coeffs1(1,2)=pTracer(prof_i+1,prof_j,k,bi,bj,1) !SE
191     tab_coeffs1(2,1)=pTracer(prof_i,prof_j+1,k,bi,bj,1) !NO
192     tab_coeffs1(2,2)=pTracer(prof_i+1,prof_j+1,k,bi,bj,1) !NZ
193     #else
194     tab_coeffs1(1,1)=0 !SO
195     tab_coeffs1(1,2)=0 !SE
196     tab_coeffs1(2,1)=0 !NO
197     tab_coeffs1(2,2)=0 !NZ
198     #endif
199     tab_coeffs3(1,1)=maskC(prof_i,prof_j,k,bi,bj) !SO
200     tab_coeffs3(1,2)=maskC(prof_i+1,prof_j,k,bi,bj) !SE
201     tab_coeffs3(2,1)=maskC(prof_i,prof_j+1,k,bi,bj) !NO
202     tab_coeffs3(2,2)=maskC(prof_i+1,prof_j+1,k,bi,bj) !NZ
203     elseif (type_cur.EQ.6) then
204     tab_coeffs1(1,1)=etan(prof_i,prof_j,bi,bj) !SO
205     tab_coeffs1(1,2)=etan(prof_i+1,prof_j,bi,bj) !SE
206     tab_coeffs1(2,1)=etan(prof_i,prof_j+1,bi,bj) !NO
207     tab_coeffs1(2,2)=etan(prof_i+1,prof_j+1,bi,bj) !NZ
208     tab_coeffs3(1,1)=maskC(prof_i,prof_j,1,bi,bj) !SO
209     tab_coeffs3(1,2)=maskC(prof_i+1,prof_j,1,bi,bj) !SE
210     tab_coeffs3(2,1)=maskC(prof_i,prof_j+1,1,bi,bj) !NO
211     tab_coeffs3(2,2)=maskC(prof_i+1,prof_j+1,1,bi,bj) !NZ
212 heimbach 1.1 else
213 gforget 1.2 tab_coeffs1(1,1)=0.
214     tab_coeffs1(2,1)=0.
215     tab_coeffs1(1,2)=0.
216     tab_coeffs1(2,2)=0.
217     tab_coeffs3(1,1)=0.
218     tab_coeffs3(2,1)=0.
219     tab_coeffs3(1,2)=0.
220     tab_coeffs3(2,2)=0.
221 heimbach 1.1 endif
222    
223 gforget 1.2 ponderations(1,1)=tab_coeffs3(1,1)*tab_coeffs2(1,1)
224     ponderations(1,2)=tab_coeffs3(1,2)*tab_coeffs2(1,2)
225     ponderations(2,1)=tab_coeffs3(2,1)*tab_coeffs2(2,1)
226     ponderations(2,2)=tab_coeffs3(2,2)*tab_coeffs2(2,2)
227     pondsSUM=ponderations(1,1)+ponderations(2,1)+ponderations(1,2)+
228     & ponderations(2,2)
229    
230     if (pondsSUM.GT.0) then
231 heimbach 1.1 tab_coeffs1(1,1)=tab_coeffs1(1,1)*ponderations(1,1)/pondsSUM
232     tab_coeffs1(1,2)=tab_coeffs1(1,2)*ponderations(1,2)/pondsSUM
233     tab_coeffs1(2,1)=tab_coeffs1(2,1)*ponderations(2,1)/pondsSUM
234     tab_coeffs1(2,2)=tab_coeffs1(2,2)*ponderations(2,2)/pondsSUM
235 gforget 1.2 traj_cur(k)=tab_coeffs1(1,1)+tab_coeffs1(2,1)+
236     & tab_coeffs1(1,2)+tab_coeffs1(2,2)
237     mask_cur(k)=1
238     else
239 heimbach 1.1 traj_cur(k)=0
240     mask_cur(k)=0
241 gforget 1.2 endif
242     enddo
243 heimbach 1.1
244     else
245 gforget 1.2 do k=1,nr
246 heimbach 1.1 traj_cur(k)=0
247     mask_cur(k)=0
248 gforget 1.2 enddo
249 heimbach 1.1 endif
250    
251     cgf vertical interpolation:
252     do kk=1,NLEVELMAX
253     traj_cur_out(kk)=0
254 gforget 1.5 prof_mask1D_cur(kk,bi,bj)=0
255 heimbach 1.1 enddo
256 gforget 1.5 do kk=1,ProfDepthNo(file_cur,bi,bj)
257 heimbach 1.1 c case 1: above first grid center=> first grid center value
258 gforget 1.5 if (prof_depth(file_cur,kk,bi,bj).LT.-rC(1)) then
259 gforget 1.2 traj_cur_out(kk)=traj_cur(1)
260 gforget 1.5 prof_mask1D_cur(kk,bi,bj)=mask_cur(1)
261 gforget 1.2 c case 2: just below last grid center=> last cell value
262 gforget 1.5 elseif (prof_depth(file_cur,kk,bi,bj).GE.-rC(nr)) then
263     if ( prof_depth(file_cur,kk,bi,bj) .LT.
264 gforget 1.2 & (-rC(nr)+drC(nr)/2) ) then
265     traj_cur_out(kk)=traj_cur(nr)
266 gforget 1.5 prof_mask1D_cur(kk,bi,bj)=mask_cur(nr)
267 gforget 1.2 endif
268 heimbach 1.1 c case 3: between two grid centers
269 gforget 1.2 else
270     kcur=0
271     do k=1,nr-1
272 gforget 1.5 if ((prof_depth(file_cur,kk,bi,bj).GE.-rC(k)).AND.
273     & (prof_depth(file_cur,kk,bi,bj).LT.-rC(k+1))) then
274 gforget 1.2 kcur=k
275     endif
276     enddo
277     if (kcur.EQ.0) then
278 heimbach 1.1 print*,"profiles_interp unexpected case: stop 1"
279     stop
280 gforget 1.2 endif
281     if (mask_cur(kcur+1).EQ.1.) then
282     c subcase 1: 2 wet points=>linear interpolation
283 gforget 1.5 tmp_coeff=(prof_depth(file_cur,kk,bi,bj)+rC(kcur))/
284 gforget 1.2 & (-rC(kcur+1)+rC(kcur))
285     traj_cur_out(kk)=(1-tmp_coeff)*traj_cur(kcur)
286     & +tmp_coeff*traj_cur(kcur+1)
287 gforget 1.5 prof_mask1D_cur(kk,bi,bj)=1
288 gforget 1.2 if (mask_cur(kcur).EQ.0.) then
289 heimbach 1.1 print*,"profiles_interp unexpected case: stop 2"
290     stop
291 gforget 1.2 endif
292 gforget 1.5 elseif (prof_depth(file_cur,kk,bi,bj).LT.-rF(kcur+1)) then
293 gforget 1.2 c subcase 2: only 1 wet point just above=>upper cell value
294     traj_cur_out(kk)=traj_cur(kcur)
295 gforget 1.5 prof_mask1D_cur(kk,bi,bj)=mask_cur(kcur)
296 gforget 1.2 endif
297     endif
298 heimbach 1.1 enddo
299    
300    
301     #endif
302    
303     end
304    

  ViewVC Help
Powered by ViewVC 1.1.22