Philippe 07/03/2019: IO bugfix: io_set_mnhversion must be called by all the processes
[MNH-git_open_source-lfs.git] / src / SURFEX / prep_grid_gauss.F90
1 !SFX_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
2 !SFX_LIC This is part of the SURFEX software governed by the CeCILL-C licence
3 !SFX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
4 !SFX_LIC for details. version 1.
5 !     #########
6       SUBROUTINE PREP_GRID_GAUSS (&
7                                   HFILETYPE,HINTERP_TYPE,KNI)
8 !     ##########################################################################
9 !
10 !!****  *PREP_GRID_GAUSS* - reads EXTERNALIZED Surface grid.
11 !!
12 !!    PURPOSE
13 !!    -------
14 !!
15 !!**  METHOD
16 !!    ------
17 !!
18 !!    EXTERNAL
19 !!    --------
20 !!
21 !!    IMPLICIT ARGUMENTS
22 !!    ------------------
23 !!
24 !!
25 !!    REFERENCE
26 !!    ---------
27 !!
28 !!
29 !!    AUTHOR
30 !!    ------
31 !!
32 !!      V. Masson
33 !!
34 !!    MODIFICATIONS
35 !!    -------------
36 !!      Original   06/2003
37 !!      M. Jidane    Nov 2013 : correct allocation of NINLO and reading of INLOPA
38 !!      F. Taillefer Dec 2013 : debug estimation of XILO2
39 !-------------------------------------------------------------------------------
40 !
41 !*      0. DECLARATIONS
42 !          ------------
43 !
44 !
45 !
46 !
47 USE MODI_READ_SURF
48 !
49 USE MODD_GRID_GAUSS, ONLY : XILA1, XILO1, XILA2, XILO2, NINLA, NINLO, NILEN, LROTPOLE, &
50                             XCOEF, XLAP, XLOP, XILATARRAY
51 !
52 USE YOMHOOK   ,ONLY : LHOOK,   DR_HOOK
53 USE PARKIND1  ,ONLY : JPRB
54 !
55 IMPLICIT NONE
56 !
57 !* 0.1. Declaration of arguments
58 !       ------------------------
59 !
60 !
61 !
62  CHARACTER(LEN=6),  INTENT(IN)    :: HFILETYPE    ! file type
63  CHARACTER(LEN=6),  INTENT(OUT)   :: HINTERP_TYPE ! Grid type
64 INTEGER,           INTENT(OUT)   :: KNI          ! number of points
65 !
66 !* 0.2 Declaration of local variables
67 !      ------------------------------
68 !
69  CHARACTER(LEN=LEN_HREC) :: YRECFM    ! Name of the article to be read
70 INTEGER           :: IRESP
71 !
72 !
73 REAL, DIMENSION(:), ALLOCATABLE :: ZLAT    ! latitudes
74 REAL, DIMENSION(:), ALLOCATABLE :: ZW ! work array
75 !
76 INTEGER :: JL, ICPT        ! loop counter
77 INTEGER :: INLATI  ! number of pseudo-latitudes
78 INTEGER :: INLATI2 ! number of half pseudo-latitudes
79 REAL    :: ZLAPO   ! latitude of the rotated pole  (deg)
80 REAL    :: ZLOPO   ! longitude of the rotated pole (deg)
81 REAL    :: ZCODIL  ! stretching factor (must be greater than or equal to 1)
82 INTEGER, DIMENSION(:), ALLOCATABLE :: INLOPA ! number of pseudo-longitudes on each
83                                              ! pseudo-latitude circle
84 REAL(KIND=JPRB) :: ZHOOK_HANDLE
85 !
86 !-----------------------------------------------------------------------
87 IF (LHOOK) CALL DR_HOOK('PREP_GRID_GAUSS',0,ZHOOK_HANDLE)
88 !-----------------------------------------------------------------------
89 !
90 !*   1 Projection
91 !      ----------
92 !
93 YRECFM = 'LAPO'
94  CALL READ_SURF(HFILETYPE,YRECFM,ZLAPO,IRESP)
95 YRECFM = 'LOPO'
96  CALL READ_SURF(HFILETYPE,YRECFM,ZLOPO,IRESP)
97 YRECFM = 'CODIL'
98  CALL READ_SURF(HFILETYPE,YRECFM,ZCODIL,IRESP)
99 !
100 !-----------------------------------------------------------------------
101 !
102 !*   2 Grid
103 !      ----
104 !
105 YRECFM = 'NLATI'
106  CALL READ_SURF(HFILETYPE,YRECFM,INLATI,IRESP)
107 !
108 IF (ALLOCATED(INLOPA)) DEALLOCATE(INLOPA)
109 ALLOCATE(INLOPA(INLATI))
110 IF (ALLOCATED(NINLO)) DEALLOCATE(NINLO)
111 ALLOCATE(NINLO(INLATI))
112 YRECFM = 'NLOPA'
113  CALL READ_SURF(HFILETYPE,YRECFM,INLOPA,IRESP,HDIR='-')
114 !
115 KNI = SUM(INLOPA)
116 !
117 ALLOCATE(ZLAT(KNI))
118 ! CALL READ_SURF(HFILETYPE,'LATGAUSS',ZLAT(:),IRESP,HDIR='-')
119  CALL READ_SURF(HFILETYPE,'LAT_G_XY',ZLAT(:),IRESP,HDIR='-')
120 !
121 IF (ALLOCATED(XILATARRAY)) DEALLOCATE(XILATARRAY)
122 ALLOCATE(XILATARRAY(INLATI))
123 XILATARRAY(1) = ZLAT(1)
124 ICPT = 1
125 DO JL = 2,KNI
126   IF (ZLAT(JL)/=ZLAT(JL-1)) THEN
127     ICPT = ICPT + 1
128     XILATARRAY(ICPT) = ZLAT(JL)
129   ENDIF
130 ENDDO
131 !
132 DEALLOCATE(ZLAT)
133 !-----------------------------------------------------------------------
134 !
135 !*   3 Computes additional quantities used in interpolation
136 !      ----------------------------------------------------
137 !
138 INLATI2  = NINT(REAL(INLATI)/2.0)
139 NINLA    = INLATI
140 NILEN    = KNI
141 XLOP     = ZLOPO
142 XLAP     = ZLAPO
143 XCOEF    = ZCODIL
144 !
145 NINLO(:) = INLOPA(:)
146 !
147 !* type of transform
148 IF (ZLAPO>89.99 .AND. ABS(ZLOPO)<0.00001) THEN
149   LROTPOLE = .FALSE.
150 ELSE
151   LROTPOLE = .TRUE.
152 ENDIF
153 !
154 !XILA1=90.0*(1.0-1.0/(REAL(INLATI)))
155 !XILA2=-90.0*(1.0-1.0/(REAL(INLATI)))
156 XILA1 = XILATARRAY(1)
157 XILA2 = XILATARRAY(INLATI)
158 XILO1=0.0
159 XILO2=360.0*(REAL(INLOPA(INLATI2))-1.0)/REAL(INLOPA(INLATI2))
160 !
161 HINTERP_TYPE = 'HORIBL'
162 !-----------------------------------------------------------------------
163 IF (LHOOK) CALL DR_HOOK('PREP_GRID_GAUSS',1,ZHOOK_HANDLE)
164 !-----------------------------------------------------------------------
165 !
166 END SUBROUTINE PREP_GRID_GAUSS