/[MITgcm]/MITgcm/pkg/dic/dic_aver_init.F
ViewVC logotype

Diff of /MITgcm/pkg/dic/dic_aver_init.F

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

revision 1.8 by dfer, Mon Apr 7 20:31:16 2008 UTC revision 1.13 by jmc, Sat Jan 2 23:07:39 2010 UTC
# Line 7  CStartOfInterFace Line 7  CStartOfInterFace
7        SUBROUTINE DIC_AVER_INIT(        SUBROUTINE DIC_AVER_INIT(
8       I           myThid)       I           myThid)
9    
10  C     /==========================================================\  C     *==========================================================*
11  C     | SUBROUTINE DIC_AVER_INIT  i                            |  C     | SUBROUTINE DIC_AVER_INIT
12  C     |==========================================================|  C     *==========================================================*
13        IMPLICIT NONE        IMPLICIT NONE
14    
15  C     == GLobal variables ==  C     == GLobal variables ==
# Line 19  C     == GLobal variables == Line 19  C     == GLobal variables ==
19  #include "PARAMS.h"  #include "PARAMS.h"
20  #include "GRID.h"  #include "GRID.h"
21  #include "PTRACERS_SIZE.h"  #include "PTRACERS_SIZE.h"
 #include "PTRACERS_PARAMS.h"  
22  #include "PTRACERS_FIELDS.h"  #include "PTRACERS_FIELDS.h"
 #include "GCHEM.h"  
23  #include "DIC_VARS.h"  #include "DIC_VARS.h"
24  #ifdef DIC_BIOTIC  #ifdef DIC_BIOTIC
25  #include "DIC_DIAGS.h"  #include "DIC_DIAGS.h"
26  #include "DIC_COST.h"  #include "DIC_COST.h"
27  #endif  #endif
 c#ifdef ALLOW_SEAICE  
 c#include "SEAICE.h"  
 c#endif  
28    
29  C     == Routine arguments ==  C     == Routine arguments ==
30        INTEGER myThid        INTEGER myThid
31    
32  #ifdef ALLOW_DIC_COST  #ifdef ALLOW_DIC_COST
33    #ifdef ALLOW_TIMEAVE
34    
35  C     == Local variables ==  C     == Local variables ==
36        INTEGER i, j, bi, bj, k, it        INTEGER i, j, bi, bj, k, it
37        _RL po4av(nR)        _RL po4av(Nr)
38        _RL o2av(nR)        _RL o2av(Nr)
39        _RL volvar(nR)        _RL volvar(Nr)
40  cswdmonth -add-  Cswdmonth -add-
41        _RL po4avm(12,4)        _RL po4avm(12,4)
42        _RL o2avm(12,4)        _RL o2avm(12,4)
43        _RL rdt  Cswddmonth -- end-
       INTEGER nForcingPeriods,Imytm,Ifprd,Ifcyc,Iftm  
44    
45  cswddmonth -- end-  C initialize to zero
46  c  cph        totcost=0. _d 0
 c initialize to zero  
         totcost=0.d0  
47          DO bj = myByLo(myThid), myByHi(myThid)          DO bj = myByLo(myThid), myByHi(myThid)
48           DO bi = myBxLo(myThid), myBxHi(myThid)           DO bi = myBxLo(myThid), myBxHi(myThid)
49            CALL TIMEAVE_RESET(PO4obs,   Nr,  bi, bj, myThid)            CALL TIMEAVE_RESET(PO4obs,   Nr,  bi, bj, myThid)
50            CALL TIMEAVE_RESET(O2obs,   Nr,  bi, bj, myThid)            CALL TIMEAVE_RESET(O2obs,   Nr,  bi, bj, myThid)
51  cswdmonth  Cswdmonth
52            CALL TIMEAVE_RESET(PO4obsl1,   Nr,  bi, bj, myThid)            CALL TIMEAVE_RESET(PO4obsl1,   Nr,  bi, bj, myThid)
53            CALL TIMEAVE_RESET(PO4obsl2,   Nr,  bi, bj, myThid)            CALL TIMEAVE_RESET(PO4obsl2,   Nr,  bi, bj, myThid)
54            CALL TIMEAVE_RESET(PO4obsl3,   Nr,  bi, bj, myThid)            CALL TIMEAVE_RESET(PO4obsl3,   Nr,  bi, bj, myThid)
# Line 64  cQQ       CALL TIMEAVE_RESET(PO4obsl4, Line 57  cQQ       CALL TIMEAVE_RESET(PO4obsl4,
57            CALL TIMEAVE_RESET(O2obsl2,   Nr,  bi, bj, myThid)            CALL TIMEAVE_RESET(O2obsl2,   Nr,  bi, bj, myThid)
58            CALL TIMEAVE_RESET(O2obsl3,   Nr,  bi, bj, myThid)            CALL TIMEAVE_RESET(O2obsl3,   Nr,  bi, bj, myThid)
59  cQQ       CALL TIMEAVE_RESET(O2obsl4,   Nr,  bi, bj, myThid)  cQQ       CALL TIMEAVE_RESET(O2obsl4,   Nr,  bi, bj, myThid)
60  cswdmonth -end-  Cswdmonth -end-
61              OBS_timetave(bi,bj) = 0. _d 0
62            do k=1,Nr            do k=1,Nr
63              OBS_Timetave(bi,bj,k)=0.d0              po4av(k)  = 0. _d 0
64              po4av(k)=0.d0              o2av(k)   = 0. _d 0
65              o2av(k)=0.d0              po4var(k) = 0. _d 0
66              po4var(k)=0.d0              o2var(k)  = 0. _d 0
67              o2var(k)=0.d0              volvar(k) = 0. _d 0
             volvar(k)=0.d0  
68            enddo            enddo
69  cswdmonth  Cswdmonth
70            do k=1,3            do k=1,3
71             do it=1,12             do it=1,12
72              OBSM_Timetave(bi,bj,it)=0.d0              OBSM_Timetave(it,bi,bj) = 0. _d 0
73              po4avm(it,k)=0.d0              po4avm(it,k) = 0. _d 0
74              o2avm(it,k)=0.d0              o2avm(it,k)  = 0. _d 0
75              po4varm(it,k)=0.d0              po4varm(it,k)= 0. _d 0
76              o2varm(it,k)=0.d0              o2varm(it,k) = 0. _d 0
77             enddo             enddo
78            enddo            enddo
79           ENDDO           ENDDO
80          ENDDO          ENDDO
         _BEGIN_MASTER( myThid )  
81          CALL READ_FLD_XYZ_RL( 'input/po4obs.bin', ' ',          CALL READ_FLD_XYZ_RL( 'input/po4obs.bin', ' ',
82       &                          po4obs, 0, myThid )       &                          po4obs, 0, myThid )
83          CALL READ_FLD_XYZ_RL( 'input/o2obs.bin', ' ',          CALL READ_FLD_XYZ_RL( 'input/o2obs.bin', ' ',
84       &                          o2obs, 0, myThid )       &                          o2obs, 0, myThid )
85  cswdmonth  Cswdmonth
86          CALL READ_FLD_XYZ_RL( 'input/po4lev1.bin', ' ',          CALL READ_FLD_XYZ_RL( 'input/po4lev1.bin', ' ',
87       &                          po4obsl1, 0, myThid )       &                          po4obsl1, 0, myThid )
88          CALL READ_FLD_XYZ_RL( 'input/po4lev2.bin', ' ',          CALL READ_FLD_XYZ_RL( 'input/po4lev2.bin', ' ',
# Line 107  cQQ  &                          po4obsl4 Line 99  cQQ  &                          po4obsl4
99       &                          o2obsl3, 0, myThid )       &                          o2obsl3, 0, myThid )
100  cQQ     CALL READ_FLD_XYZ_RL( 'input/o2lev4.bin', ' ',  cQQ     CALL READ_FLD_XYZ_RL( 'input/o2lev4.bin', ' ',
101  cQQ  &                          o2obsl4, 0, myThid )  cQQ  &                          o2obsl4, 0, myThid )
102  cswdmonth -end-  Cswdmonth -end-
103         _END_MASTER(myThid)         _EXCH_XYZ_RL(po4obs  , myThid )
104         _EXCH_XYZ_R8(po4obs  , myThid )         _EXCH_XYZ_RL(o2obs  , myThid )
105         _EXCH_XYZ_R8(o2obs  , myThid )  Cswdmonth -add-
106  cswdmonth -add-         _EXCH_XYZ_RL(po4obsl1  , myThid )
107         _EXCH_XYZ_R8(po4obsl1  , myThid )         _EXCH_XYZ_RL(po4obsl2  , myThid )
108         _EXCH_XYZ_R8(po4obsl2  , myThid )         _EXCH_XYZ_RL(po4obsl3  , myThid )
109         _EXCH_XYZ_R8(po4obsl3  , myThid )  cQQ    _EXCH_XYZ_RL(po4obsl4  , myThid )
110  cQQ    _EXCH_XYZ_R8(po4obsl4  , myThid )         _EXCH_XYZ_RL(o2obsl1  , myThid )
111         _EXCH_XYZ_R8(o2obsl1  , myThid )         _EXCH_XYZ_RL(o2obsl2  , myThid )
112         _EXCH_XYZ_R8(o2obsl2  , myThid )         _EXCH_XYZ_RL(o2obsl3  , myThid )
113         _EXCH_XYZ_R8(o2obsl3  , myThid )  cQQ    _EXCH_XYZ_RL(o2obsl4  , myThid )
114  cQQ    _EXCH_XYZ_R8(o2obsl4  , myThid )  Cswdmonth -end-
115  cswdmonth -end-  
116    C calculate layer means
         _BARRIER  
 c calculate layer means  
         _BEGIN_MASTER( mythid )  
117          do k=1,Nr          do k=1,Nr
118           call tracer_meanarea(myThid,po4obs, k,           call tracer_meanarea(po4obs, k,
119       &                    po4av(k))       &                    po4av(k),myThid)
120           call tracer_meanarea(myThid,o2obs, k,           call tracer_meanarea(o2obs, k,
121       &                    o2av(k))       &                    o2av(k),myThid)
122  c        print*,po4av(k), o2av(k)  c        print*,po4av(k), o2av(k)
123          enddo          enddo
124  cswdmonth -add-  Cswdmonth -add-
125          do it=1,12          do it=1,MIN(12,Nr)
126           call tracer_meanarea(myThid,po4obsl1,it,           call tracer_meanarea(po4obsl1,it,
127       &                    po4avm(it,1))       &                    po4avm(it,1),myThid)
128           call tracer_meanarea(myThid,po4obsl2,it,           call tracer_meanarea(po4obsl2,it,
129       &                    po4avm(it,2))       &                    po4avm(it,2),myThid)
130           call tracer_meanarea(myThid,po4obsl3,it,           call tracer_meanarea(po4obsl3,it,
131       &                    po4avm(it,3))       &                    po4avm(it,3),myThid)
132  cQQ      call tracer_meanarea(myThid,po4obsl4,it,  cQQ      call tracer_meanarea(po4obsl4,it,
133  cQQ  &                    po4avm(it,4))  cQQ  &                    po4avm(it,4),myThid)
134           call tracer_meanarea(myThid,o2obsl1,it,           call tracer_meanarea(o2obsl1,it,
135       &                    o2avm(it,1))       &                    o2avm(it,1),myThid)
136           call tracer_meanarea(myThid,o2obsl2,it,           call tracer_meanarea(o2obsl2,it,
137       &                    o2avm(it,2))       &                    o2avm(it,2),myThid)
138           call tracer_meanarea(myThid,o2obsl3,it,           call tracer_meanarea(o2obsl3,it,
139       &                    o2avm(it,3))       &                    o2avm(it,3),myThid)
140  cQQ      call tracer_meanarea(myThid,o2obsl4,it,  cQQ      call tracer_meanarea(o2obsl4,it,
141  cQQ  &                    o2avm(it,4))  cQQ  &                    o2avm(it,4),myThid)
142    
143          enddo          enddo
         _END_MASTER(myThid)  
144  c calculate layer variance  c calculate layer variance
         _BEGIN_MASTER( mythid )  
145          DO bj = myByLo(myThid), myByHi(myThid)          DO bj = myByLo(myThid), myByHi(myThid)
146          DO bi = myBxLo(myThid), myBxHi(myThid)          DO bi = myBxLo(myThid), myBxHi(myThid)
147           DO j=1-OLy,sNy+OLy           DO j=1,sNy
148           DO i=1-OLx,sNx+OLx           DO i=1,sNx
149            DO k=1,Nr            DO k=1,Nr
150              volvar(k)=volvar(k)+              volvar(k)=volvar(k)+
151       &                rA(i,j,bi,bj)*drF(k)*maskC(i,j,k,bi,bj)       &                rA(i,j,bi,bj)*drF(k)*maskC(i,j,k,bi,bj)
# Line 167  c calculate layer variance Line 154  c calculate layer variance
154              o2var(k)=o2var(k)+(o2obs(i,j,k,bi,bj)-o2av(k))**2              o2var(k)=o2var(k)+(o2obs(i,j,k,bi,bj)-o2av(k))**2
155       &                *rA(i,j,bi,bj)*drF(k)*maskC(i,j,k,bi,bj)       &                *rA(i,j,bi,bj)*drF(k)*maskC(i,j,k,bi,bj)
156            ENDDO            ENDDO
157  cswdmonth -add-  Cswdmonth -add-
158            DO it=1,12            DO it=1,12
159             po4varm(it,1)=po4varm(it,1)+             po4varm(it,1)=po4varm(it,1)+
160       &                (po4obsl1(i,j,it,bi,bj)-po4avm(it,1))**2       &                (po4obsl1(i,j,it,bi,bj)-po4avm(it,1))**2
# Line 204  cQQ  &                *rA(i,j,bi,bj)*drF Line 191  cQQ  &                *rA(i,j,bi,bj)*drF
191              o2var(k)=o2var(k)/volvar(k)              o2var(k)=o2var(k)/volvar(k)
192  cQQ         print*,po4var(k),o2var(k)  cQQ         print*,po4var(k),o2var(k)
193          ENDDO          ENDDO
194  cswdmonth- add-  Cswdmonth- add-
195          DO k=1,3          DO k=1,3
196           Do it=1,12           Do it=1,12
197             po4varm(it,k)=po4varm(it,k)/volvar(k)             po4varm(it,k)=po4varm(it,k)/volvar(k)
198             o2varm(it,k)=o2varm(it,k)/volvar(k)             o2varm(it,k)=o2varm(it,k)/volvar(k)
199           ENDDO           ENDDO
200          ENDDO          ENDDO
201  cswdmonth -end-  Cswdmonth -end-
202          _END_MASTER(myThid)  
 C  
203  C Reset averages to zero  C Reset averages to zero
204         print*,'QQ dic_diags, set to zero, gchem_init'         print*,'QQ dic_aver_init, set to zero'
205         DO bj = myByLo(myThid), myByHi(myThid)         DO bj = myByLo(myThid), myByHi(myThid)
206          DO bi = myBxLo(myThid), myBxHi(myThid)          DO bi = myBxLo(myThid), myBxHi(myThid)
207           CALL TIMEAVE_RESET(PO4ann,Nr,bi,bj,myThid)           CALL TIMEAVE_RESET(PO4ann,Nr,bi,bj,myThid)
# Line 229  cQQ       CALL TIMEAVE_RESET(PO4lev4, Line 215  cQQ       CALL TIMEAVE_RESET(PO4lev4,
215           CALL TIMEAVE_RESET(O2lev3,    12,  bi, bj, myThid)           CALL TIMEAVE_RESET(O2lev3,    12,  bi, bj, myThid)
216  cQQ       CALL TIMEAVE_RESET(O2lev4,    12,  bi, bj, myThid)  cQQ       CALL TIMEAVE_RESET(O2lev4,    12,  bi, bj, myThid)
217    
218           do k=1,Nr           OBS_timetave(bi,bj) = 0. _d 0
            OBS_Timetave(bi,bj,k)=0.d0  
          enddo  
219           do it=1,12           do it=1,12
220             OBSM_Timetave(bi,bj,it)=0.d0             OBSM_Timetave(it,bi,bj) = 0. _d 0
221           enddo           enddo
222          ENDDO          ENDDO
223         ENDDO         ENDDO
224  c  
225  #endif  #endif /* ALLOW_TIMEAVE */
226  c  #endif /* ALLOW_DIC_COST */
227    
228        RETURN        RETURN
229        END        END
 cswd -- end added subroutine --  
Legend:
Removed from v.1.8  
changed lines
  Added in v.1.13
  ViewVC Help
Powered by ViewVC 1.1.22