/[MITgcm]/MITgcm_contrib/SOSE/code_ad/ctrl_map_ini_gen.F
ViewVC logotype

Annotation of /MITgcm_contrib/SOSE/code_ad/ctrl_map_ini_gen.F

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


Revision 1.1 - (hide annotations) (download)
Fri Apr 23 19:55:11 2010 UTC (15 years, 3 months ago) by mmazloff
Branch: MAIN
CVS Tags: HEAD
original files

1 mmazloff 1.1 C $Header: /u/gcmpack/MITgcm/pkg/ctrl/ctrl_map_ini_gen.F,v 1.5 2009/04/28 18:09:28 jmc Exp $
2     C $Name: $
3    
4     #include "CTRL_CPPOPTIONS.h"
5    
6    
7     subroutine ctrl_map_ini_gen3D(xxFileCur, wFileCur, xxDummyCur,
8     & boundsVec, paramFld3d, maskFld3d, paramSmooth, mythid )
9    
10     c ==================================================================
11     c SUBROUTINE ctrl_map_ini_gen3D
12     c ==================================================================
13     c
14     c started: Gael Forget gforget@mit.edu 8-Feb-2008
15     c
16     c - Generetic routine for an individual 3D control term
17     c (to be called from ctrl_map_ini in a loop e.g.)
18     c
19     c ==================================================================
20     c SUBROUTINE ctrl_map_ini_gen3D
21     c ==================================================================
22    
23     implicit none
24    
25     c == global variables ==
26    
27     #include "EEPARAMS.h"
28     #include "SIZE.h"
29     #include "PARAMS.h"
30     #include "GRID.h"
31     #include "DYNVARS.h"
32     #include "FFIELDS.h"
33    
34     #include "ctrl.h"
35     #include "ctrl_dummy.h"
36     #include "optim.h"
37     #ifdef ALLOW_ECCO
38     #include "ecco_cost.h"
39     #endif
40    
41     #ifdef ALLOW_AUTODIFF_TAMC
42     #include "tamc.h"
43     #include "tamc_keys.h"
44     #endif /* ALLOW_AUTODIFF_TAMC */
45    
46     c == routine arguments ==
47    
48     integer mythid
49     character*(*) wFileCur,xxFileCur
50     _RL boundsVec(5),tmpMax,xxDummyCur
51    
52     _RL wFld3d(1-olx:snx+olx,1-oly:sny+oly,nr,nsx,nsy)
53     _RL xxFld3d(1-olx:snx+olx,1-oly:sny+oly,nr,nsx,nsy)
54     _RL paramFld3d(1-olx:snx+olx,1-oly:sny+oly,nr,nsx,nsy)
55     _RL maskFld3d(1-olx:snx+olx,1-oly:sny+oly,nr,nsx,nsy)
56     integer paramSmooth
57    
58     c == local variables ==
59    
60     integer bi,bj
61     integer i,j,k
62     integer itlo,ithi
63     integer jtlo,jthi
64     integer jmin,jmax
65     integer imin,imax
66     integer il
67    
68     logical doglobalread
69     logical ladinit
70    
71     character*( 80) fnamegeneric
72    
73     c == external ==
74    
75     integer ilnblnk
76     external ilnblnk
77    
78     c == end of interface ==
79    
80     #ifdef ALLOW_AUTODIFF_TAMC
81     act3 = myThid - 1
82     max3 = nTx*nTy
83     act4 = 0
84     ikey = (act3 + 1) + act4*max3
85     #endif /* ALLOW_AUTODIFF_TAMC */
86    
87     jtlo = mybylo(mythid)
88     jthi = mybyhi(mythid)
89     itlo = mybxlo(mythid)
90     ithi = mybxhi(mythid)
91     c-- only do interior, and exchange outside
92     jmin = 1
93     jmax = sny
94     imin = 1
95     imax = snx
96    
97     doglobalread = .false.
98     ladinit = .false.
99    
100    
101     call mdsreadfield(wFileCur,32,'RL',nR,wFld3d,1,mythid)
102     _EXCH_XYZ_RL( wFld3d, mythid )
103    
104     il=ilnblnk( xxFileCur )
105     write(fnamegeneric(1:80),'(2a,i10.10)')
106     & xxFileCur(1:il),'.',optimcycle
107     call active_read_xyz( fnamegeneric, xxFld3d, 1,
108     & doglobalread, ladinit, optimcycle, mythid, xxDummyCur )
109    
110    
111     #ifndef ALLOW_SMOOTH_CORREL3D
112     c avoid xx larger than boundsVec(5) X uncertainty
113     #ifdef CMM_dont_want
114     if ( boundsVec(5).GT.0.) then
115     do bj = jtlo,jthi
116     do bi = itlo,ithi
117     do k = 1,nr
118     do j = jmin,jmax
119     do i = imin,imax
120     if ( (maskFld3d(i,j,k,bi,bj).NE.0.).AND.
121     & (wFld3d(i,j,k,bi,bj).GT.0.) ) then
122     tmpMax=boundsVec(5)/sqrt(wFld3d(i,j,k,bi,bj))
123     if ( abs(xxFld3d(i,j,k,bi,bj)).GT.tmpMax ) then
124     xxFld3d(i,j,k,bi,bj)=sign(tmpMax,xxFld3d(i,j,k,bi,bj))
125     else
126     xxFld3d(i,j,k,bi,bj)=xxFld3d(i,j,k,bi,bj)
127     endif
128     endif
129     enddo
130     enddo
131     enddo
132     enddo
133     enddo
134     endif
135     #endif
136     #else
137     c apply Weaver And Courtier correlation operator
138     CMM( my ctrls are dimensional so need to remove this
139     c scale param adjustment
140     do bj = jtlo,jthi
141     do bi = itlo,ithi
142     do k = 1,nr
143     do j = jmin,jmax
144     do i = imin,imax
145     if ( (maskFld3d(i,j,k,bi,bj).NE.0.)
146     & .AND. (wFld3d(i,j,k,bi,bj).GT.0.) ) then
147     xxFld3d(i,j,k,bi,bj)=xxFld3d(i,j,k,bi,bj)
148     & *sqrt( wFld3d(i,j,k,bi,bj) )
149     else
150     xxFld3d(i,j,k,bi,bj)=0.
151     endif
152     enddo
153     enddo
154     enddo
155     enddo
156     enddo
157     CMM)
158     if (paramSmooth.NE.0) then
159     call smooth_correl3D(xxFld3d,paramSmooth,mythid)
160     endif
161     #endif
162    
163     do bj = jtlo,jthi
164     do bi = itlo,ithi
165     do k = 1,nr
166     do j = jmin,jmax
167     do i = imin,imax
168     #ifdef ALLOW_SMOOTH_CORREL3D
169     c scale param adjustment
170     if ( (maskFld3d(i,j,k,bi,bj).NE.0.)
171     & .AND. (wFld3d(i,j,k,bi,bj).GT.0.) ) then
172     xxFld3d(i,j,k,bi,bj)=xxFld3d(i,j,k,bi,bj)
173     & /sqrt( wFld3d(i,j,k,bi,bj) )
174     else
175     xxFld3d(i,j,k,bi,bj)=0.
176     endif
177     #endif
178     paramFld3d(i,j,k,bi,bj) = paramFld3d(i,j,k,bi,bj)
179     & + xxFld3d(i,j,k,bi,bj)
180     enddo
181     enddo
182     enddo
183     enddo
184     enddo
185    
186     c avoid param out of [boundsVec(1) boundsVec(4)]
187     CMM no need CALL CTRL_BOUND_3D(paramFld3d,maskFld3d,boundsVec,myThid)
188    
189     #ifdef ALLOW_SMOOTH_CORREL3D
190     write(fnamegeneric(1:80),'(2a,i10.10)')
191     & xxFileCur(1:il),'.effective.',optimcycle
192     call mdswritefield(fnamegeneric,ctrlprec,.FALSE.,'RL',
193     & nr, paramFld3d, 1, optimcycle, mythid)
194     #endif
195    
196    
197     end
198    
199    
200    
201     subroutine ctrl_map_ini_gen2D(xxFileCur, wFileCur, xxDummyCur,
202     & boundsVec, paramFld2d, maskFld3d, paramSmooth, mythid )
203    
204     c ==================================================================
205     c SUBROUTINE ctrl_map_ini_gen2D
206     c ==================================================================
207     c
208     c started: Gael Forget gforget@mit.edu 8-Feb-2008
209     c
210     c - Generetic routine for an individual 2D control term
211     c (to be called from ctrl_map_ini in a loop e.g.)
212     c
213     c ==================================================================
214     c SUBROUTINE ctrl_map_ini_gen3D
215     c ==================================================================
216    
217     implicit none
218    
219     c == global variables ==
220    
221     #include "EEPARAMS.h"
222     #include "SIZE.h"
223     #include "PARAMS.h"
224     #include "GRID.h"
225     #include "DYNVARS.h"
226     #include "FFIELDS.h"
227    
228     #include "ctrl.h"
229     #include "ctrl_dummy.h"
230     #include "optim.h"
231     #ifdef ALLOW_ECCO
232     #include "ecco_cost.h"
233     #endif
234    
235     #ifdef ALLOW_AUTODIFF_TAMC
236     #include "tamc.h"
237     #include "tamc_keys.h"
238     #endif /* ALLOW_AUTODIFF_TAMC */
239    
240     c == routine arguments ==
241    
242     integer mythid
243     character*(*) wFileCur,xxFileCur
244     _RL boundsVec(5),tmpMax,xxDummyCur
245    
246     _RL wFld2d(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy)
247     _RL xxFld2d(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy)
248     _RL paramFld2d(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy)
249     _RL maskFld3d(1-olx:snx+olx,1-oly:sny+oly,nr,nsx,nsy)
250     integer paramSmooth
251    
252     c == local variables ==
253    
254     integer bi,bj
255     integer i,j,k
256     integer itlo,ithi
257     integer jtlo,jthi
258     integer jmin,jmax
259     integer imin,imax
260     integer il
261    
262     logical doglobalread
263     logical ladinit
264    
265     character*( 80) fnamegeneric
266    
267     c == external ==
268    
269     integer ilnblnk
270     external ilnblnk
271    
272     c == end of interface ==
273    
274     #ifdef ALLOW_AUTODIFF_TAMC
275     act3 = myThid - 1
276     max3 = nTx*nTy
277     act4 = 0
278     ikey = (act3 + 1) + act4*max3
279     #endif /* ALLOW_AUTODIFF_TAMC */
280    
281     jtlo = mybylo(mythid)
282     jthi = mybyhi(mythid)
283     itlo = mybxlo(mythid)
284     ithi = mybxhi(mythid)
285     c-- only do interior, and exchange outside
286     jmin = 1
287     jmax = sny
288     imin = 1
289     imax = snx
290    
291     doglobalread = .false.
292     ladinit = .false.
293    
294    
295     call mdsreadfield(wFileCur,32,'RL',1,wFld2d,1,mythid)
296     _EXCH_XY_RL( wFld2d, mythid )
297    
298     il=ilnblnk( xxFileCur )
299     write(fnamegeneric(1:80),'(2a,i10.10)')
300     & xxFileCur(1:il),'.',optimcycle
301     call active_read_xy( fnamegeneric, xxFld2d, 1,
302     & doglobalread, ladinit, optimcycle, mythid, xxDummyCur )
303    
304     CMM(#ifndef ALLOW_SMOOTH_CORREL3D
305     #ifndef ALLOW_SMOOTH_CORREL2D
306     CMM)
307     c avoid xx larger than boundsVec(5) X uncertainty
308     if ( boundsVec(5).GT.0.) then
309     do bj = jtlo,jthi
310     do bi = itlo,ithi
311     do j = jmin,jmax
312     do i = imin,imax
313     if ( (maskFld3d(i,j,1,bi,bj).NE.0.).AND.
314     & (wFld2d(i,j,bi,bj).GT.0.) ) then
315     tmpMax=boundsVec(5)/sqrt(wFld2d(i,j,bi,bj))
316     if ( abs(xxFld2d(i,j,bi,bj)).GT.tmpMax ) then
317     xxFld2d(i,j,bi,bj)=sign(tmpMax,xxFld2d(i,j,bi,bj))
318     else
319     xxFld2d(i,j,bi,bj)=xxFld2d(i,j,bi,bj)
320     endif
321     endif
322     enddo
323     enddo
324     enddo
325     enddo
326     endif
327     #else
328     c apply Weaver And Courtier correlation operator
329     CMM( my controls are dimensional so get rid of that
330     do bj = jtlo,jthi
331     do bi = itlo,ithi
332     do j = jmin,jmax
333     do i = imin,imax
334     if ( (maskFld3d(i,j,1,bi,bj).NE.0.)
335     & .AND. (wFld2d(i,j,bi,bj).GT.0.) ) then
336     xxFld2d(i,j,bi,bj)=xxFld2d(i,j,bi,bj)
337     & *sqrt( wFld2d(i,j,bi,bj) )
338     else
339     xxFld2d(i,j,bi,bj)=0.
340     endif
341     enddo
342     enddo
343     enddo
344     enddo
345     CMM)
346     if (paramSmooth.NE.0) then
347     call smooth_correl2d(xxFld2d,maskFld3d,paramSmooth,mythid)
348     endif
349     #endif
350    
351     do bj = jtlo,jthi
352     do bi = itlo,ithi
353     do j = jmin,jmax
354     do i = imin,imax
355     #ifdef ALLOW_SMOOTH_CORREL2D
356     c scale param adjustment
357     if ( (maskFld3d(i,j,1,bi,bj).NE.0.)
358     & .AND. (wFld2d(i,j,bi,bj).GT.0.) ) then
359     xxFld2d(i,j,bi,bj)=xxFld2d(i,j,bi,bj)
360     & /sqrt( wFld2d(i,j,bi,bj) )
361     else
362     xxFld2d(i,j,bi,bj)=0.
363     endif
364     #endif
365     paramFld2d(i,j,bi,bj) = paramFld2d(i,j,bi,bj)
366     & + xxFld2d(i,j,bi,bj)
367     enddo
368     enddo
369     enddo
370     enddo
371    
372     CALL CTRL_BOUND_2D(paramFld2d,maskFld3d,boundsVec,myThid)
373    
374     #ifdef ALLOW_SMOOTH_CORREL2D
375     write(fnamegeneric(1:80),'(2a,i10.10)')
376     & xxFileCur(1:il),'.effective.',optimcycle
377     CMM( make consistent with ctrl_get_gen.F
378     call mdswritefield(fnamegeneric,ctrlprec,.FALSE.,'RL',
379     & 1, xxFld2d, 1, optimcycle, mythid)
380     CMM & 1, paramFld2d, 1, optimcycle, mythid)
381     CMM)
382     #endif
383    
384    
385     end
386    
387    
388    
389    
390    
391    

  ViewVC Help
Powered by ViewVC 1.1.22