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

Annotation of /MITgcm_contrib/SOSE/code_ad/ctrl_get_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_get_gen.F,v 1.13 2007/06/19 03:42:30 gforget Exp $
2     C $Name: $
3    
4     #include "CTRL_CPPOPTIONS.h"
5    
6    
7     subroutine ctrl_get_gen(
8     I xx_gen_file, xx_genstartdate, xx_genperiod,
9     I genmask, genfld, xx_gen0, xx_gen1, xx_gen_dummy,
10     I xx_gen_remo_intercept, xx_gen_remo_slope,
11     I mytime, myiter, mythid
12     & )
13    
14     c ==================================================================
15     c SUBROUTINE ctrl_get_gen
16     c ==================================================================
17     c
18     c o new generic routine for reading time dependent control variables
19     c heimbach@mit.edu 12-Jun-2003
20     c
21     c ==================================================================
22     c SUBROUTINE ctrl_get_gen
23     c ==================================================================
24    
25     implicit none
26    
27     c == global variables ==
28    
29     #include "EEPARAMS.h"
30     #include "SIZE.h"
31     #include "PARAMS.h"
32     #include "GRID.h"
33    
34     #include "ctrl.h"
35     #include "ctrl_dummy.h"
36     #include "optim.h"
37     #ifdef ALLOW_EXF
38     # include "EXF_FIELDS.h"
39     #endif
40    
41     c == routine arguments ==
42    
43     character*(80) fnamegeneric
44     character*(MAX_LEN_FNAM) xx_gen_file
45     integer xx_genstartdate(4)
46     _RL xx_genperiod
47     _RL genmask(1-olx:snx+olx,1-oly:sny+oly,nr,nsx,nsy)
48     _RL genfld(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy)
49     _RL xx_gen0(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy)
50     _RL xx_gen1(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy)
51     _RL xx_gen_dummy
52     _RL xx_gen_remo_intercept
53     _RL xx_gen_remo_slope
54    
55     _RL mytime
56     integer myiter
57     integer mythid
58    
59     c == local variables ==
60    
61     #ifdef ALLOW_EXF
62    
63     integer bi,bj
64     integer i,j,k
65     integer itlo,ithi
66     integer jtlo,jthi
67     integer jmin,jmax
68     integer imin,imax
69     integer ilgen
70    
71     _RL gensign
72     _RL genfac
73     logical doCtrlUpdate
74     logical genfirst
75     logical genchanged
76     integer gencount0
77     integer gencount1
78    
79     logical doglobalread
80     logical ladinit
81    
82     character*(80) fnamegen
83    
84     c == external functions ==
85    
86     integer ilnblnk
87     external ilnblnk
88    
89    
90     c == end of interface ==
91    
92     jtlo = mybylo(mythid)
93     jthi = mybyhi(mythid)
94     itlo = mybxlo(mythid)
95     ithi = mybxhi(mythid)
96     jmin = 1-oly
97     jmax = sny+oly
98     imin = 1-olx
99     imax = snx+olx
100    
101     c-- Now, read the control vector.
102     doglobalread = .false.
103     ladinit = .false.
104    
105     if (optimcycle .ge. 0) then
106     ilgen=ilnblnk( xx_gen_file )
107     write(fnamegen(1:80),'(2a,i10.10)')
108     & xx_gen_file(1:ilgen), '.', optimcycle
109     endif
110    
111     # ifdef ALLOW_CAL
112     if ( xx_genperiod .EQ. 0 ) then
113     c record numbers are assumed 1 to 12 corresponding to
114     c Jan. through Dec.
115     call cal_GetMonthsRec(
116     O genfac, genfirst, genchanged,
117     O gencount0, gencount1,
118     I mytime, myiter, mythid
119     & )
120     else
121     c-- Get the counters, flags, and the interpolation factor.
122     call ctrl_get_gen_rec(
123     I xx_genstartdate, xx_genperiod,
124     O genfac, genfirst, genchanged,
125     O gencount0,gencount1,
126     I mytime, myiter, mythid )
127     endif
128     # else
129     c-- Get the counters, flags, and the interpolation factor.
130     call ctrl_get_gen_rec(
131     I xx_genstartdate, xx_genperiod,
132     O genfac, genfirst, genchanged,
133     O gencount0,gencount1,
134     I mytime, myiter, mythid )
135     # endif
136    
137     if ( genfirst ) then
138     call active_read_xy( fnamegen, xx_gen1, gencount0,
139     & doglobalread, ladinit, optimcycle,
140     & mythid, xx_gen_dummy )
141     #ifdef ALLOW_CTRL_SMOOTH
142     if ( xx_gen_file .EQ. xx_tauu_file .OR.
143     & xx_gen_file .EQ. xx_tauv_file )
144     & call ctrl_smooth(xx_gen1,genmask)
145     #endif
146     #ifdef ALLOW_SMOOTH_CORREL2D
147     CMM(
148     call smooth_correl2Dw(xx_gen1,genmask,xx_gen_file,1,mythid)
149     call smooth_correl2D(xx_gen1,genmask,1,mythid)
150     call smooth_correl2Dw(xx_gen1,genmask,xx_gen_file,2,mythid)
151     CMM)
152     write(fnamegeneric(1:80),'(2a,i10.10)')
153     & xx_gen_file(1:ilgen),'.effective.',optimcycle
154     call mdswritefield(fnamegeneric,ctrlprec,.FALSE.,'RL',
155     & 1, xx_gen1, gencount0, optimcycle, mythid)
156    
157     CMM(
158     CMM call smooth_correl2Dw(xx_gen1,genmask,xx_gen_file,2,mythid)
159     CMM)
160     #endif
161     endif
162     if (( genfirst ) .or. ( genchanged )) then
163     call exf_SwapFFields( xx_gen0, xx_gen1, mythid )
164    
165     call active_read_xy( fnamegen, xx_gen1 , gencount1,
166     & doglobalread, ladinit, optimcycle,
167     & mythid, xx_gen_dummy )
168     #ifdef ALLOW_CTRL_SMOOTH
169     if ( xx_gen_file .EQ. xx_tauu_file .OR.
170     & xx_gen_file .EQ. xx_tauv_file )
171     & call ctrl_smooth(xx_gen1,genmask)
172     #endif
173     #ifdef ALLOW_SMOOTH_CORREL2D
174     CMM(
175     call smooth_correl2Dw(xx_gen1,genmask,xx_gen_file,1,mythid)
176     call smooth_correl2D(xx_gen1,genmask,1,mythid)
177     call smooth_correl2Dw(xx_gen1,genmask,xx_gen_file,2,mythid)
178     CMM)
179     write(fnamegeneric(1:80),'(2a,i10.10)')
180     & xx_gen_file(1:ilgen),'.effective.',optimcycle
181     call mdswritefield(fnamegeneric,ctrlprec,.FALSE.,'RL',
182     & 1, xx_gen1, gencount1, optimcycle, mythid)
183    
184     CMM(
185     CMM call smooth_correl2Dw(xx_gen1,genmask,xx_gen_file,2,mythid)
186     CMM)
187     #endif
188     endif
189    
190     c-- Add control to model variable.
191     cph(
192     cph this flag ported from the SIO code
193     cph Initial wind stress adjustments are too vigorous.
194     if ( gencount0 .LE. 2 .AND.
195     & ( xx_gen_file .EQ. xx_tauu_file .OR.
196     & xx_gen_file .EQ. xx_tauv_file ) .AND.
197     & ( xx_genperiod .NE. 0 ) ) then
198     doCtrlUpdate = .FALSE.
199     else
200     doCtrlUpdate = .TRUE.
201     endif
202     if ( xx_gen_file .EQ. xx_tauu_file .OR.
203     & xx_gen_file .EQ. xx_tauv_file ) then
204     gensign = -1.
205     else
206     gensign = 1.
207     endif
208     c
209     cph since the above is ECCO specific, we undo it here:
210     cph doCtrlUpdate = .TRUE.
211     c
212     if ( doCtrlUpdate ) then
213     cph)
214     do bj = jtlo,jthi
215     do bi = itlo,ithi
216     c-- Calculate mask for tracer cells (0 => land, 1 => water).
217     k = 1
218     do j = 1,sny
219     do i = 1,snx
220     genfld(i,j,bi,bj) = genfld (i,j,bi,bj)
221     & + gensign*genfac *xx_gen0(i,j,bi,bj)
222     & + gensign*(1. _d 0 - genfac)*xx_gen1(i,j,bi,bj)
223     genfld(i,j,bi,bj) =
224     & genmask(i,j,k,bi,bj)*( genfld (i,j,bi,bj) -
225     & ( xx_gen_remo_intercept +
226     & xx_gen_remo_slope*(mytime-starttime) ) )
227     enddo
228     enddo
229     enddo
230     enddo
231     cph(
232     endif
233     cph)
234    
235     #endif /* ALLOW_EXF */
236    
237     end
238    

  ViewVC Help
Powered by ViewVC 1.1.22