Changeset 2214 in lmdz_wrf

Timestamp:

Nov 5, 2018, 7:16:03 PM (6 years ago)

Author:

lfita

Message:

Introducing multiple ranges in `ranges_faces'

Location:

trunk/tools

Files:

• diag_tools.py (modified) (3 diffs)
• diagnostics.py (modified) (2 diffs)
• module_ForDiagnostics.f90 (modified) (6 diffs)
• module_ForDiagnosticsVars.f90 (modified) (5 diffs)

trunk/tools/diag_tools.py

@@ -1394 +1394 @@
         dy = hgt.shape[0]
 
-        hgtmaxt, pthgtmaxt, dhgtt, peakst, valleyst, ofacest, ffacest =              \
+        hgtmaxt, pthgtmaxt, dhgtt, peakst, valleyst, ofacest, ffacest, rnghgtmaxt,   \
+          ptrnghgtmaxt =                                                             \
           fdin.module_fordiagnostics.compute_range_faces(lon=lon[:].transpose(),     \
           lat=lat[:].transpose(), hgt=hgt[:].transpose(), face=face, nfilt=Nfilt,    \
@@ -1405 +1406 @@
         origfaces = ofacest.transpose()
         filtfaces = ffacest.transpose()
+        rangeshgtmax = rnghgtmaxt.transpose()
+        ptrangeshgtmax = ptrnghgtmaxt.transpose()
     else:
         print errormsg
@@ -1412 +1415 @@
 
     return hgtmax, pthgtmax, dhgt, peaks, valleys, origfaces, filtfaces, vardims,    \
-      varvdims
+      varvdims, rangeshgtmax, ptrangeshgtmax
 
 ####### ###### ##### #### ### ## # END Fortran diagnostics

trunk/tools/diagnostics.py

@@ -704 +704 @@
 
         orogmax, ptorogmax, dhgt, peaks, valleys, origfaces, diagout, diagoutd,      \
-          diagoutvd = diag.Forcompute_range_faces(lon, lat, hgt, face, Nptfilt,      \
-          Nptnewrange, dnamesvar, dvnamesvar)
+          diagoutvd, rngorogmax, ptrngorogmax = diag.Forcompute_range_faces(lon, lat,\
+          hgt, face, Nptfilt, Nptnewrange, dnamesvar, dvnamesvar)
 
         # Removing the nonChecking variable-dimensions from the initial list
@@ -717 +717 @@
         if face == 'WE': axis = 'lon'
         elif face == 'SN': axis = 'lat'
+
+        ncvar.insert_variable(ncobj, 'orogmax', rngorogmax, diagoutd, diagoutvd,     \
+          newnc)
+        newnc.renameVariable('orogmax', 'rangeorogmax')
+        ovar = newnc.variables['rangeorogmax']
+        ncvar.set_attribute(ovar, 'deriv', axis)
+        stdn = ovar.standard_name
+        ncvar.set_attribute(ovar, 'standard_name', 'range_' + stdn)
+        Ln = ovar.long_name
+        ncvar.set_attribute(ovar, 'long_name', 'range ' + stdn)
+
+        ncvar.insert_variable(ncobj, 'ptorogmax', ptrngorogmax, diagoutd, diagoutvd, \
+          newnc)
+        newnc.renameVariable('ptorogmax', 'rangeptorogmax')
+        ovar = newnc.variables['rangeptorogmax']
+        ncvar.set_attribute(ovar, 'deriv', axis)
+        stdn = ovar.standard_name
+        ncvar.set_attribute(ovar, 'standard_name', 'range_' + stdn)
+        Ln = ovar.long_name
+        ncvar.set_attribute(ovar, 'long_name', 'range ' + stdn)
 
         ncvar.insert_variable(ncobj, 'orogmax', orogmax, diagoutd, diagoutvd,        \

trunk/tools/module_ForDiagnostics.f90

@@ -950 +950 @@
 
   SUBROUTINE compute_range_faces(d1, d2, lon, lat, hgt, face, Nfilt, newrange, hgtmax, pthgtmax,      &
-    derivhgt, peaks, valleys, origfaces, filtfaces)
+    derivhgt, peaks, valleys, origfaces, filtfaces, rangeshgtmax, ptrangeshgtmax)
 ! Subroutine to compute faces [uphill, valleys, downhill] of a mountain range along a given face
 
@@ -958 +958 @@
     REAL(r_k), DIMENSION(d1,d2), INTENT(in)              :: lon, lat, hgt
     CHARACTER(len=*)                                     :: face
-    REAL(r_k), DIMENSION(d1,d2), INTENT(out)             :: derivhgt, hgtmax
+    REAL(r_k), DIMENSION(d1,d2), INTENT(out)             :: derivhgt, hgtmax, rangeshgtmax
     INTEGER, DIMENSION(d1,d2), INTENT(out)               :: pthgtmax, origfaces, filtfaces, peaks,    &
-      valleys
- 
+      valleys, ptrangeshgtmax
 ! Local
     INTEGER                                              :: i, j
-    INTEGER                                              :: pthgtmax1, Npeaks, Nvalleys
+    INTEGER                                              :: pthgtmax1, Npeaks, Nvalleys, Nranges
     REAL(r_k)                                            :: hgtmax1
-    INTEGER, DIMENSION(d1)                               :: ipeaks1, ivalleys1
-    INTEGER, DIMENSION(d2)                               :: ipeaks2, ivalleys2
+    INTEGER, DIMENSION(d1)                               :: ipeaks1, ivalleys1, irangeshgtmax1
+    INTEGER, DIMENSION(d2)                               :: ipeaks2, ivalleys2, irangeshgtmax2
+    REAL(r_k), DIMENSION(d1)                             :: rangeshgtmax1
+    REAL(r_k), DIMENSION(d2)                             :: rangeshgtmax2
+    INTEGER, DIMENSION(2,d1)                             :: ranges1
+    INTEGER, DIMENSION(2,d2)                             :: ranges2
 
 !!!!!!! Variables
@@ -983 +986 @@
 ! origfaces: original faces (-1, downhill; 0: valley; 1: uphill)
 ! filtfaces: filtered faces (-1, downhill; 0: valley; 1: uphill)
+! rangeshgtmax: maximum height for each individual range [m]
+! ptrangeshgtmax: grid point maximum height for each individual range [1]
 
     fname = 'compute_range_faces'
@@ -992 +997 @@
       DO j=1, d2
         CALL var_range_faces(d1, lon(:,j), lat(:,j), hgt(:,j), Nfilt, newrange, hgtmax1, pthgtmax1,   &
-          derivhgt(:,j), Npeaks, ipeaks1, Nvalleys, ivalleys1, origfaces(:,j), filtfaces(:,j))
+          derivhgt(:,j), Npeaks, ipeaks1, Nvalleys, ivalleys1, origfaces(:,j), filtfaces(:,j),        &
+          Nranges, ranges1, rangeshgtmax1, irangeshgtmax1)
         hgtmax(:,j) = hgtmax1
         pthgtmax(pthgtmax1,j) = 1
@@ -1001 +1007 @@
           valleys(ivalleys1(i),j) = 1
         END DO
+        DO i=1, Nranges
+          rangeshgtmax(ranges1(1,i):ranges1(2,i),j) = rangeshgtmax1(i)
+          ptrangeshgtmax(irangeshgtmax1(i),j) = 1
+        END DO
       END DO
     ELSE IF (TRIM(face) == 'SN') THEN
       DO i=1, d1
         CALL var_range_faces(d2, lon(i,:), lat(i,:), hgt(i,:), Nfilt, newrange, hgtmax1, pthgtmax1,   &
-          derivhgt(i,:), Npeaks, ipeaks2, Nvalleys, ivalleys2, origfaces(i,:), filtfaces(i,:))
+          derivhgt(i,:), Npeaks, ipeaks2, Nvalleys, ivalleys2, origfaces(i,:), filtfaces(i,:),        &
+          Nranges, ranges2, rangeshgtmax2, irangeshgtmax2)
         hgtmax(i,:) = hgtmax1
         pthgtmax(i,pthgtmax1) = 1
@@ -1013 +1024 @@
         DO j=1, Nvalleys
           valleys(i,ivalleys2(j)) = 1
+        END DO
+        DO j=1, Nranges
+          rangeshgtmax(i,ranges2(1,j):ranges2(2,j)) = rangeshgtmax2(j)
+          ptrangeshgtmax(i,irangeshgtmax2(j)) = 1
         END DO
       END DO

trunk/tools/module_ForDiagnosticsVars.f90

@@ -1581 +1581 @@
 
   SUBROUTINE var_range_faces(d, lon, lat, hgt, filt, newrange, hgtmax, ihgtmax, dhgt, Npeaks, ipeaks, &
-    Nvalleys, ivalleys, faces0, faces)
+    Nvalleys, ivalleys, faces0, faces, Nranges, ranges, rangeshgtmax, irangeshgtmax)
 ! Subroutine to compute faces [uphill, valleys, downhill] of a monuntain range along a face
 
@@ -1588 +1588 @@
     INTEGER, INTENT(in)                                  :: d, filt, newrange
     REAL(r_k), DIMENSION(d), INTENT(in)                  :: lon, lat, hgt
-    INTEGER, INTENT(out)                                 :: ihgtmax
+    INTEGER, INTENT(out)                                 :: ihgtmax, Npeaks, Nvalleys, Nranges
+    INTEGER, DIMENSION(d), INTENT(out)                   :: ipeaks, ivalleys, faces0, faces,          &
+      irangeshgtmax
+    INTEGER, DIMENSION(2,d), INTENT(out)                 :: ranges
     REAL(r_k), INTENT(out)                               :: hgtmax
     REAL(r_k), DIMENSION(d), INTENT(out)                 :: dhgt, rangeshgtmax
-    INTEGER, DIMENSION(d), INTENT(out)                   :: ipeaks, ivalleys, faces0, faces,          &
-      irangeshgtmax
-    INTEGER, INTENT(out)                                 :: Npeaks, Nvalleys, Nranges
 
 ! Local
     INTEGER                                              :: i, j, j1, k, l, m 
-    INTEGER                                              :: iface, Nranges
+    INTEGER                                              :: iface
     INTEGER                                              :: Nfaces, Nfaces1, Npeaks1, Nvalleys1
     INTEGER                                              :: fbeg, fend
@@ -1604 +1604 @@
     INTEGER, DIMENSION(d)                                :: faces1, Ndhgt1, ipeaks1, ivalleys1
     REAL(r_k), DIMENSION(d)                              :: dLl, peaks, valleys, peaks1, valleys1
-    LOGICAL                                              :: peakwithin, valleywithin
+    REAL(r_k), DIMENSION(d)                              :: sortedpeaks
+    INTEGER, DIMENSION(d)                                :: isortedpeaks
+    LOGICAL                                              :: rangewithin, peakwithin, valleywithin
 
 !!!!!!! Variables
@@ -1682 +1684 @@
     !PRINT *, Nfaces, ' length faces:', Ndhgt(1:Nfaces)
 
-    ! Defining quantitiy of ranges within the face
+    ! Defining quantitiy of ranges within the face 
+    !   sorting peaks within the face and defining ranges as maximum peaks distanced > newrage
+    sortedpeaks = peaks
+    CALL SortR_K(sortedpeaks, Npeaks)
+    isortedpeaks = 0
+    DO i=1, Npeaks
+      DO j=1, Npeaks
+        IF (peaks(j) == sortedpeaks(i)) isortedpeaks(i) = ipeaks(j)
+      END DO
+    END DO
+ 
+    ranges = -1
     Nranges = 1
-    IF (Npeaks > 2) THEN
-      DO i=1, Npeaks-1
-        IF (ipeaks(j+1)-ipeaks(j) < newrange)) THEN
-         AQUI
+    ranges(1,Nranges) = ipeaks(1)
+    ranges(2,Nranges) = ipeaks(1)
+    irangeshgtmax(Nranges) = ipeaks(1)
+    rangeshgtmax(Nranges) = peaks(1)
+    IF (Npeaks >= 2) THEN
+      DO j=2, Npeaks
+        rangewithin = .FALSE.
+        IF (ranges(1,Nranges) - isortedpeaks(j) < newrange) THEN
+          ranges(2,Nranges) = ipeaks(j)
+          IF (peaks(j) > rangeshgtmax(Nranges)) THEN
+            rangeshgtmax(Nranges) = peaks(j)
+            irangeshgtmax(Nranges) = ipeaks(j)
+          END IF
+        ELSE
+          Nranges = Nranges + 1
+          ranges(1,Nranges) = ipeaks(j)
+          ranges(2,Nranges) = ipeaks(j)
+          rangeshgtmax(Nranges) = peaks(j)
+          irangeshgtmax(Nranges) = ipeaks(j)
         END IF
       END DO
     END IF
-
-    ! Defining valleys as that consecutive grid points below surroungind peaks from the same range
+    !PRINT *, Nranges, ' _______'
+    !DO i=, Nranges
+    !   PRINT *,i,':', ranges(:,i),' max:', rangeshgtmax(i), ' ,', irangeshgtmax(i)
+    !END DO
+
+    ! Defining valleys as that consecutive grid points below surrounding peaks from the same range and
+    !   below range max
+    k = 1
     IF (Npeaks > 1) THEN
       j = 1
       j1 = 2
       DO i=2, d
-        IF (i >= ipeaks(j) .AND. i < ipeaks(j+1) .AND. (ipeaks(j+1)-ipeaks(j) < newrange)) THEN
+        IF (i >= ipeaks(j) .AND. i < ipeaks(j+1) .AND. (ipeaks(j+1) < ranges(2,k))) THEN
           IF (hgt(i) < peaks(j) .AND. hgt(i) < peaks(j+1)) faces0(i) = 0
-          IF (i == ipeaks(j) .AND. hgt(i) < hgtmax) faces0(i) = 0
+          IF (i == ipeaks(j) .AND. hgt(i) < rangeshgtmax(k)) faces0(i) = 0
         ELSE IF (i == ipeaks(j+1)) THEN
           j = j1
@@ -1706 +1740 @@
         END IF
       END DO
+      IF (ipeaks(j+1) > ranges(2,k)) k = k+1
     END IF