Back to home page

MITgcm
 
 

    

File indexing completed on 2018-03-02 18:39:05 UTC

view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
b614ad87ad Jean*0001 #include "DIAG_OPTIONS.h"
                0002 
                0003 CBOP
                0004 C     !ROUTINE: DIAGSTATS_CALC
                0005 C     !INTERFACE:
                0006       SUBROUTINE DIAGSTATS_CALC(
                0007      O                  statArr,
                0008      I                  inpArr, frcArr, scaleFact, power, useFract,
0ca9894620 Jean*0009      I                  useReg, regMskVal,
b614ad87ad Jean*0010      I                  nStats,sizI1,sizI2,sizJ1,sizJ2, iRun,jRun,
                0011      I                  regMask, arrMask, arrhFac, arrArea,
                0012      I                  arrDr, specialVal, exclSpVal, useWeight,
                0013      I                  myThid )
                0014 
                0015 C     !DESCRIPTION:
                0016 C     Compute statistics for this tile, level, region
                0017 
                0018 C     !USES:
                0019       IMPLICIT NONE
                0020 
                0021 #include "EEPARAMS.h"
                0022 #include "SIZE.h"
                0023 
                0024 C     !INPUT/OUTPUT PARAMETERS:
                0025 C     == Routine Arguments ==
                0026 C     statArr     :: output statistics array
                0027 C     inpArr      :: input field array to process (compute stats & add to statFld)
                0028 C     frcArr      :: fraction used for weighted-average diagnostics
                0029 C     scaleFact   :: scaling factor
                0030 C     power       :: option to fill-in with the field square (power=2)
                0031 C     useFract    :: if True, use fraction-weight
0ca9894620 Jean*0032 C     useReg      :: how to use region-mask: =0 : not used ;
                0033 C                    =1 : grid-center location ; =2 : U location ; =3 : V location
b614ad87ad Jean*0034 C     regMskVal   :: region-mask identificator value
0ca9894620 Jean*0035 C                    (point i,j belong to region <=> regMask(i,j) = regMskVal)
b614ad87ad Jean*0036 C     nStats      :: size of output array: statArr
                0037 C     sizI1,sizI2 :: size of inpArr array: 1rst index range (min,max)
                0038 C     sizJ1,sizJ2 :: size of inpArr array: 2nd  index range (min,max)
                0039 C     iRun,jRun   :: range of 1rst & 2nd index to process
                0040 C     regMask     :: regional mask
                0041 C     arrMask     :: mask for this input array
                0042 C     arrhFac     :: weight factor (horizontally varying)
                0043 C     arrArea     :: Area weighting factor
                0044 C     arrDr       :: uniform weighting factor
                0045 C     specialVal  :: special value in input array (to exclude if exclSpVal=T)
                0046 C     exclSpVal   :: if T, exclude "specialVal" in input array
                0047 C     useWeight   :: use weight factor "arrhFac"
                0048 Cc    k,bi,bj     :: level and tile indices used for weighting (mask,area ...)
                0049 Cc    parsFld     :: parser field with characteristics of the diagnostics
                0050 C     myThid      :: my Thread Id number
                0051       INTEGER nStats,sizI1,sizI2,sizJ1,sizJ2
                0052       INTEGER iRun, jRun
                0053       _RL statArr(0:nStats)
                0054       _RL inpArr (sizI1:sizI2,sizJ1:sizJ2)
                0055       _RL frcArr (sizI1:sizI2,sizJ1:sizJ2)
                0056       _RL scaleFact
                0057       INTEGER power
                0058       LOGICAL useFract
0ca9894620 Jean*0059       INTEGER useReg
b614ad87ad Jean*0060       _RS regMskVal
                0061       _RS regMask(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0062       _RS arrMask(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0063       _RS arrhFac(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0064       _RS arrArea(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0065       _RL arrDr
                0066       _RL specialVal
                0067       LOGICAL exclSpVal
                0068       LOGICAL useWeight
                0069       INTEGER myThid
                0070 CEOP
                0071 
                0072 C     !LOCAL VARIABLES:
                0073 C     i,j    :: loop indices
                0074       INTEGER i, j, n
                0075       INTEGER im, ix
1079a52e0a Mart*0076 #ifndef TARGET_NEC_SX
0ca9894620 Jean*0077       LOGICAL inside(sNx+1,sNy+1)
                0078       _RL     tmpFld(sNx+1,sNy+1)
                0079       _RL     tmpVol(sNx+1,sNy+1)
1079a52e0a Mart*0080 #else
                0081 C     Extra variables and fields to support vectorization.
                0082 C     This code also uses the intrinsic F90 routines SUM, MINVAL, MAXVAL
                0083 C     and thus will not compile with a F77 compiler.
                0084       _RL     arrMaskL(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0085       _RL     tmpFld  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0086       _RL     tmpVol  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0087 #endif
b614ad87ad Jean*0088 
                0089       im = nStats - 1
                0090       ix = nStats
                0091       DO n=0,nStats
                0092         statArr(n) = 0.
                0093       ENDDO
                0094 
1079a52e0a Mart*0095 #ifndef TARGET_NEC_SX
b614ad87ad Jean*0096 
0ca9894620 Jean*0097 C-    Apply Scaling Factor and power option to Input Field (-> tmpFld):
                0098       IF ( power.EQ.2 ) THEN
b614ad87ad Jean*0099        DO j = 1,jRun
                0100         DO i = 1,iRun
0ca9894620 Jean*0101           tmpFld(i,j) = scaleFact*inpArr(i,j)
                0102           tmpFld(i,j) = tmpFld(i,j)*tmpFld(i,j)
b614ad87ad Jean*0103         ENDDO
                0104        ENDDO
0ca9894620 Jean*0105       ELSE
b614ad87ad Jean*0106        DO j = 1,jRun
                0107         DO i = 1,iRun
0ca9894620 Jean*0108           tmpFld(i,j) = scaleFact*inpArr(i,j)
b614ad87ad Jean*0109         ENDDO
                0110        ENDDO
0ca9894620 Jean*0111       ENDIF
b614ad87ad Jean*0112 
0ca9894620 Jean*0113 C-    Set weight factor "tmpVol" (area and hFac and/or fraction field)
                0114 C     and part of domain (=inside) where to compute stats
                0115       IF ( useFract .AND. useWeight ) THEN
b614ad87ad Jean*0116        DO j = 1,jRun
                0117         DO i = 1,iRun
0ca9894620 Jean*0118           inside(i,j) = arrMask(i,j).NE.0.
                0119      &            .AND. arrhFac(i,j).NE.0.
                0120      &            .AND. frcArr(i,j) .NE.0.
                0121           tmpVol(i,j) = arrArea(i,j)*arrhFac(i,j)*frcArr(i,j)
b614ad87ad Jean*0122         ENDDO
                0123        ENDDO
0ca9894620 Jean*0124       ELSEIF ( useFract ) THEN
b614ad87ad Jean*0125        DO j = 1,jRun
                0126         DO i = 1,iRun
0ca9894620 Jean*0127           inside(i,j) = arrMask(i,j).NE.0.
                0128      &            .AND. arrhFac(i,j).NE.0.
                0129      &            .AND. frcArr(i,j) .NE.0.
                0130           tmpVol(i,j) = arrArea(i,j)*frcArr(i,j)
b614ad87ad Jean*0131         ENDDO
                0132        ENDDO
0ca9894620 Jean*0133       ELSEIF ( useWeight ) THEN
b614ad87ad Jean*0134        DO j = 1,jRun
                0135         DO i = 1,iRun
0ca9894620 Jean*0136           inside(i,j) = arrMask(i,j).NE.0.
                0137      &            .AND. arrhFac(i,j).NE.0.
                0138           tmpVol(i,j) = arrArea(i,j)*arrhFac(i,j)
b614ad87ad Jean*0139         ENDDO
                0140        ENDDO
0ca9894620 Jean*0141       ELSE
b614ad87ad Jean*0142        DO j = 1,jRun
                0143         DO i = 1,iRun
0ca9894620 Jean*0144           inside(i,j) = arrMask(i,j).NE.0.
                0145      &            .AND. arrhFac(i,j).NE.0.
                0146           tmpVol(i,j) = arrArea(i,j)
b614ad87ad Jean*0147         ENDDO
                0148        ENDDO
0ca9894620 Jean*0149       ENDIF
b614ad87ad Jean*0150 
0ca9894620 Jean*0151 C-    Exclude (setting inside=F) Special Value:
                0152       IF ( exclSpVal ) THEN
b614ad87ad Jean*0153        DO j = 1,jRun
                0154         DO i = 1,iRun
0ca9894620 Jean*0155           inside(i,j) = inside(i,j) .AND. inpArr(i,j).NE.specialVal
b614ad87ad Jean*0156         ENDDO
                0157        ENDDO
0ca9894620 Jean*0158       ENDIF
                0159 C-    Account for Region-mask (refine "inside"):
                0160       IF ( useReg.EQ.1 ) THEN
b614ad87ad Jean*0161        DO j = 1,jRun
                0162         DO i = 1,iRun
0ca9894620 Jean*0163           inside(i,j) = inside(i,j) .AND. regMask(i,j).EQ.regMskVal
b614ad87ad Jean*0164         ENDDO
                0165        ENDDO
0ca9894620 Jean*0166       ELSEIF ( useReg.EQ.2 ) THEN
1079a52e0a Mart*0167        DO j = 1,jRun
                0168         DO i = 1,iRun
0ca9894620 Jean*0169           inside(i,j) = inside(i,j) .AND.( regMask(i,j).EQ.regMskVal
                0170      &                              .OR. regMask(i-1,j).EQ.regMskVal )
1079a52e0a Mart*0171         ENDDO
                0172        ENDDO
0ca9894620 Jean*0173       ELSEIF ( useReg.EQ.3 ) THEN
1079a52e0a Mart*0174        DO j = 1,jRun
                0175         DO i = 1,iRun
0ca9894620 Jean*0176           inside(i,j) = inside(i,j) .AND.( regMask(i,j).EQ.regMskVal
                0177      &                              .OR. regMask(i,j-1).EQ.regMskVal )
1079a52e0a Mart*0178         ENDDO
                0179        ENDDO
0ca9894620 Jean*0180       ENDIF
1079a52e0a Mart*0181 
0ca9894620 Jean*0182 C-    Calculate Stats
                0183       DO j = 1,jRun
                0184        DO i = 1,iRun
                0185         IF ( inside(i,j) ) THEN
                0186           statArr(im) = tmpFld(i,j)
                0187           statArr(0) = statArr(0) + tmpVol(i,j)
                0188           statArr(1) = statArr(1) + tmpVol(i,j)*tmpFld(i,j)
                0189           statArr(2) = statArr(2) + tmpVol(i,j)*tmpFld(i,j)*tmpFld(i,j)
                0190         ENDIF
1079a52e0a Mart*0191        ENDDO
0ca9894620 Jean*0192       ENDDO
                0193       statArr(ix) = statArr(im)
                0194       DO j = 1,jRun
                0195        DO i = 1,iRun
                0196         IF ( inside(i,j) ) THEN
                0197           statArr(im) = MIN(tmpFld(i,j),statArr(im))
                0198           statArr(ix) = MAX(tmpFld(i,j),statArr(ix))
                0199         ENDIF
                0200        ENDDO
                0201       ENDDO
                0202       statArr(0) = statArr(0)*arrDr
                0203       statArr(1) = statArr(1)*arrDr
                0204       statArr(2) = statArr(2)*arrDr
1079a52e0a Mart*0205 
0ca9894620 Jean*0206 #else /* TARGET_NEC_SX defined */
1079a52e0a Mart*0207 
0ca9894620 Jean*0208       arrMaskL = 0. _d 0
1079a52e0a Mart*0209 
0ca9894620 Jean*0210       IF ( useFract .AND. exclSpVal ) THEN
1079a52e0a Mart*0211 
                0212        DO j = 1,jRun
                0213         DO i = 1,iRun
                0214          IF ( arrMask(i,j).NE.0. .AND. frcArr(i,j).NE.0.
0ca9894620 Jean*0215      &             .AND. arrhFac(i,j).NE.0.
                0216      &             .AND. inpArr(i,j).NE.specialVal )
1079a52e0a Mart*0217      &        arrMaskL(i,j) = 1. _d 0
                0218         ENDDO
                0219        ENDDO
                0220        IF ( useWeight ) THEN
3a14199d51 Mart*0221         tmpVol = arrhFac*arrArea*frcArr
1079a52e0a Mart*0222        ELSE
3a14199d51 Mart*0223         tmpVol = arrArea*frcArr
1079a52e0a Mart*0224        ENDIF
                0225 
                0226       ELSEIF ( useFract ) THEN
                0227 
                0228        DO j = 1,jRun
                0229         DO i = 1,iRun
                0230          IF ( arrMask(i,j).NE.0. .AND. frcArr(i,j).NE.0.
0ca9894620 Jean*0231      &             .AND. arrhFac(i,j).NE.0. )
1079a52e0a Mart*0232      &        arrMaskL(i,j) = 1. _d 0
                0233         ENDDO
                0234        ENDDO
                0235        IF ( useWeight ) THEN
3a14199d51 Mart*0236         tmpVol = arrhFac*arrArea*frcArr
1079a52e0a Mart*0237        ELSE
3a14199d51 Mart*0238         tmpVol = arrArea*frcArr
1079a52e0a Mart*0239        ENDIF
                0240 
                0241       ELSEIF ( exclSpVal ) THEN
                0242 
                0243        DO j = 1,jRun
                0244         DO i = 1,iRun
                0245          IF ( arrMask(i,j).NE.0.
0ca9894620 Jean*0246      &             .AND. arrhFac(i,j).NE.0.
                0247      &             .AND. inpArr(i,j).NE.specialVal )
1079a52e0a Mart*0248      &        arrMaskL(i,j) = 1. _d 0
                0249         ENDDO
                0250        ENDDO
                0251        IF ( useWeight ) THEN
3a14199d51 Mart*0252         tmpVol = arrhFac*arrArea
1079a52e0a Mart*0253        ELSE
3a14199d51 Mart*0254         tmpVol = arrArea
1079a52e0a Mart*0255        ENDIF
                0256 
                0257       ELSE
                0258 
                0259        DO j = 1,jRun
                0260         DO i = 1,iRun
                0261          IF ( arrMask(i,j).NE.0.
0ca9894620 Jean*0262      &             .AND. arrhFac(i,j).NE.0. )
1079a52e0a Mart*0263      &        arrMaskL(i,j) = 1. _d 0
                0264         ENDDO
                0265        ENDDO
                0266        IF ( useWeight ) THEN
3a14199d51 Mart*0267         tmpVol = arrhFac*arrArea
1079a52e0a Mart*0268        ELSE
3a14199d51 Mart*0269         tmpVol = arrArea
1079a52e0a Mart*0270        ENDIF
                0271 
                0272       ENDIF
0ca9894620 Jean*0273 
                0274 C-    Account for Region-mask:
                0275       IF ( useReg.EQ.1 ) THEN
                0276        DO j = 1,jRun
                0277         DO i = 1,iRun
                0278           IF ( regMask(i,j).NE.regMskVal ) arrMaskL(i,j) = 0. _d 0
                0279         ENDDO
                0280        ENDDO
                0281       ELSEIF ( useReg.EQ.2 ) THEN
                0282        DO j = 1,jRun
                0283         DO i = 1,iRun
                0284           IF ( regMask(i,j).NE.regMskVal .AND.
                0285      &       regMask(i-1,j).NE.regMskVal ) arrMaskL(i,j) = 0. _d 0
                0286         ENDDO
                0287        ENDDO
                0288       ELSEIF ( useReg.EQ.3 ) THEN
                0289        DO j = 1,jRun
                0290         DO i = 1,iRun
                0291           IF ( regMask(i,j).NE.regMskVal .AND.
                0292      &       regMask(i,j-1).NE.regMskVal ) arrMaskL(i,j) = 0. _d 0
                0293         ENDDO
                0294        ENDDO
                0295       ENDIF
                0296 
53f764b624 Mart*0297 C     inpArr can be undefined/non-initialised in overlaps, so we need
                0298 C     to clean this fields first by copying the defined range to tmpFld
                0299       tmpFld = 0. _d 0
                0300       DO j = 1,jRun
                0301        DO i = 1,iRun
0ca9894620 Jean*0302         tmpFld(i,j) = inpArr(i,j)*scaleFact
53f764b624 Mart*0303        ENDDO
                0304       ENDDO
1079a52e0a Mart*0305       IF ( power.EQ.2) THEN
53f764b624 Mart*0306        tmpFld = tmpFld*tmpFld
1079a52e0a Mart*0307       ENDIF
                0308 C     sum up the volume
                0309       tmpVol = tmpVol*arrMaskL
3a14199d51 Mart*0310       statArr(0)  = SUM(tmpVol)*arrDr
1079a52e0a Mart*0311 C     compute and sum up volume*field
                0312       tmpVol = tmpVol*tmpFld
3ad53f77c9 Mart*0313       statArr(1)  = SUM(tmpVol)*arrDr
1079a52e0a Mart*0314 C     compute and sum up volume*field**2
                0315       tmpVol = tmpVol*tmpFld
3ad53f77c9 Mart*0316       statArr(2)  = SUM(tmpVol)*arrDr
                0317       statArr(im) = MINVAL(tmpFld, MASK = arrMaskL>0.)
                0318       statArr(ix) = MAXVAL(tmpFld, MASK = arrMaskL>0.)
1079a52e0a Mart*0319 
                0320 #endif /* TARGET_NEC_SX */
                0321 
b614ad87ad Jean*0322       RETURN
                0323       END
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated from https://github.com/MITgcm/MITgcm by the 2.2.1-MITgcm-0.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!