Changeset 2208 in lmdz_wrf for trunk/tools

Timestamp:

Oct 31, 2018, 8:23:47 PM (6 years ago)

Author:

lfita

Message:

Adding:

• `range_faces': Function to compute faces [uphill, valley, downhill] of sections of a mountain rage, along a given face

Location:

trunk/tools

Files:

• diag_tools.py (modified) (1 diff)
• diagnostics.py (modified) (1 diff)
• module_ForDiagnostics.f90 (modified) (1 diff)
• module_ForDiagnosticsVars.f90 (modified) (1 diff)

trunk/tools/diag_tools.py

@@ -1367 +1367 @@
 
     return var1, var2, vardims, varvdims
+
+def Forcompute_range_faces(lon, lat, hgt, face, dimns, dimvns):
+    """ Function to compute faces [uphill, valley, downhill] of sections of a mountain
+      rage, along a given face
+    Forcompute_range_faces(lon, lat, hgt, face, dimns, dimvns)
+      [lon]= longitude values (assuming [y,x]) [degrees east]
+      [lat]= latitude values (assuming [y,x]) [degrees north]
+      [hgt]= height values (assuming [y,x]) [m]
+      [dimns]= list of the name of the dimensions of [smois]
+      [dimvns]= list of the name of the variables with the values of the 
+        dimensions of [smois]
+    """
+    fname = 'Forcompute_range_faces'
+
+    vardims = dimns[:]
+    varvdims = dimvns[:]
+
+    if len(var.shape) == 2:
+        faces = np.zeros((var.shape[0],var.shape[1]), dtype=np.float)
+        dx = var.shape[1]
+        dy = var.shape[0]
+
+        facest=fdin.module_fordiagnostics.compute_range_faces(lon=lon[:].transpose(),\
+          lat=lat[:].transpose(), hgt=hgt[:].transpose(), face=face, d1=dx, d2=dy)
+        faces = facest.transpose()
+    else:
+        print errormsg
+        print '  ' + fname + ': rank', len(var.shape), 'not ready !!'
+        print '  it only computes 2D [y,x] rank values'
+        quit(-1)
+
+    return faces, vardims, varvdims
 
 ####### ###### ##### #### ### ## # END Fortran diagnostics

trunk/tools/diagnostics.py

@@ -676 +676 @@
         ncvar.insert_variable(ncobj, 'hurs', diagout, diagoutd, diagoutvd, newnc)
 
+# range_faces: LON, LAT, HGT, 'face,[WE],[SN]'
+    elif diagn == 'range_faces':
+            
+        var0 = ncobj.variables[depvars[0]][:]
+        var1 = ncobj.variables[depvars[1]][:]
+        var2 = ncobj.variables[depvars[2]][:]
+        var3 = depvars[3].split(',')[1]
+
+        dnamesvar = list(ncobj.variables[depvars[0]].dimensions)
+        dvnamesvar = ncvar.var_dim_dimv(dnamesvar,dnames,dvnames)
+
+        diagout, diagoutd, diagoutvd = diag.Forcompute_range_faces(var0, var1, var2, \
+          var3, dnamesvar, dvnamesvar)
+
+        # Removing the nonChecking variable-dimensions from the initial list
+        varsadd = []
+        diagoutvd = list(dvnames)
+        for nonvd in NONchkvardims:
+            if gen.searchInlist(dvnames,nonvd): diagoutvd.remove(nonvd)
+            varsadd.append(nonvd)
+        ncvar.insert_variable(ncobj, 'range_faces', diagout, diagoutd, diagoutvd,    \
+          newnc)
+
 # mrso: total soil moisture SMOIS, DZS
     elif diagn == 'WRFmrso':

trunk/tools/module_ForDiagnostics.f90

@@ -948 +948 @@
   END SUBROUTINE compute_fog_FRAML50
 
+  SUBROUTINE compute_range_faces(d1, d2, lon, lat, hgt, face, faces)
+! Subroutine to compute faces [uphill, valleys, downhill] of a mountain range along a given face
+
+    IMPLICIT NONE
+
+    INTEGER, INTENT(in)                                  :: d1, d2
+    REAL(r_k), DIMENSION(d1,d2), INTENT(in)              :: lon, lat, hgt
+    CHARACTER(len=*)                                     :: face
+    INTEGER, DIMENSION(d1,d2), INTENT(out)               :: faces
+ 
+! Local
+    INTEGER                                              :: i, j
+
+!!!!!!! Variables
+! lon: longitude [degrees east]
+! lat: latitude [degrees north]
+! hgt: topograpical height [m]
+
+    fname = 'compute_range_faces'
+
+    IF (TRIM(face) == 'WE') THEN
+      DO j=1, d2
+        CALL var_range_faces(d1, lon(:,j), lat(:,j), hgt(:,j), faces(:,j))
+      END DO
+    ELSE IF (TRIM(face) == 'SN') THEN
+      DO i=1, d1
+        CALL var_range_faces(d2, lon(i,:), lat(i,:), hgt(i,:), faces(i,:))
+      END DO
+    ELSE
+      PRINT *,TRIM(ErrWarnMsg('err'))
+      PRINT *,'  ' // TRIM(fname) // ": wrong face: '" // TRIM(face) // "' !!"
+      PRINT *,'    accepted ones: WE, SN'
+      STOP
+    END IF
+
+    RETURN
+
+  END SUBROUTINE compute_range_faces
+
 END MODULE module_ForDiagnostics

trunk/tools/module_ForDiagnosticsVars.f90

@@ -1578 +1578 @@
   END SUBROUTINE var_fog_FRAML50
 
+  SUBROUTINE var_range_faces(d, lon, lat, hgt, faces)
+! Subroutine to compute faces of a monuntain range along a face
+
+    IMPLICIT NONE
+
+    INTEGER, INTENT(in)                                  :: d
+    REAL, DIMENSION(d), INTENT(in)                       :: lon, lat, hgt
+    INTEGER, DIMENSION(d), INTENT(out)                   :: faces
+
+! Local
+    INTEGER                                              :: i, Nfaces, Npeaks
+    INTEGER, DIMENSION(1)                                :: ihmax
+    REAL                                                 :: hgtmax
+    INTEGER, DIMENSION(d)                                :: ddhgt, Ndhgt, ipeaks
+    REAL, DIMENSION(d)                                   :: dhgt, peaks
+
+!!!!!!! Variables
+! lon: longitude [degrees east]
+! lat: latitude [degrees north]
+! hgt: topograpical height [m]
+
+    fname = 'var_range_faces'
+
+    PRINT *, 'heights:', hgt
+
+    ! Looking for the maximum height
+    hgtmax = MAXVAL(hgt)
+    ihmax = MAXLOC(hgt)
+
+    PRINT *, 'height max:', hgtmax, 'location:', ihmax
+
+    ! range slope
+    dhgt(1:d) = hgt(2:d) - hgt(1:d-1)
+
+    PRINT *, 'slope:', dhgt
+
+    ! Classification
+    Npeaks = 0
+    DO i=1, d-1
+      IF (dhgt(i) > 0.) THEN
+        faces(i) = 1
+      ELSE 
+        faces(i) = -1
+      END IF
+      ! peaks
+      IF (dhgt(i) > 0. .AND. dhgt(i+1) < 0.) THEN
+        Npeaks = Npeaks + 1
+        peaks(Npeaks) = hgt(i+1)
+        ipeaks(Npeaks) = i+1
+      END IF
+    END DO
+
+    PRINT *, 'faces:', faces
+    PRINT *, Npeaks, ' peaks:', peaks(1:Npeaks), ' ipeak:', ipeaks(1:Npeaks)
+
+    ! tendency of the slope
+    ddhgt(1) = 1 
+    Nfaces = 1
+    Ndhgt(Nfaces) = 1
+    DO i=2, d-1
+      IF (faces(i) /= faces(i-1)) THEN
+        ddhgt(i) = ddhgt(i-1) + 1
+        Nfaces = Nfaces + 1
+        Ndhgt(Nfaces) = 1
+      ELSE
+        Ndhgt(Nfaces) = Ndhgt(Nfaces) + 1
+      END IF
+    END DO
+
+    PRINT *, 'ddhgt:', ddhgt
+    PRINT *, Nfaces, ' length faces:', Ndhgt(1:Nfaces)
+
+    RETURN
+
+  END SUBROUTINE var_range_faces
+
   SUBROUTINE var_hur(t, press, qv, hur)
 ! Subroutine to compute relative humidity using August-Roche-Magnus approximation [1]