Changeset 3605 for LMDZ6/branches/Ocean_skin/libf/phylmd/orografi_strato.F90

Timestamp:

Nov 21, 2019, 4:43:45 PM (4 years ago)

Author:

lguez

Message:

Merge revisions 3427:3600 of trunk into branch Ocean_skin

Location:

LMDZ6/branches/Ocean_skin

Files:

• . (modified) (1 prop)
• libf/phylmd/orografi_strato.F90 (modified) (3 diffs)

LMDZ6/branches/Ocean_skin/libf/phylmd/orografi_strato.F90

@@ -1819 +1819 @@
   USE mod_phys_lmdz_para
   USE mod_grid_phy_lmdz
+  USE geometry_mod
   IMPLICIT NONE
 
@@ -1832 +1833 @@
   INTEGER jk
   REAL zpr, ztop, zsigt, zpm1r
-  REAL :: pplay_glo(klon_glo, nlev)
-  REAL :: paprs_glo(klon_glo, nlev+1)
+  INTEGER :: cell,ij,nstra_tmp,nktopg_tmp
+  REAL :: current_dist, dist_min,dist_min_glo
 
   ! *       1.    SET THE VALUES OF THE PARAMETERS
@@ -1848 +1849 @@
   ! old  ZSIGT=0.85
 
-  CALL gather(pplay, pplay_glo)
-  CALL bcast(pplay_glo)
-  CALL gather(paprs, paprs_glo)
-  CALL bcast(paprs_glo)
-
-  DO jk = 1, nlev
-    zpm1r = pplay_glo(klon_glo/2+1, jk)/paprs_glo(klon_glo/2+1, 1)
-    IF (zpm1r>=zsigt) THEN
-      nktopg = jk
-    END IF
-    zpm1r = pplay_glo(klon_glo/2+1, jk)/paprs_glo(klon_glo/2+1, 1)
-    IF (zpm1r>=ztop) THEN
-      nstra = jk
-    END IF
-  END DO
-
+
+!ym Take the point at equator close to (0,0) coordinates.
+  dist_min=360
+  dist_min_glo=360.
+  cell=-1
+  DO ij=1,klon
+    current_dist=sqrt(longitude_deg(ij)**2+latitude_deg(ij)**2)
+    current_dist=current_dist*(1+(1e-10*ind_cell_glo(ij))/klon_glo) ! For point unicity
+    IF (dist_min>current_dist) THEN
+      dist_min=current_dist
+      cell=ij    
+    ENDIF  
+  ENDDO
+  
+  !PRINT *, 'SUGWD distmin cell=', dist_min,cell
+  CALL reduce_min(dist_min,dist_min_glo)
+  CALL bcast(dist_min_glo)
+  IF (dist_min/=dist_min_glo) cell=-1
+!ym in future find the point at equator close to (0,0) coordinates.
+  PRINT *, 'SUGWD distmin dist_min_glo cell=', dist_min,dist_min_glo,cell
+
+  nktopg_tmp=nktopg
+  nstra_tmp=nstra
+  
+  IF (cell/=-1) THEN
+
+    !print*,'SUGWD shape ',shape(pplay),cell+1
+
+    DO jk = 1, nlev
+      !zpm1r = pplay(cell+1, jk)/paprs(cell+1, 1)
+      zpm1r = pplay(cell, jk)/paprs(cell, 1)
+      IF (zpm1r>=zsigt) THEN
+        nktopg_tmp = jk
+      END IF
+      IF (zpm1r>=ztop) THEN
+        nstra_tmp = jk
+      END IF
+    END DO
+  ELSE
+    nktopg_tmp=0
+    nstra_tmp=0
+  ENDIF
+  
+  CALL reduce_sum(nktopg_tmp,nktopg)
+  CALL bcast(nktopg)
+  CALL reduce_sum(nstra_tmp,nstra)
+  CALL bcast(nstra)
+  
   ! inversion car dans orodrag on compte les niveaux a l'envers
   nktopg = nlev - nktopg + 1