[2136] | 1 | !************************************************************************ |
---|
| 2 | ! Fortran 95 OPA Nesting tools * |
---|
| 3 | ! * |
---|
| 4 | ! Copyright (C) 2005 Florian Lemarié (Florian.Lemarie@imag.fr) * |
---|
| 5 | ! * |
---|
| 6 | !************************************************************************ |
---|
| 7 | ! |
---|
| 8 | MODULE agrif_readwrite |
---|
| 9 | ! |
---|
| 10 | USE agrif_types |
---|
| 11 | ! |
---|
| 12 | IMPLICIT NONE |
---|
| 13 | ! |
---|
| 14 | CONTAINS |
---|
| 15 | ! |
---|
| 16 | !************************************************************************ |
---|
| 17 | ! * |
---|
| 18 | ! MODULE AGRIF_READWRITE * |
---|
| 19 | ! * |
---|
| 20 | ! module containing subroutine used for : * |
---|
[2455] | 21 | ! - Coordinates files reading/writing * |
---|
[2136] | 22 | ! - Bathymetry files reading/writing (meter and levels) * |
---|
| 23 | ! - Naming of child grid files * |
---|
| 24 | ! * |
---|
| 25 | !************************************************************************ |
---|
| 26 | ! |
---|
| 27 | !***************************************************** |
---|
| 28 | ! function Read_Coordinates(name,Grid) |
---|
| 29 | !***************************************************** |
---|
| 30 | |
---|
| 31 | INTEGER FUNCTION Read_Coordinates(name,Grid,Pacifique) |
---|
| 32 | ! |
---|
| 33 | USE io_netcdf |
---|
| 34 | ! |
---|
| 35 | ! file name to open |
---|
| 36 | ! |
---|
| 37 | CHARACTER(*) name |
---|
| 38 | LOGICAL,OPTIONAL :: Pacifique |
---|
| 39 | ! |
---|
| 40 | TYPE(Coordinates) :: Grid |
---|
| 41 | ! |
---|
[10025] | 42 | CALL read_ncdf_var('glamt',name,Grid%glamt) |
---|
| 43 | CALL read_ncdf_var('glamu',name,Grid%glamu) |
---|
| 44 | CALL read_ncdf_var('glamv',name,Grid%glamv) |
---|
| 45 | CALL read_ncdf_var('glamf',name,Grid%glamf) |
---|
| 46 | CALL read_ncdf_var('gphit',name,Grid%gphit) |
---|
| 47 | CALL read_ncdf_var('gphiu',name,Grid%gphiu) |
---|
| 48 | CALL read_ncdf_var('gphiv',name,Grid%gphiv) |
---|
| 49 | CALL read_ncdf_var('gphif',name,Grid%gphif) |
---|
| 50 | CALL read_ncdf_var('e1t',name,Grid%e1t) |
---|
| 51 | CALL read_ncdf_var('e1u',name,Grid%e1u) |
---|
| 52 | CALL read_ncdf_var('e1v',name,Grid%e1v) |
---|
| 53 | CALL read_ncdf_var('e1f',name,Grid%e1f) |
---|
| 54 | CALL read_ncdf_var('e2t',name,Grid%e2t) |
---|
| 55 | CALL read_ncdf_var('e2u',name,Grid%e2u) |
---|
| 56 | CALL read_ncdf_var('e2v',name,Grid%e2v) |
---|
| 57 | CALL read_ncdf_var('e2f',name,Grid%e2f) |
---|
| 58 | CALL read_ncdf_var('nav_lon',name,Grid%nav_lon) |
---|
| 59 | CALL read_ncdf_var('nav_lat',name,Grid%nav_lat) |
---|
[2136] | 60 | ! |
---|
| 61 | IF( PRESENT(Pacifique) )THEN |
---|
| 62 | IF ( Grid%glamt(1,1) > Grid%glamt(nxfin,nyfin) ) THEN |
---|
| 63 | Pacifique = .TRUE. |
---|
| 64 | WHERE ( Grid%glamt < 0 ) |
---|
| 65 | Grid%glamt = Grid%glamt + 360. |
---|
| 66 | END WHERE |
---|
| 67 | WHERE ( Grid%glamf < 0 ) |
---|
| 68 | Grid%glamf = Grid%glamf + 360. |
---|
| 69 | END WHERE |
---|
| 70 | WHERE ( Grid%glamu < 0 ) |
---|
| 71 | Grid%glamu = Grid%glamu + 360. |
---|
| 72 | END WHERE |
---|
| 73 | WHERE ( Grid%glamv < 0 ) |
---|
| 74 | Grid%glamv = Grid%glamv + 360. |
---|
| 75 | END WHERE |
---|
| 76 | WHERE ( Grid%nav_lon < 0 ) |
---|
| 77 | Grid%nav_lon = Grid%nav_lon + 360. |
---|
| 78 | END WHERE |
---|
| 79 | ENDIF |
---|
| 80 | ENDIF |
---|
| 81 | ! |
---|
| 82 | WRITE(*,*) ' ' |
---|
| 83 | WRITE(*,*) 'Reading coordinates file: ',name |
---|
| 84 | WRITE(*,*) ' ' |
---|
| 85 | ! |
---|
| 86 | Read_Coordinates = 1 |
---|
| 87 | ! |
---|
| 88 | END FUNCTION Read_Coordinates |
---|
| 89 | |
---|
| 90 | !***************************************************** |
---|
| 91 | ! function Read_Coordinates(name,Grid) |
---|
| 92 | !***************************************************** |
---|
| 93 | |
---|
| 94 | INTEGER FUNCTION Read_Local_Coordinates(name,Grid,strt,cnt) |
---|
| 95 | ! |
---|
| 96 | USE io_netcdf |
---|
| 97 | ! |
---|
| 98 | ! file name to open |
---|
| 99 | ! |
---|
| 100 | CHARACTER(*) name |
---|
| 101 | INTEGER, DIMENSION(2) :: strt,cnt |
---|
| 102 | ! |
---|
| 103 | TYPE(Coordinates) :: Grid |
---|
| 104 | ! |
---|
[10025] | 105 | CALL read_ncdf_var('glamt',name,Grid%glamt,strt,cnt) |
---|
| 106 | CALL read_ncdf_var('glamu',name,Grid%glamu,strt,cnt) |
---|
| 107 | CALL read_ncdf_var('glamv',name,Grid%glamv,strt,cnt) |
---|
| 108 | CALL read_ncdf_var('glamf',name,Grid%glamf,strt,cnt) |
---|
| 109 | CALL read_ncdf_var('gphit',name,Grid%gphit,strt,cnt) |
---|
| 110 | CALL read_ncdf_var('gphiu',name,Grid%gphiu,strt,cnt) |
---|
| 111 | CALL read_ncdf_var('gphiv',name,Grid%gphiv,strt,cnt) |
---|
| 112 | CALL read_ncdf_var('gphif',name,Grid%gphif,strt,cnt) |
---|
| 113 | CALL read_ncdf_var('e1t',name,Grid%e1t,strt,cnt) |
---|
| 114 | CALL read_ncdf_var('e1u',name,Grid%e1u,strt,cnt) |
---|
| 115 | CALL read_ncdf_var('e1v',name,Grid%e1v,strt,cnt) |
---|
| 116 | CALL read_ncdf_var('e1f',name,Grid%e1f,strt,cnt) |
---|
| 117 | CALL read_ncdf_var('e2t',name,Grid%e2t,strt,cnt) |
---|
| 118 | CALL read_ncdf_var('e2u',name,Grid%e2u,strt,cnt) |
---|
| 119 | CALL read_ncdf_var('e2v',name,Grid%e2v,strt,cnt) |
---|
| 120 | CALL read_ncdf_var('e2f',name,Grid%e2f,strt,cnt) |
---|
| 121 | CALL read_ncdf_var('nav_lon',name,Grid%nav_lon,strt,cnt) |
---|
| 122 | CALL read_ncdf_var('nav_lat',name,Grid%nav_lat,strt,cnt) |
---|
[2136] | 123 | ! |
---|
| 124 | WRITE(*,*) ' ' |
---|
| 125 | WRITE(*,*) 'Reading coordinates file: ',name |
---|
| 126 | WRITE(*,*) ' ' |
---|
| 127 | ! |
---|
| 128 | Read_Local_Coordinates = 1 |
---|
| 129 | ! |
---|
| 130 | END FUNCTION Read_Local_Coordinates |
---|
| 131 | |
---|
| 132 | !***************************************************** |
---|
| 133 | ! function Write_Coordinates(name,Grid) |
---|
| 134 | !***************************************************** |
---|
| 135 | |
---|
| 136 | INTEGER FUNCTION Write_Coordinates(name,Grid) |
---|
| 137 | ! |
---|
| 138 | USE io_netcdf |
---|
| 139 | CHARACTER(*) name |
---|
| 140 | TYPE(Coordinates) :: Grid |
---|
[2455] | 141 | INTEGER :: status,ncid |
---|
| 142 | CHARACTER(len=1),DIMENSION(2) :: dimnames |
---|
[2136] | 143 | ! |
---|
| 144 | status = nf90_create(name,NF90_WRITE,ncid) |
---|
| 145 | status = nf90_close(ncid) |
---|
| 146 | ! |
---|
[2455] | 147 | dimnames = (/ 'x','y' /) |
---|
[10025] | 148 | CALL write_ncdf_dim(dimnames(1),name,nxfin) |
---|
| 149 | CALL write_ncdf_dim(dimnames(2),name,nyfin) |
---|
[2136] | 150 | ! |
---|
[10025] | 151 | CALL write_ncdf_var('nav_lon',dimnames,name,Grid%nav_lon,'float') |
---|
| 152 | CALL write_ncdf_var('nav_lat',dimnames,name,Grid%nav_lat,'float') |
---|
[2136] | 153 | ! |
---|
[10025] | 154 | CALL write_ncdf_var('glamt',dimnames,name,Grid%glamt,'double') |
---|
| 155 | CALL write_ncdf_var('glamu',dimnames,name,Grid%glamu,'double') |
---|
| 156 | CALL write_ncdf_var('glamv',dimnames,name,Grid%glamv,'double') |
---|
| 157 | CALL write_ncdf_var('glamf',dimnames,name,Grid%glamf,'double') |
---|
| 158 | CALL write_ncdf_var('gphit',dimnames,name,Grid%gphit,'double') |
---|
| 159 | CALL write_ncdf_var('gphiu',dimnames,name,Grid%gphiu,'double') |
---|
| 160 | CALL write_ncdf_var('gphiv',dimnames,name,Grid%gphiv,'double') |
---|
| 161 | CALL write_ncdf_var('gphif',dimnames,name,Grid%gphif,'double') |
---|
| 162 | CALL write_ncdf_var('e1t',dimnames,name,Grid%e1t,'double') |
---|
| 163 | CALL write_ncdf_var('e1u',dimnames,name,Grid%e1u,'double') |
---|
| 164 | CALL write_ncdf_var('e1v',dimnames,name,Grid%e1v,'double') |
---|
| 165 | CALL write_ncdf_var('e1f',dimnames,name,Grid%e1f,'double') |
---|
| 166 | CALL write_ncdf_var('e2t',dimnames,name,Grid%e2t,'double') |
---|
| 167 | CALL write_ncdf_var('e2u',dimnames,name,Grid%e2u,'double') |
---|
| 168 | CALL write_ncdf_var('e2v',dimnames,name,Grid%e2v,'double') |
---|
| 169 | CALL write_ncdf_var('e2f',dimnames,name,Grid%e2f,'double') |
---|
[2136] | 170 | ! |
---|
[10025] | 171 | CALL copy_ncdf_att('nav_lon',TRIM(parent_coordinate_file),name,MINVAL(Grid%nav_lon),MAXVAL(Grid%nav_lon)) |
---|
| 172 | CALL copy_ncdf_att('nav_lat',TRIM(parent_coordinate_file),name,MINVAL(Grid%nav_lat),MAXVAL(Grid%nav_lat)) |
---|
| 173 | CALL copy_ncdf_att('glamt',TRIM(parent_coordinate_file),name) |
---|
| 174 | CALL copy_ncdf_att('glamu',TRIM(parent_coordinate_file),name) |
---|
| 175 | CALL copy_ncdf_att('glamv',TRIM(parent_coordinate_file),name) |
---|
| 176 | CALL copy_ncdf_att('glamf',TRIM(parent_coordinate_file),name) |
---|
| 177 | CALL copy_ncdf_att('gphit',TRIM(parent_coordinate_file),name) |
---|
| 178 | CALL copy_ncdf_att('gphiu',TRIM(parent_coordinate_file),name) |
---|
| 179 | CALL copy_ncdf_att('gphiv',TRIM(parent_coordinate_file),name) |
---|
| 180 | CALL copy_ncdf_att('gphif',TRIM(parent_coordinate_file),name) |
---|
| 181 | CALL copy_ncdf_att('e1t',TRIM(parent_coordinate_file),name) |
---|
| 182 | CALL copy_ncdf_att('e1u',TRIM(parent_coordinate_file),name) |
---|
| 183 | CALL copy_ncdf_att('e1v',TRIM(parent_coordinate_file),name) |
---|
| 184 | CALL copy_ncdf_att('e1f',TRIM(parent_coordinate_file),name) |
---|
| 185 | CALL copy_ncdf_att('e2t',TRIM(parent_coordinate_file),name) |
---|
| 186 | CALL copy_ncdf_att('e2u',TRIM(parent_coordinate_file),name) |
---|
| 187 | CALL copy_ncdf_att('e2v',TRIM(parent_coordinate_file),name) |
---|
| 188 | CALL copy_ncdf_att('e2f',TRIM(parent_coordinate_file),name) |
---|
[2136] | 189 | ! |
---|
| 190 | WRITE(*,*) ' ' |
---|
| 191 | WRITE(*,*) 'Writing coordinates file: ',name |
---|
| 192 | WRITE(*,*) ' ' |
---|
| 193 | ! |
---|
| 194 | Write_Coordinates = 1 |
---|
| 195 | ! |
---|
| 196 | END FUNCTION Write_Coordinates |
---|
| 197 | ! |
---|
| 198 | ! |
---|
| 199 | ! |
---|
| 200 | !***************************************************** |
---|
| 201 | ! function Read_Bathy_level(name,Grid) |
---|
| 202 | !***************************************************** |
---|
| 203 | ! |
---|
| 204 | INTEGER FUNCTION Read_Bathy_level(name,Grid) |
---|
| 205 | ! |
---|
| 206 | USE io_netcdf |
---|
| 207 | ! |
---|
| 208 | CHARACTER(*) name |
---|
| 209 | TYPE(Coordinates) :: Grid |
---|
| 210 | ! |
---|
[10381] | 211 | CALL read_ncdf_var(parent_level_name,name,Grid%Bathy_level) |
---|
[2136] | 212 | ! |
---|
| 213 | WRITE(*,*) ' ' |
---|
| 214 | WRITE(*,*) 'Reading bathymetry file: ',name |
---|
| 215 | WRITE(*,*) ' ' |
---|
| 216 | ! |
---|
| 217 | Read_Bathy_level = 1 |
---|
| 218 | ! |
---|
| 219 | END FUNCTION Read_Bathy_level |
---|
| 220 | ! |
---|
| 221 | !***************************************************** |
---|
| 222 | ! function Write_Bathy_level(name,Grid) |
---|
| 223 | !***************************************************** |
---|
| 224 | ! |
---|
| 225 | INTEGER FUNCTION Write_Bathy_level(name,Grid) |
---|
| 226 | ! |
---|
| 227 | USE io_netcdf |
---|
| 228 | ! |
---|
| 229 | CHARACTER(*) name |
---|
| 230 | TYPE(Coordinates) :: Grid |
---|
[2455] | 231 | INTEGER :: status,ncid |
---|
| 232 | CHARACTER(len=1),DIMENSION(2) :: dimnames |
---|
[2136] | 233 | ! |
---|
[10025] | 234 | status = nf90_create(name,NF90_WRITE,ncid) |
---|
[2136] | 235 | status = nf90_close(ncid) |
---|
| 236 | ! |
---|
[2455] | 237 | dimnames = (/ 'x','y' /) |
---|
[10025] | 238 | CALL write_ncdf_dim(dimnames(1),name,nxfin) |
---|
| 239 | CALL write_ncdf_dim(dimnames(2),name,nyfin) |
---|
[2136] | 240 | ! |
---|
[10381] | 241 | CALL write_ncdf_var('nav_lon' ,dimnames,name,Grid%nav_lon ,'float') |
---|
| 242 | CALL write_ncdf_var('nav_lat' ,dimnames,name,Grid%nav_lat ,'float') |
---|
[12253] | 243 | CALL write_ncdf_var(parent_level_name,dimnames,name,NINT(Grid%bathy_level),'integer') |
---|
[2136] | 244 | ! |
---|
[10381] | 245 | CALL copy_ncdf_att('nav_lon' ,TRIM(parent_bathy_level),name,MINVAL(Grid%nav_lon),MAXVAL(Grid%nav_lon)) |
---|
| 246 | CALL copy_ncdf_att('nav_lat' ,TRIM(parent_bathy_level),name,MINVAL(Grid%nav_lat),MAXVAL(Grid%nav_lat)) |
---|
| 247 | CALL copy_ncdf_att(parent_level_name,TRIM(parent_bathy_level),name) |
---|
[2136] | 248 | ! |
---|
| 249 | WRITE(*,*) ' ' |
---|
[10093] | 250 | WRITE(*,*) 'Writing bathymetry (levels) in: ',name |
---|
[2136] | 251 | WRITE(*,*) ' ' |
---|
| 252 | ! |
---|
| 253 | Write_Bathy_level = 1 |
---|
| 254 | ! |
---|
| 255 | END FUNCTION Write_Bathy_level |
---|
| 256 | ! |
---|
| 257 | !***************************************************** |
---|
[10381] | 258 | ! function read_bathy_coord(name,CoarseGrid,ChildGrid) |
---|
[2136] | 259 | !***************************************************** |
---|
| 260 | ! |
---|
[10381] | 261 | INTEGER FUNCTION read_bathy_coord(name,CoarseGrid,ChildGrid,Pacifique) |
---|
[2136] | 262 | ! |
---|
| 263 | USE io_netcdf |
---|
| 264 | CHARACTER(*) name |
---|
| 265 | INTEGER :: i,j,tabdim1,tabdim2 |
---|
| 266 | INTEGER, DIMENSION(1) :: i_min,i_max,j_min,j_max |
---|
| 267 | REAL*8,POINTER,DIMENSION(:) :: topo_lon,topo_lat |
---|
| 268 | INTEGER :: status,ncid,varid |
---|
| 269 | LOGICAL,OPTIONAL :: Pacifique |
---|
| 270 | TYPE(Coordinates) :: CoarseGrid,ChildGrid |
---|
[9149] | 271 | REAL*8 :: zdel |
---|
| 272 | zdel = 0.5 ! Offset in degrees to extend extraction of bathymetry data |
---|
[2136] | 273 | ! |
---|
[2455] | 274 | IF( Dims_Existence('lon',name) .AND. Dims_Existence('lat',name) ) THEN |
---|
[2136] | 275 | WRITE(*,*) '****' |
---|
| 276 | WRITE(*,*) ' etopo format for external high resolution database ' |
---|
| 277 | WRITE(*,*) '****' |
---|
[10025] | 278 | CALL read_ncdf_var('lon',name,topo_lon) |
---|
| 279 | CALL read_ncdf_var('lat',name,topo_lat) |
---|
[2455] | 280 | ELSE IF( Dims_Existence('x',name) .AND. Dims_Existence('y',name) ) THEN |
---|
[2136] | 281 | WRITE(*,*) '****' |
---|
| 282 | WRITE(*,*) ' OPA format for external high resolution database ' |
---|
| 283 | WRITE(*,*) '****' |
---|
[10025] | 284 | CALL read_ncdf_var('nav_lon',name,CoarseGrid%nav_lon) |
---|
| 285 | CALL read_ncdf_var('nav_lat',name,CoarseGrid%nav_lat) |
---|
[10381] | 286 | CALL read_ncdf_var(parent_meter_name,name,CoarseGrid%Bathy_meter) |
---|
[2136] | 287 | ! |
---|
| 288 | IF ( PRESENT(Pacifique) ) THEN |
---|
| 289 | IF(Pacifique) THEN |
---|
| 290 | WHERE(CoarseGrid%nav_lon < 0.001) |
---|
| 291 | CoarseGrid%nav_lon = CoarseGrid%nav_lon + 360. |
---|
| 292 | END WHERE |
---|
| 293 | ENDIF |
---|
| 294 | ENDIF |
---|
| 295 | ! |
---|
[10381] | 296 | read_bathy_coord = 1 |
---|
[2136] | 297 | RETURN |
---|
| 298 | ELSE |
---|
| 299 | WRITE(*,*) '****' |
---|
| 300 | WRITE(*,*) '*** ERROR Bad format for external high resolution database' |
---|
| 301 | WRITE(*,*) '****' |
---|
| 302 | STOP |
---|
| 303 | ENDIF |
---|
| 304 | ! |
---|
[10025] | 305 | IF( MAXVAL(ChildGrid%glamt) > 180. ) THEN |
---|
[2136] | 306 | ! |
---|
[10025] | 307 | WHERE( topo_lon < 0. ) topo_lon = topo_lon + 360. |
---|
[2136] | 308 | ! |
---|
[9149] | 309 | i_min = MAXLOC(topo_lon,mask = topo_lon < MINVAL(ChildGrid%nav_lon)-zdel) |
---|
| 310 | i_max = MINLOC(topo_lon,mask = topo_lon > MAXVAL(ChildGrid%nav_lon)+zdel) |
---|
| 311 | j_min = MAXLOC(topo_lat,mask = topo_lat < MINVAL(ChildGrid%nav_lat)-zdel) |
---|
| 312 | j_max = MINLOC(topo_lat,mask = topo_lat > MAXVAL(ChildGrid%nav_lat)+zdel) |
---|
[2136] | 313 | ! |
---|
| 314 | tabdim1 = ( SIZE(topo_lon) - i_min(1) + 1 ) + i_max(1) |
---|
| 315 | ! |
---|
[10025] | 316 | IF( ln_agrif_domain ) THEN |
---|
| 317 | IF(j_min(1)-2 >= 1 .AND. j_max(1)+3 <= SIZE(topo_lat,1) ) THEN |
---|
| 318 | j_min(1) = j_min(1)-2 |
---|
| 319 | j_max(1) = j_max(1)+3 |
---|
| 320 | ENDIF |
---|
[2136] | 321 | ENDIF |
---|
[2455] | 322 | tabdim2 = j_max(1) - j_min(1) + 1 |
---|
[2136] | 323 | ! |
---|
| 324 | ALLOCATE(CoarseGrid%nav_lon(tabdim1,tabdim2)) |
---|
| 325 | ALLOCATE(CoarseGrid%nav_lat(tabdim1,tabdim2)) |
---|
| 326 | ALLOCATE(CoarseGrid%Bathy_meter(tabdim1,tabdim2)) |
---|
| 327 | ! |
---|
[2455] | 328 | DO i = 1,tabdim1 |
---|
| 329 | CoarseGrid%nav_lat(i,:) = topo_lat(j_min(1):j_max(1)) |
---|
[2136] | 330 | END DO |
---|
| 331 | ! |
---|
[2455] | 332 | DO j = 1, tabdim2 |
---|
[2136] | 333 | ! |
---|
[2455] | 334 | CoarseGrid%nav_lon(1:SIZE(topo_lon)-i_min(1)+1 ,j) = topo_lon(i_min(1):SIZE(topo_lon)) |
---|
| 335 | CoarseGrid%nav_lon(2+SIZE(topo_lon)-i_min(1):tabdim1,j) = topo_lon(1:i_max(1)) |
---|
[2136] | 336 | ! |
---|
| 337 | END DO |
---|
| 338 | status = nf90_open(name,NF90_NOWRITE,ncid) |
---|
[2455] | 339 | status = nf90_inq_varid(ncid,elevation_name,varid) |
---|
[9149] | 340 | ! |
---|
| 341 | IF (status/=nf90_noerr) THEN |
---|
| 342 | WRITE(*,*)"Can't find variable: ", elevation_name |
---|
| 343 | STOP |
---|
| 344 | ENDIF |
---|
[2136] | 345 | ! |
---|
[2455] | 346 | status=nf90_get_var(ncid,varid,CoarseGrid%Bathy_meter(1:SIZE(topo_lon)-i_min(1)+1,:), & |
---|
| 347 | start=(/i_min(1),j_min(1)/),count=(/SIZE(topo_lon)-i_min(1),tabdim2/)) |
---|
[2136] | 348 | |
---|
[2455] | 349 | status=nf90_get_var(ncid,varid,CoarseGrid%Bathy_meter(2+SIZE(topo_lon)-i_min(1):tabdim1,:), & |
---|
| 350 | start=(/1,j_min(1)/),count=(/i_max(1),tabdim2/)) |
---|
[2136] | 351 | ! |
---|
| 352 | ELSE |
---|
| 353 | ! |
---|
[10025] | 354 | WHERE( topo_lon > 180. ) topo_lon = topo_lon - 360. |
---|
| 355 | ! |
---|
[9149] | 356 | i_min = MAXLOC(topo_lon,mask = topo_lon < MINVAL(ChildGrid%nav_lon)-zdel) |
---|
| 357 | i_max = MINLOC(topo_lon,mask = topo_lon > MAXVAL(ChildGrid%nav_lon)+zdel) |
---|
| 358 | j_min = MAXLOC(topo_lat,mask = topo_lat < MINVAL(ChildGrid%nav_lat)-zdel) |
---|
| 359 | j_max = MINLOC(topo_lat,mask = topo_lat > MAXVAL(ChildGrid%nav_lat)+zdel) |
---|
[2136] | 360 | ! |
---|
[10025] | 361 | IF( ln_agrif_domain ) THEN |
---|
| 362 | IF(i_min(1)-2 >= 1 .AND. i_max(1)+3 <= SIZE(topo_lon,1) ) THEN |
---|
| 363 | i_min(1) = i_min(1)-2 |
---|
| 364 | i_max(1) = i_max(1)+3 |
---|
| 365 | ENDIF |
---|
[2136] | 366 | ENDIF |
---|
[2455] | 367 | tabdim1 = i_max(1) - i_min(1) + 1 |
---|
[2136] | 368 | ! |
---|
[10025] | 369 | IF( ln_agrif_domain ) THEN |
---|
| 370 | IF(j_min(1)-2 >= 1 .AND. j_max(1)+3 <= SIZE(topo_lat,1) ) THEN |
---|
| 371 | j_min(1) = j_min(1)-2 |
---|
| 372 | j_max(1) = j_max(1)+3 |
---|
| 373 | ENDIF |
---|
[2136] | 374 | ENDIF |
---|
[2455] | 375 | tabdim2 = j_max(1) - j_min(1) + 1 |
---|
[2136] | 376 | ! |
---|
| 377 | WRITE(*,*) ' ' |
---|
| 378 | WRITE(*,*) 'Reading bathymetry file: ',name |
---|
| 379 | WRITE(*,*) ' ' |
---|
| 380 | ! |
---|
| 381 | ALLOCATE(CoarseGrid%nav_lon(tabdim1,tabdim2)) |
---|
| 382 | ALLOCATE(CoarseGrid%nav_lat(tabdim1,tabdim2)) |
---|
| 383 | ALLOCATE(CoarseGrid%Bathy_meter(tabdim1,tabdim2)) |
---|
| 384 | ! |
---|
[2455] | 385 | DO j = 1,tabdim2 |
---|
| 386 | CoarseGrid%nav_lon(:,j) = topo_lon(i_min(1):i_max(1)) |
---|
[2136] | 387 | END DO |
---|
| 388 | ! |
---|
[2455] | 389 | DO i = 1,tabdim1 |
---|
| 390 | CoarseGrid%nav_lat(i,:) = topo_lat(j_min(1):j_max(1)) |
---|
[2136] | 391 | END DO |
---|
| 392 | ! |
---|
| 393 | status = nf90_open(name,NF90_NOWRITE,ncid) |
---|
[2455] | 394 | status = nf90_inq_varid(ncid,elevation_name,varid) |
---|
[9149] | 395 | IF (status/=nf90_noerr) THEN |
---|
| 396 | WRITE(*,*)"Can't find variable: ", elevation_name |
---|
| 397 | STOP |
---|
| 398 | ENDIF |
---|
| 399 | |
---|
[2455] | 400 | status = nf90_get_var(ncid,varid,CoarseGrid%Bathy_meter, & |
---|
| 401 | & start=(/i_min(1),j_min(1)/),count=(/tabdim1,tabdim2/)) |
---|
[2136] | 402 | ! |
---|
| 403 | ENDIF |
---|
| 404 | ! |
---|
| 405 | status = nf90_close(ncid) |
---|
| 406 | ! |
---|
| 407 | WHERE(CoarseGrid%Bathy_meter.GE.0) |
---|
| 408 | CoarseGrid%Bathy_meter = 0.0 |
---|
| 409 | END WHERE |
---|
| 410 | ! |
---|
[2455] | 411 | CoarseGrid%Bathy_meter(:,:) = -1.0 * CoarseGrid%Bathy_meter(:,:) |
---|
[2136] | 412 | ! |
---|
[10381] | 413 | read_bathy_coord = 1 |
---|
[2136] | 414 | RETURN |
---|
| 415 | ! |
---|
[10381] | 416 | END FUNCTION read_bathy_coord |
---|
[2136] | 417 | ! |
---|
| 418 | ! |
---|
| 419 | !***************************************************** |
---|
[10381] | 420 | ! function read_bathy_meter(name,CoarseGrid,ChildGrid) |
---|
[2136] | 421 | !***************************************************** |
---|
| 422 | ! |
---|
[10381] | 423 | INTEGER FUNCTION read_bathy_meter(name,Grid) |
---|
[2136] | 424 | ! |
---|
| 425 | ! |
---|
| 426 | USE io_netcdf |
---|
| 427 | ! |
---|
| 428 | CHARACTER(*) name |
---|
| 429 | TYPE(Coordinates) :: Grid |
---|
| 430 | ! |
---|
[10381] | 431 | CALL read_ncdf_var(parent_meter_name,name,Grid%Bathy_meter) |
---|
[2136] | 432 | ! |
---|
| 433 | WRITE(*,*) ' ' |
---|
| 434 | WRITE(*,*) 'Reading bathymetry file: ',name |
---|
| 435 | WRITE(*,*) ' ' |
---|
| 436 | ! |
---|
[10381] | 437 | read_bathy_meter = 1 |
---|
[2136] | 438 | ! |
---|
[10381] | 439 | END FUNCTION read_bathy_meter |
---|
[2136] | 440 | ! |
---|
| 441 | !***************************************************** |
---|
| 442 | ! function Write_Bathy_meter(name,Grid) |
---|
| 443 | !***************************************************** |
---|
| 444 | ! |
---|
| 445 | INTEGER FUNCTION Write_Bathy_meter(name,Grid) |
---|
| 446 | ! |
---|
| 447 | USE io_netcdf |
---|
| 448 | ! |
---|
| 449 | CHARACTER(*) name |
---|
| 450 | TYPE(Coordinates) :: Grid |
---|
[2455] | 451 | INTEGER :: status,ncid |
---|
| 452 | CHARACTER(len=1),DIMENSION(2) :: dimnames |
---|
[2136] | 453 | INTEGER :: nx,ny |
---|
| 454 | ! |
---|
| 455 | status = nf90_create(name,NF90_WRITE,ncid) |
---|
| 456 | status = nf90_close(ncid) |
---|
| 457 | ! |
---|
| 458 | nx = SIZE(Grid%bathy_meter,1) |
---|
| 459 | ny = SIZE(Grid%bathy_meter,2) |
---|
[2455] | 460 | dimnames = (/ 'x','y' /) |
---|
[2136] | 461 | |
---|
[10025] | 462 | CALL write_ncdf_dim(dimnames(1),name,nx) |
---|
| 463 | CALL write_ncdf_dim(dimnames(2),name,ny) |
---|
[2136] | 464 | ! |
---|
[10381] | 465 | CALL write_ncdf_var('nav_lon' ,dimnames,name,Grid%nav_lon ,'float') |
---|
| 466 | CALL write_ncdf_var('nav_lat' ,dimnames,name,Grid%nav_lat ,'float') |
---|
| 467 | CALL write_ncdf_var(parent_meter_name,dimnames,name,Grid%bathy_meter,'float') |
---|
| 468 | CALL write_ncdf_var('weight' ,dimnames,name,Grid%wgt ,'float') |
---|
[2136] | 469 | ! |
---|
[10381] | 470 | CALL copy_ncdf_att('nav_lon' ,TRIM(parent_bathy_meter),name,MINVAL(Grid%nav_lon),MAXVAL(Grid%nav_lon)) |
---|
| 471 | CALL copy_ncdf_att('nav_lat' ,TRIM(parent_bathy_meter),name,MINVAL(Grid%nav_lat),MAXVAL(Grid%nav_lat)) |
---|
| 472 | CALL copy_ncdf_att(parent_meter_name,TRIM(parent_bathy_meter),name) |
---|
[2136] | 473 | ! |
---|
| 474 | WRITE(*,*) ' ' |
---|
[10093] | 475 | WRITE(*,*) 'Writing bathymetry (meters) in: ',name |
---|
[2136] | 476 | WRITE(*,*) ' ' |
---|
| 477 | ! |
---|
| 478 | Write_Bathy_meter = 1 |
---|
| 479 | ! |
---|
| 480 | END FUNCTION Write_Bathy_meter |
---|
| 481 | ! |
---|
| 482 | !***************************************************** |
---|
[10025] | 483 | ! function write_domcfg(name,Grid) |
---|
| 484 | !***************************************************** |
---|
| 485 | |
---|
| 486 | INTEGER FUNCTION write_domcfg(name,Grid) |
---|
| 487 | !----------------------------------------- |
---|
| 488 | ! It creates a domain_cfg.nc used in NEMO4 |
---|
| 489 | !----------------------------------------- |
---|
| 490 | ! |
---|
| 491 | USE io_netcdf |
---|
| 492 | ! |
---|
| 493 | CHARACTER(*) name |
---|
| 494 | TYPE(Coordinates) :: Grid |
---|
| 495 | ! |
---|
| 496 | INTEGER :: status, ncid |
---|
| 497 | INTEGER :: nx, ny, jk |
---|
| 498 | INTEGER :: ln_sco, ln_isfcav, ln_zco, ln_zps, jperio |
---|
| 499 | REAL*8 :: rpi, rad, rday, rsiyea, rsiday, omega |
---|
| 500 | ! |
---|
| 501 | CHARACTER(len=1), DIMENSION(3) :: dimnames |
---|
[10381] | 502 | REAL*8 , DIMENSION(N) :: e3t_1d, e3w_1d, gdept_1d, gdepw_1d |
---|
| 503 | REAL*8 , ALLOCATABLE, DIMENSION(:,:) :: ff_t, ff_f, zbathy |
---|
| 504 | INTEGER, ALLOCATABLE, DIMENSION(:,:) :: bottom_level, top_level, mbathy |
---|
[10025] | 505 | REAL*8 , ALLOCATABLE, DIMENSION(:,:,:) :: e3t_0, e3u_0, e3v_0, e3f_0, e3w_0, e3uw_0, e3vw_0 |
---|
| 506 | ! |
---|
| 507 | ! size of the Grid |
---|
| 508 | nx = SIZE(Grid%bathy_meter,1) |
---|
| 509 | ny = SIZE(Grid%bathy_meter,2) |
---|
| 510 | |
---|
| 511 | ! allocate needed arrays for domain_cfg |
---|
[10381] | 512 | ALLOCATE( ff_t(nx,ny), ff_f(nx,ny), zbathy(nx,ny) ) |
---|
| 513 | ALLOCATE( bottom_level(nx,ny), top_level(nx,ny), mbathy(nx,ny) ) |
---|
[10025] | 514 | ALLOCATE( e3t_0(nx,ny,N), e3u_0 (nx,ny,N), e3v_0 (nx,ny,N), e3f_0(nx,ny,N), & |
---|
| 515 | & e3w_0(nx,ny,N), e3uw_0(nx,ny,N), e3vw_0(nx,ny,N) ) |
---|
| 516 | |
---|
| 517 | ! some physical parameters |
---|
| 518 | rpi = 3.141592653589793 |
---|
| 519 | rad = 3.141592653589793 / 180. |
---|
| 520 | rday = 24.*60.*60. |
---|
| 521 | rsiyea = 365.25 * rday * 2. * rpi / 6.283076 |
---|
| 522 | rsiday = rday / ( 1. + rday / rsiyea ) |
---|
| 523 | omega = 2. * rpi / rsiday |
---|
| 524 | |
---|
| 525 | ! Coriolis |
---|
| 526 | ff_f(:,:) = 2. * omega * SIN( rad * Grid%gphif(:,:) ) ! compute it on the sphere at f-point |
---|
| 527 | ff_t(:,:) = 2. * omega * SIN( rad * Grid%gphit(:,:) ) ! - - - at t-point |
---|
| 528 | |
---|
[10381] | 529 | ! bathy in meters |
---|
| 530 | zbathy(:,:) = Grid%bathy_meter |
---|
[10025] | 531 | |
---|
| 532 | ! vertical scale factors |
---|
[10381] | 533 | CALL zgr_z( e3t_1d, e3w_1d, gdept_1d, gdepw_1d ) |
---|
| 534 | ! DO jk = 1, N |
---|
| 535 | ! e3t_0 (:,:,jk) = e3t_1d (jk) |
---|
| 536 | ! e3u_0 (:,:,jk) = e3t_1d (jk) |
---|
| 537 | ! e3v_0 (:,:,jk) = e3t_1d (jk) |
---|
| 538 | ! e3f_0 (:,:,jk) = e3t_1d (jk) |
---|
| 539 | ! e3w_0 (:,:,jk) = e3w_1d (jk) |
---|
| 540 | ! e3uw_0 (:,:,jk) = e3w_1d (jk) |
---|
| 541 | ! e3vw_0 (:,:,jk) = e3w_1d (jk) |
---|
| 542 | ! END DO |
---|
[10025] | 543 | |
---|
| 544 | ! logicals and others |
---|
| 545 | ln_sco = 0 |
---|
| 546 | ln_isfcav = 0 |
---|
[12264] | 547 | !!$ IF( partial_steps ) THEN |
---|
[10025] | 548 | ln_zps = 1 |
---|
| 549 | ln_zco = 0 |
---|
[12264] | 550 | !!$ ELSE |
---|
| 551 | !!$ ln_zps = 0 |
---|
| 552 | !!$ ln_zco = 1 |
---|
| 553 | !!$ ENDIF |
---|
[10025] | 554 | |
---|
[12253] | 555 | CALL zgr_zps( nx, ny, gdept_1d, gdepw_1d, e3t_1d, e3w_1d, zbathy, & ! input |
---|
| 556 | & mbathy, e3t_0, e3u_0, e3v_0, e3f_0, e3w_0, e3uw_0, e3vw_0 ) ! output |
---|
| 557 | |
---|
[10381] | 558 | ! top/bottom levels |
---|
| 559 | bottom_level(:,:) = mbathy(:,:) |
---|
| 560 | top_level(:,:) = MIN( 1, mbathy(:,:) ) |
---|
| 561 | |
---|
[11206] | 562 | ! closed domain for child grid and namelist defined domain for parent grid |
---|
| 563 | IF( TRIM(name) == TRIM(parent_domcfg_updated) .OR. TRIM(name) == TRIM(parent_domcfg_out) ) THEN |
---|
| 564 | jperio = parent_jperio |
---|
| 565 | ELSE |
---|
| 566 | jperio = 0 |
---|
| 567 | ENDIF |
---|
[10025] | 568 | |
---|
[11231] | 569 | !IF( .NOT.ln_agrif_domain ) THEN |
---|
| 570 | ! bottom_level(1:nx,1 ) = 0 |
---|
| 571 | ! bottom_level(1:nx, ny) = 0 |
---|
| 572 | ! bottom_level(1 ,1:ny) = 0 |
---|
| 573 | ! bottom_level( nx,1:ny) = 0 |
---|
[10029] | 574 | |
---|
[11231] | 575 | ! top_level(1:nx,1 ) = 0 |
---|
| 576 | ! top_level(1:nx, ny) = 0 |
---|
| 577 | ! top_level(1 ,1:ny) = 0 |
---|
| 578 | ! top_level( nx,1:ny) = 0 |
---|
[10029] | 579 | |
---|
[11231] | 580 | ! zbathy(1:nx,1 ) = 0. |
---|
| 581 | ! zbathy(1:nx, ny) = 0. |
---|
| 582 | ! zbathy(1 ,1:ny) = 0. |
---|
| 583 | ! zbathy( nx,1:ny) = 0. |
---|
| 584 | !ENDIF |
---|
[10160] | 585 | |
---|
[10025] | 586 | ! ------------------- |
---|
| 587 | ! write domain_cfg.nc |
---|
| 588 | ! ------------------- |
---|
[10381] | 589 | status = nf90_create(name,IOR(NF90_64BIT_OFFSET,NF90_WRITE),ncid) |
---|
[10025] | 590 | status = nf90_close(ncid) |
---|
| 591 | ! |
---|
| 592 | ! dimensions |
---|
| 593 | dimnames = (/'x','y','z'/) |
---|
| 594 | CALL write_ncdf_dim(dimnames(1),name,nx) |
---|
| 595 | CALL write_ncdf_dim(dimnames(2),name,ny) |
---|
| 596 | CALL write_ncdf_dim(dimnames(3),name,N) |
---|
| 597 | ! |
---|
| 598 | ! variables |
---|
| 599 | CALL write_ncdf_var('nav_lon',dimnames(1:2),name,Grid%nav_lon,'float') |
---|
| 600 | CALL write_ncdf_var('nav_lat',dimnames(1:2),name,Grid%nav_lat,'float') |
---|
| 601 | CALL write_ncdf_var('nav_lev',dimnames(3) ,name,gdept_1d ,'float') |
---|
| 602 | ! |
---|
| 603 | CALL write_ncdf_var('jpiglo',name,nx ,'integer') |
---|
| 604 | CALL write_ncdf_var('jpjglo',name,ny ,'integer') |
---|
| 605 | CALL write_ncdf_var('jpkglo',name,N ,'integer') |
---|
| 606 | CALL write_ncdf_var('jperio',name,jperio,'integer') |
---|
| 607 | ! |
---|
| 608 | CALL write_ncdf_var('ln_zco' ,name,ln_zco ,'integer') |
---|
| 609 | CALL write_ncdf_var('ln_zps' ,name,ln_zps ,'integer') |
---|
| 610 | CALL write_ncdf_var('ln_sco' ,name,ln_sco ,'integer') |
---|
| 611 | CALL write_ncdf_var('ln_isfcav',name,ln_isfcav,'integer') |
---|
| 612 | |
---|
| 613 | CALL write_ncdf_var('glamt',dimnames(1:2),name,Grid%glamt,'double') |
---|
| 614 | CALL write_ncdf_var('glamu',dimnames(1:2),name,Grid%glamu,'double') |
---|
| 615 | CALL write_ncdf_var('glamv',dimnames(1:2),name,Grid%glamv,'double') |
---|
| 616 | CALL write_ncdf_var('glamf',dimnames(1:2),name,Grid%glamf,'double') |
---|
| 617 | CALL write_ncdf_var('gphit',dimnames(1:2),name,Grid%gphit,'double') |
---|
| 618 | CALL write_ncdf_var('gphiu',dimnames(1:2),name,Grid%gphiu,'double') |
---|
| 619 | CALL write_ncdf_var('gphiv',dimnames(1:2),name,Grid%gphiv,'double') |
---|
| 620 | CALL write_ncdf_var('gphif',dimnames(1:2),name,Grid%gphif,'double') |
---|
| 621 | |
---|
| 622 | CALL write_ncdf_var('e1t',dimnames(1:2),name,Grid%e1t,'double') |
---|
| 623 | CALL write_ncdf_var('e1u',dimnames(1:2),name,Grid%e1u,'double') |
---|
| 624 | CALL write_ncdf_var('e1v',dimnames(1:2),name,Grid%e1v,'double') |
---|
| 625 | CALL write_ncdf_var('e1f',dimnames(1:2),name,Grid%e1f,'double') |
---|
| 626 | CALL write_ncdf_var('e2t',dimnames(1:2),name,Grid%e2t,'double') |
---|
| 627 | CALL write_ncdf_var('e2u',dimnames(1:2),name,Grid%e2u,'double') |
---|
| 628 | CALL write_ncdf_var('e2v',dimnames(1:2),name,Grid%e2v,'double') |
---|
| 629 | CALL write_ncdf_var('e2f',dimnames(1:2),name,Grid%e2f,'double') |
---|
| 630 | |
---|
| 631 | CALL write_ncdf_var('ff_f',dimnames(1:2),name,ff_f,'double') |
---|
| 632 | CALL write_ncdf_var('ff_t',dimnames(1:2),name,ff_t,'double') |
---|
| 633 | |
---|
| 634 | CALL write_ncdf_var('e3t_1d',dimnames(3),name,e3t_1d,'double') |
---|
| 635 | CALL write_ncdf_var('e3w_1d',dimnames(3),name,e3w_1d,'double') |
---|
| 636 | |
---|
| 637 | CALL write_ncdf_var('e3t_0' ,dimnames(:),name,e3t_0 ,'double') |
---|
| 638 | CALL write_ncdf_var('e3w_0' ,dimnames(:),name,e3w_0 ,'double') |
---|
| 639 | CALL write_ncdf_var('e3u_0' ,dimnames(:),name,e3u_0 ,'double') |
---|
| 640 | CALL write_ncdf_var('e3v_0' ,dimnames(:),name,e3v_0 ,'double') |
---|
| 641 | CALL write_ncdf_var('e3f_0' ,dimnames(:),name,e3f_0 ,'double') |
---|
| 642 | CALL write_ncdf_var('e3w_0' ,dimnames(:),name,e3w_0 ,'double') |
---|
| 643 | CALL write_ncdf_var('e3uw_0',dimnames(:),name,e3uw_0,'double') |
---|
| 644 | CALL write_ncdf_var('e3vw_0',dimnames(:),name,e3vw_0,'double') |
---|
| 645 | |
---|
| 646 | CALL write_ncdf_var('bottom_level',dimnames(1:2),name,bottom_level,'integer') |
---|
| 647 | CALL write_ncdf_var('top_level' ,dimnames(1:2),name,top_level ,'integer') |
---|
[10160] | 648 | CALL write_ncdf_var('bathy_meter' ,dimnames(1:2),name,zbathy ,'float') |
---|
[10025] | 649 | |
---|
| 650 | ! some attributes |
---|
| 651 | CALL copy_ncdf_att('nav_lon',TRIM(parent_coordinate_file),name,MINVAL(Grid%nav_lon),MAXVAL(Grid%nav_lon)) |
---|
| 652 | CALL copy_ncdf_att('nav_lat',TRIM(parent_coordinate_file),name,MINVAL(Grid%nav_lat),MAXVAL(Grid%nav_lat)) |
---|
| 653 | |
---|
| 654 | CALL copy_ncdf_att('glamt',TRIM(parent_coordinate_file),name) |
---|
| 655 | CALL copy_ncdf_att('glamu',TRIM(parent_coordinate_file),name) |
---|
| 656 | CALL copy_ncdf_att('glamv',TRIM(parent_coordinate_file),name) |
---|
| 657 | CALL copy_ncdf_att('glamf',TRIM(parent_coordinate_file),name) |
---|
| 658 | CALL copy_ncdf_att('gphit',TRIM(parent_coordinate_file),name) |
---|
| 659 | CALL copy_ncdf_att('gphiu',TRIM(parent_coordinate_file),name) |
---|
| 660 | CALL copy_ncdf_att('gphiv',TRIM(parent_coordinate_file),name) |
---|
| 661 | CALL copy_ncdf_att('gphif',TRIM(parent_coordinate_file),name) |
---|
| 662 | |
---|
| 663 | CALL copy_ncdf_att('e1t',TRIM(parent_coordinate_file),name) |
---|
| 664 | CALL copy_ncdf_att('e1u',TRIM(parent_coordinate_file),name) |
---|
| 665 | CALL copy_ncdf_att('e1v',TRIM(parent_coordinate_file),name) |
---|
| 666 | CALL copy_ncdf_att('e1f',TRIM(parent_coordinate_file),name) |
---|
| 667 | CALL copy_ncdf_att('e2t',TRIM(parent_coordinate_file),name) |
---|
| 668 | CALL copy_ncdf_att('e2u',TRIM(parent_coordinate_file),name) |
---|
| 669 | CALL copy_ncdf_att('e2v',TRIM(parent_coordinate_file),name) |
---|
| 670 | CALL copy_ncdf_att('e2f',TRIM(parent_coordinate_file),name) |
---|
| 671 | ! |
---|
| 672 | ! control print |
---|
| 673 | WRITE(*,*) ' ' |
---|
| 674 | WRITE(*,*) 'Writing domcfg file: ',name |
---|
| 675 | WRITE(*,*) ' ' |
---|
| 676 | ! |
---|
[10381] | 677 | DEALLOCATE( ff_t, ff_f, zbathy ) |
---|
| 678 | DEALLOCATE( bottom_level, top_level, mbathy ) |
---|
[10025] | 679 | DEALLOCATE( e3t_0, e3u_0, e3v_0, e3f_0, e3w_0, e3uw_0, e3vw_0 ) |
---|
| 680 | ! |
---|
| 681 | write_domcfg = 1 |
---|
| 682 | |
---|
| 683 | END FUNCTION write_domcfg |
---|
| 684 | ! |
---|
[10381] | 685 | SUBROUTINE zgr_z(e3t_1d, e3w_1d, gdept_1d, gdepw_1d ) |
---|
[10025] | 686 | !!---------------------------------------------------------------------- |
---|
| 687 | !! *** ROUTINE zgr_z (from NEMO4) *** |
---|
| 688 | !! |
---|
| 689 | !! ** Purpose : set the depth of model levels and the resulting |
---|
| 690 | !! vertical scale factors. |
---|
| 691 | !! |
---|
| 692 | !! ** Method : z-coordinate system (use in all type of coordinate) |
---|
| 693 | !! The depth of model levels is defined from an analytical |
---|
| 694 | !! function the derivative of which gives the scale factors. |
---|
| 695 | !! both depth and scale factors only depend on k (1d arrays). |
---|
| 696 | !! w-level: gdepw_1d = gdep(k) |
---|
| 697 | !! e3w_1d(k) = dk(gdep)(k) = e3(k) |
---|
| 698 | !! t-level: gdept_1d = gdep(k+0.5) |
---|
| 699 | !! e3t_1d(k) = dk(gdep)(k+0.5) = e3(k+0.5) |
---|
| 700 | !! |
---|
| 701 | !! ** Action : - gdept_1d, gdepw_1d : depth of T- and W-point (m) |
---|
| 702 | !! - e3t_1d , e3w_1d : scale factors at T- and W-levels (m) |
---|
| 703 | !! |
---|
| 704 | !! Reference : Marti, Madec & Delecluse, 1992, JGR, 97, No8, 12,763-12,766. |
---|
| 705 | !!---------------------------------------------------------------------- |
---|
| 706 | INTEGER :: jk ! dummy loop indices |
---|
| 707 | REAL*8 :: zt, zw ! temporary scalars |
---|
| 708 | REAL*8 :: zsur, za0, za1, zkth ! Values set from parameters in |
---|
| 709 | REAL*8 :: zacr, zdzmin, zhmax ! par_CONFIG_Rxx.h90 |
---|
| 710 | REAL*8 :: za2, zkth2, zacr2 ! Values for optional double tanh function set from parameters |
---|
| 711 | |
---|
[10381] | 712 | REAL*8, DIMENSION(:), INTENT(out) :: e3t_1d, e3w_1d, gdept_1d, gdepw_1d |
---|
[10025] | 713 | !!---------------------------------------------------------------------- |
---|
| 714 | ! |
---|
| 715 | ! |
---|
| 716 | ! Set variables from parameters |
---|
| 717 | ! ------------------------------ |
---|
| 718 | zkth = ppkth ; zacr = ppacr |
---|
| 719 | zdzmin = ppdzmin ; zhmax = pphmax |
---|
| 720 | zkth2 = ppkth2 ; zacr2 = ppacr2 ! optional (ldbletanh=T) double tanh parameters |
---|
| 721 | |
---|
| 722 | ! If pa1 and pa0 and psur are et to pp_to_be_computed |
---|
| 723 | ! za0, za1, zsur are computed from ppdzmin , pphmax, ppkth, ppacr |
---|
| 724 | IF ( ( pa0 == 0 .OR. pa1 == 0 .OR. psur == 0 ) & |
---|
| 725 | .AND. ppdzmin.NE.0 .AND. pphmax.NE.0 ) THEN |
---|
| 726 | ! |
---|
| 727 | za1 = ( ppdzmin - pphmax / FLOAT(N-1) ) & |
---|
| 728 | & / ( TANH((1-ppkth)/ppacr) - ppacr/FLOAT(N-1) * ( LOG( COSH( (N - ppkth) / ppacr) ) & |
---|
| 729 | & - LOG( COSH( ( 1 - ppkth) / ppacr) ) ) ) |
---|
| 730 | za0 = ppdzmin - za1 * TANH( (1-ppkth) / ppacr ) |
---|
| 731 | zsur = - za0 - za1 * ppacr * LOG( COSH( (1-ppkth) / ppacr ) ) |
---|
| 732 | ELSE |
---|
| 733 | za1 = pa1 ; za0 = pa0 ; zsur = psur |
---|
| 734 | za2 = pa2 ! optional (ldbletanh=T) double tanh parameter |
---|
| 735 | ENDIF |
---|
| 736 | |
---|
| 737 | ! Reference z-coordinate (depth - scale factor at T- and W-points) |
---|
| 738 | ! ====================== |
---|
| 739 | IF( ppkth == 0. ) THEN ! uniform vertical grid |
---|
| 740 | |
---|
| 741 | za1 = zhmax / FLOAT(N-1) |
---|
| 742 | |
---|
| 743 | DO jk = 1, N |
---|
| 744 | zw = FLOAT( jk ) |
---|
| 745 | zt = FLOAT( jk ) + 0.5 |
---|
| 746 | gdepw_1d(jk) = ( zw - 1 ) * za1 |
---|
| 747 | gdept_1d(jk) = ( zt - 1 ) * za1 |
---|
| 748 | e3w_1d (jk) = za1 |
---|
| 749 | e3t_1d (jk) = za1 |
---|
| 750 | END DO |
---|
| 751 | ELSE ! Madec & Imbard 1996 function |
---|
| 752 | IF( .NOT. ldbletanh ) THEN |
---|
| 753 | DO jk = 1, N |
---|
| 754 | zw = REAL( jk ) |
---|
| 755 | zt = REAL( jk ) + 0.5 |
---|
| 756 | gdepw_1d(jk) = ( zsur + za0 * zw + za1 * zacr * LOG ( COSH( (zw-zkth) / zacr ) ) ) |
---|
| 757 | gdept_1d(jk) = ( zsur + za0 * zt + za1 * zacr * LOG ( COSH( (zt-zkth) / zacr ) ) ) |
---|
| 758 | e3w_1d (jk) = za0 + za1 * TANH( (zw-zkth) / zacr ) |
---|
| 759 | e3t_1d (jk) = za0 + za1 * TANH( (zt-zkth) / zacr ) |
---|
| 760 | END DO |
---|
| 761 | ELSE |
---|
| 762 | DO jk = 1, N |
---|
| 763 | zw = FLOAT( jk ) |
---|
| 764 | zt = FLOAT( jk ) + 0.5 |
---|
| 765 | ! Double tanh function |
---|
| 766 | gdepw_1d(jk) = ( zsur + za0 * zw + za1 * zacr * LOG ( COSH( (zw-zkth ) / zacr ) ) & |
---|
| 767 | & + za2 * zacr2* LOG ( COSH( (zw-zkth2) / zacr2 ) ) ) |
---|
| 768 | gdept_1d(jk) = ( zsur + za0 * zt + za1 * zacr * LOG ( COSH( (zt-zkth ) / zacr ) ) & |
---|
| 769 | & + za2 * zacr2* LOG ( COSH( (zt-zkth2) / zacr2 ) ) ) |
---|
| 770 | e3w_1d (jk) = za0 + za1 * TANH( (zw-zkth ) / zacr ) & |
---|
| 771 | & + za2 * TANH( (zw-zkth2) / zacr2 ) |
---|
| 772 | e3t_1d (jk) = za0 + za1 * TANH( (zt-zkth ) / zacr ) & |
---|
| 773 | & + za2 * TANH( (zt-zkth2) / zacr2 ) |
---|
| 774 | END DO |
---|
| 775 | ENDIF |
---|
| 776 | gdepw_1d(1) = 0. ! force first w-level to be exactly at zero |
---|
| 777 | ENDIF |
---|
| 778 | |
---|
| 779 | IF ( ln_e3_dep ) THEN ! e3. = dk[gdep] |
---|
| 780 | ! |
---|
| 781 | !==>>> need to be like this to compute the pressure gradient with ISF. |
---|
| 782 | ! If not, level beneath the ISF are not aligned (sum(e3t) /= depth) |
---|
| 783 | ! define e3t_0 and e3w_0 as the differences between gdept and gdepw respectively |
---|
| 784 | ! |
---|
| 785 | DO jk = 1, N-1 |
---|
| 786 | e3t_1d(jk) = gdepw_1d(jk+1)-gdepw_1d(jk) |
---|
| 787 | END DO |
---|
| 788 | e3t_1d(N) = e3t_1d(N-1) ! we don't care because this level is masked in NEMO |
---|
| 789 | |
---|
| 790 | DO jk = 2, N |
---|
| 791 | e3w_1d(jk) = gdept_1d(jk) - gdept_1d(jk-1) |
---|
| 792 | END DO |
---|
| 793 | e3w_1d(1 ) = 2. * (gdept_1d(1) - gdepw_1d(1)) |
---|
| 794 | END IF |
---|
| 795 | |
---|
| 796 | ! |
---|
| 797 | END SUBROUTINE zgr_z |
---|
| 798 | |
---|
[10381] | 799 | SUBROUTINE zgr_zps( nx, ny, gdept_1d, gdepw_1d, e3t_1d, e3w_1d, zbathy, & |
---|
| 800 | & mbathy, e3t_0, e3u_0, e3v_0, e3f_0, e3w_0, e3uw_0, e3vw_0 ) |
---|
| 801 | !!---------------------------------------------------------------------- |
---|
| 802 | !! *** ROUTINE zgr_zps (from NEMO4) *** |
---|
| 803 | !! |
---|
| 804 | !! ** Purpose : the depth and vertical scale factor in partial step |
---|
| 805 | !! reference z-coordinate case |
---|
| 806 | !! |
---|
| 807 | !! ** Method : Partial steps : computes the 3D vertical scale factors |
---|
| 808 | !! of T-, U-, V-, W-, UW-, VW and F-points that are associated with |
---|
| 809 | !! a partial step representation of bottom topography. |
---|
| 810 | !! |
---|
| 811 | !! The reference depth of model levels is defined from an analytical |
---|
| 812 | !! function the derivative of which gives the reference vertical |
---|
| 813 | !! scale factors. |
---|
| 814 | !! From depth and scale factors reference, we compute there new value |
---|
| 815 | !! with partial steps on 3d arrays ( i, j, k ). |
---|
| 816 | !! |
---|
| 817 | !! w-level: gdepw_0(i,j,k) = gdep(k) |
---|
| 818 | !! e3w_0(i,j,k) = dk(gdep)(k) = e3(i,j,k) |
---|
| 819 | !! t-level: gdept_0(i,j,k) = gdep(k+0.5) |
---|
| 820 | !! e3t_0(i,j,k) = dk(gdep)(k+0.5) = e3(i,j,k+0.5) |
---|
| 821 | !! |
---|
| 822 | !! With the help of the bathymetric file ( bathymetry_depth_ORCA_R2.nc), |
---|
| 823 | !! we find the mbathy index of the depth at each grid point. |
---|
| 824 | !! This leads us to three cases: |
---|
| 825 | !! |
---|
| 826 | !! - bathy = 0 => mbathy = 0 |
---|
| 827 | !! - 1 < mbathy < jpkm1 |
---|
| 828 | !! - bathy > gdepw_0(jpk) => mbathy = jpkm1 |
---|
| 829 | !! |
---|
| 830 | !! Then, for each case, we find the new depth at t- and w- levels |
---|
| 831 | !! and the new vertical scale factors at t-, u-, v-, w-, uw-, vw- |
---|
| 832 | !! and f-points. |
---|
| 833 | !! |
---|
| 834 | !! This routine is given as an example, it must be modified |
---|
| 835 | !! following the user s desiderata. nevertheless, the output as |
---|
| 836 | !! well as the way to compute the model levels and scale factors |
---|
| 837 | !! must be respected in order to insure second order accuracy |
---|
| 838 | !! schemes. |
---|
| 839 | !! |
---|
| 840 | !! c a u t i o n : gdept_1d, gdepw_1d and e3._1d are positives |
---|
| 841 | !! - - - - - - - gdept_0, gdepw_0 and e3. are positives |
---|
| 842 | !! |
---|
| 843 | !! Reference : Pacanowsky & Gnanadesikan 1997, Mon. Wea. Rev., 126, 3248-3270. |
---|
| 844 | !!---------------------------------------------------------------------- |
---|
| 845 | INTEGER :: ji, jj, jk ! dummy loop indices |
---|
| 846 | INTEGER :: ik |
---|
| 847 | INTEGER :: ikb, ikt ! temporary integers |
---|
| 848 | REAL*8 :: ze3tp , ze3wp ! Last ocean level thickness at T- and W-points |
---|
| 849 | REAL*8 :: zdepwp, zdepth ! Ajusted ocean depth to avoid too small e3t |
---|
| 850 | REAL*8 :: zmax ! temporary scalar |
---|
| 851 | ! |
---|
| 852 | REAL*8, ALLOCATABLE, DIMENSION(:,:,:) :: gdept_0, gdepw_0, zprt |
---|
| 853 | ! |
---|
| 854 | INTEGER, INTENT(in ) :: nx, ny |
---|
| 855 | REAL*8 , DIMENSION(:) , INTENT(in ) :: gdept_1d, gdepw_1d, e3t_1d, e3w_1d |
---|
| 856 | REAL*8 , DIMENSION(:,:) , INTENT(inout) :: zbathy |
---|
| 857 | INTEGER, DIMENSION(:,:) , INTENT( out) :: mbathy |
---|
| 858 | REAL*8 , DIMENSION(:,:,:), INTENT( out) :: e3t_0, e3u_0, e3v_0, e3f_0, e3w_0, e3uw_0, e3vw_0 |
---|
| 859 | ! |
---|
| 860 | !!--------------------------------------------------------------------- |
---|
| 861 | |
---|
| 862 | ALLOCATE( zprt(nx,ny,N), gdept_0(nx,ny,N), gdepw_0(nx,ny,N) ) |
---|
| 863 | ! |
---|
| 864 | ! bathymetry in level (from bathy_meter) |
---|
| 865 | ! =================== |
---|
| 866 | zmax = gdepw_1d(N) + e3t_1d(N) ! maximum depth (i.e. the last ocean level thickness <= 2*e3t_1d(N-1) ) |
---|
| 867 | zbathy(:,:) = MIN( zmax , zbathy(:,:) ) ! bounded value of bathy (min already set at the end of zgr_bat) |
---|
[12264] | 868 | WHERE( zbathy(:,:) == 0. ) ; mbathy(:,:) = 0 ! land : set mbathy to 0 |
---|
[10381] | 869 | ELSEWHERE ; mbathy(:,:) = N-1 ! ocean : initialize mbathy to the max ocean level |
---|
| 870 | END WHERE |
---|
[10025] | 871 | |
---|
[12264] | 872 | !!$ IF( partial_steps ) THEN |
---|
[12253] | 873 | ! Compute mbathy for ocean points (i.e. the number of ocean levels) |
---|
| 874 | ! find the number of ocean levels such that the last level thickness |
---|
| 875 | ! is larger than the minimum of e3zps_min and e3zps_rat * e3t_1d (where |
---|
| 876 | ! e3t_1d is the reference level thickness |
---|
| 877 | DO jk = N-1, 1, -1 |
---|
| 878 | zdepth = gdepw_1d(jk) + MIN( e3zps_min, e3t_1d(jk)*e3zps_rat ) |
---|
| 879 | WHERE( 0. < zbathy(:,:) .AND. zbathy(:,:) <= zdepth ) mbathy(:,:) = jk-1 |
---|
| 880 | END DO |
---|
[12264] | 881 | !!$ ELSE |
---|
| 882 | !!$ DO jk = 1, N |
---|
| 883 | !!$ WHERE( 0. < zbathy(:,:) .AND. zbathy(:,:) >= gdepw_1d(jk) ) mbathy(:,:) = jk-1 |
---|
| 884 | !!$ END DO |
---|
| 885 | !!$ ENDIF |
---|
[12253] | 886 | |
---|
[10381] | 887 | ! Scale factors and depth at T- and W-points |
---|
| 888 | DO jk = 1, N ! intitialization to the reference z-coordinate |
---|
| 889 | gdept_0(:,:,jk) = gdept_1d(jk) |
---|
| 890 | gdepw_0(:,:,jk) = gdepw_1d(jk) |
---|
| 891 | e3t_0 (:,:,jk) = e3t_1d (jk) |
---|
| 892 | e3w_0 (:,:,jk) = e3w_1d (jk) |
---|
| 893 | END DO |
---|
| 894 | |
---|
| 895 | ! Bathy, iceshelf draft, scale factor and depth at T- and W- points in case of isf |
---|
| 896 | !!clem: not implemented yet |
---|
| 897 | !! IF ( ln_isfcav == 1 ) CALL zgr_isf |
---|
| 898 | ! |
---|
[12264] | 899 | ! Scale factors and depth at T- and W-points |
---|
| 900 | ! IF ( ln_isfcav == 0 ) THEN |
---|
| 901 | DO jj = 1, ny |
---|
| 902 | DO ji = 1, nx |
---|
| 903 | ik = mbathy(ji,jj) |
---|
| 904 | IF( ik > 0 ) THEN ! ocean point only |
---|
| 905 | ! max ocean level case |
---|
| 906 | IF( ik == N-1 ) THEN |
---|
| 907 | zdepwp = zbathy(ji,jj) |
---|
| 908 | ze3tp = zbathy(ji,jj) - gdepw_1d(ik) |
---|
| 909 | ze3wp = 0.5 * e3w_1d(ik) * ( 1. + ( ze3tp/e3t_1d(ik) ) ) |
---|
| 910 | e3t_0(ji,jj,ik ) = ze3tp |
---|
| 911 | e3t_0(ji,jj,ik+1) = ze3tp |
---|
| 912 | e3w_0(ji,jj,ik ) = ze3wp |
---|
| 913 | e3w_0(ji,jj,ik+1) = ze3tp |
---|
| 914 | gdepw_0(ji,jj,ik+1) = zdepwp |
---|
| 915 | gdept_0(ji,jj,ik ) = gdept_1d(ik-1) + ze3wp |
---|
| 916 | gdept_0(ji,jj,ik+1) = gdept_0(ji,jj,ik) + ze3tp |
---|
| 917 | ! |
---|
| 918 | ELSE ! standard case |
---|
| 919 | IF( zbathy(ji,jj) <= gdepw_1d(ik+1) ) THEN ; gdepw_0(ji,jj,ik+1) = zbathy(ji,jj) |
---|
| 920 | ELSE ; gdepw_0(ji,jj,ik+1) = gdepw_1d(ik+1) |
---|
[10381] | 921 | ENDIF |
---|
[12264] | 922 | !gm Bug? check the gdepw_1d |
---|
| 923 | ! ... on ik |
---|
| 924 | gdept_0(ji,jj,ik) = gdepw_1d(ik) + ( gdepw_0(ji,jj,ik+1) - gdepw_1d(ik) ) & |
---|
| 925 | & * ((gdept_1d( ik ) - gdepw_1d(ik) ) & |
---|
| 926 | & / ( gdepw_1d( ik+1) - gdepw_1d(ik) )) |
---|
| 927 | e3t_0 (ji,jj,ik) = e3t_1d (ik) * ( gdepw_0 (ji,jj,ik+1) - gdepw_1d(ik) ) & |
---|
| 928 | & / ( gdepw_1d( ik+1) - gdepw_1d(ik) ) |
---|
| 929 | e3w_0(ji,jj,ik) = 0.5 * ( gdepw_0(ji,jj,ik+1) + gdepw_1d(ik+1) - 2. * gdepw_1d(ik) ) & |
---|
| 930 | & * ( e3w_1d(ik) / ( gdepw_1d(ik+1) - gdepw_1d(ik) ) ) |
---|
| 931 | ! ... on ik+1 |
---|
| 932 | e3w_0 (ji,jj,ik+1) = e3t_0 (ji,jj,ik) |
---|
| 933 | e3t_0 (ji,jj,ik+1) = e3t_0 (ji,jj,ik) |
---|
| 934 | gdept_0(ji,jj,ik+1) = gdept_0(ji,jj,ik) + e3t_0(ji,jj,ik) |
---|
[10381] | 935 | ENDIF |
---|
[12264] | 936 | ENDIF |
---|
[10381] | 937 | END DO |
---|
[12264] | 938 | END DO |
---|
[10381] | 939 | ! |
---|
[12264] | 940 | ! DO jj = 1, ny |
---|
| 941 | ! DO ji = 1, nx |
---|
| 942 | ! ik = mbathy(ji,jj) |
---|
| 943 | ! IF( ik > 0 ) THEN ! ocean point only |
---|
| 944 | ! e3tp (ji,jj) = e3t_0(ji,jj,ik) |
---|
| 945 | ! e3wp (ji,jj) = e3w_0(ji,jj,ik) |
---|
| 946 | ! ENDIF |
---|
| 947 | ! END DO |
---|
| 948 | ! END DO |
---|
| 949 | ! END IF |
---|
| 950 | ! |
---|
[10381] | 951 | ! Scale factors and depth at U-, V-, UW and VW-points |
---|
| 952 | DO jk = 1, N ! initialisation to z-scale factors |
---|
| 953 | e3u_0 (:,:,jk) = e3t_1d(jk) |
---|
| 954 | e3v_0 (:,:,jk) = e3t_1d(jk) |
---|
| 955 | e3uw_0(:,:,jk) = e3w_1d(jk) |
---|
| 956 | e3vw_0(:,:,jk) = e3w_1d(jk) |
---|
| 957 | END DO |
---|
| 958 | |
---|
| 959 | DO jk = 1,N ! Computed as the minimum of neighbooring scale factors |
---|
| 960 | DO jj = 1, ny-1 |
---|
| 961 | DO ji = 1, nx-1 ! vector opt. |
---|
| 962 | e3u_0 (ji,jj,jk) = MIN( e3t_0(ji,jj,jk), e3t_0(ji+1,jj,jk) ) |
---|
| 963 | e3v_0 (ji,jj,jk) = MIN( e3t_0(ji,jj,jk), e3t_0(ji,jj+1,jk) ) |
---|
| 964 | e3uw_0(ji,jj,jk) = MIN( e3w_0(ji,jj,jk), e3w_0(ji+1,jj,jk) ) |
---|
| 965 | e3vw_0(ji,jj,jk) = MIN( e3w_0(ji,jj,jk), e3w_0(ji,jj+1,jk) ) |
---|
| 966 | END DO |
---|
| 967 | END DO |
---|
| 968 | END DO |
---|
| 969 | |
---|
| 970 | ! IF ( ln_isfcav == 1 ) THEN |
---|
| 971 | ! ! (ISF) define e3uw (adapted for 2 cells in the water column) |
---|
| 972 | ! DO jj = 2, ny-1 |
---|
| 973 | ! DO ji = 2, nx-1 ! vector opt. |
---|
| 974 | ! ikb = MAX(mbathy (ji,jj),mbathy (ji+1,jj)) |
---|
| 975 | ! ikt = MAX(misfdep(ji,jj),misfdep(ji+1,jj)) |
---|
| 976 | ! IF (ikb == ikt+1) e3uw_0(ji,jj,ikb) = MIN( gdept_0(ji,jj,ikb ), gdept_0(ji+1,jj ,ikb ) ) & |
---|
| 977 | ! & - MAX( gdept_0(ji,jj,ikb-1), gdept_0(ji+1,jj ,ikb-1) ) |
---|
| 978 | ! ikb = MAX(mbathy (ji,jj),mbathy (ji,jj+1)) |
---|
| 979 | ! ikt = MAX(misfdep(ji,jj),misfdep(ji,jj+1)) |
---|
| 980 | ! IF (ikb == ikt+1) e3vw_0(ji,jj,ikb) = MIN( gdept_0(ji,jj,ikb ), gdept_0(ji ,jj+1,ikb ) ) & |
---|
| 981 | ! & - MAX( gdept_0(ji,jj,ikb-1), gdept_0(ji ,jj+1,ikb-1) ) |
---|
| 982 | ! END DO |
---|
| 983 | ! END DO |
---|
| 984 | ! END IF |
---|
| 985 | ! |
---|
| 986 | |
---|
| 987 | DO jk = 1, N ! set to z-scale factor if zero (i.e. along closed boundaries) |
---|
| 988 | WHERE( e3u_0 (:,:,jk) == 0. ) e3u_0 (:,:,jk) = e3t_1d(jk) |
---|
| 989 | WHERE( e3v_0 (:,:,jk) == 0. ) e3v_0 (:,:,jk) = e3t_1d(jk) |
---|
| 990 | WHERE( e3uw_0(:,:,jk) == 0. ) e3uw_0(:,:,jk) = e3w_1d(jk) |
---|
| 991 | WHERE( e3vw_0(:,:,jk) == 0. ) e3vw_0(:,:,jk) = e3w_1d(jk) |
---|
| 992 | END DO |
---|
| 993 | |
---|
| 994 | ! Scale factor at F-point |
---|
| 995 | DO jk = 1, N ! initialisation to z-scale factors |
---|
| 996 | e3f_0(:,:,jk) = e3t_1d(jk) |
---|
| 997 | END DO |
---|
| 998 | DO jk = 1, N ! Computed as the minimum of neighbooring V-scale factors |
---|
| 999 | DO jj = 1, ny-1 |
---|
| 1000 | DO ji = 1, nx-1 ! vector opt. |
---|
| 1001 | e3f_0(ji,jj,jk) = MIN( e3v_0(ji,jj,jk), e3v_0(ji+1,jj,jk) ) |
---|
| 1002 | END DO |
---|
| 1003 | END DO |
---|
| 1004 | END DO |
---|
| 1005 | ! |
---|
| 1006 | DO jk = 1, N ! set to z-scale factor if zero (i.e. along closed boundaries) |
---|
| 1007 | WHERE( e3f_0(:,:,jk) == 0. ) e3f_0(:,:,jk) = e3t_1d(jk) |
---|
| 1008 | END DO |
---|
| 1009 | |
---|
| 1010 | ! Control of the sign |
---|
| 1011 | IF( MINVAL( e3t_0 (:,:,:) ) <= 0. ) STOP ' zgr_zps : e r r o r e3t_0 <= 0' |
---|
| 1012 | IF( MINVAL( e3w_0 (:,:,:) ) <= 0. ) STOP ' zgr_zps : e r r o r e3w_0 <= 0' |
---|
| 1013 | IF( MINVAL( gdept_0(:,:,:) ) < 0. ) STOP ' zgr_zps : e r r o r gdept_0 < 0' |
---|
| 1014 | IF( MINVAL( gdepw_0(:,:,:) ) < 0. ) STOP ' zgr_zps : e r r o r gdepw_0 < 0' |
---|
| 1015 | |
---|
| 1016 | ! Compute gde3w_0 (vertical sum of e3w) |
---|
| 1017 | ! IF ( ln_isfcav ==1 ) THEN ! if cavity |
---|
| 1018 | ! WHERE( misfdep == 0 ) misfdep = 1 |
---|
| 1019 | ! DO jj = 1,ny |
---|
| 1020 | ! DO ji = 1,nx |
---|
| 1021 | ! gde3w_0(ji,jj,1) = 0.5 * e3w_0(ji,jj,1) |
---|
| 1022 | ! DO jk = 2, misfdep(ji,jj) |
---|
| 1023 | ! gde3w_0(ji,jj,jk) = gde3w_0(ji,jj,jk-1) + e3w_0(ji,jj,jk) |
---|
| 1024 | ! END DO |
---|
| 1025 | ! IF( misfdep(ji,jj) >= 2 ) gde3w_0(ji,jj,misfdep(ji,jj)) = risfdep(ji,jj) + 0.5 * e3w_0(ji,jj,misfdep(ji,jj)) |
---|
| 1026 | ! DO jk = misfdep(ji,jj) + 1, N |
---|
| 1027 | ! gde3w_0(ji,jj,jk) = gde3w_0(ji,jj,jk-1) + e3w_0(ji,jj,jk) |
---|
| 1028 | ! END DO |
---|
| 1029 | ! END DO |
---|
| 1030 | ! END DO |
---|
| 1031 | ! ELSE ! no cavity |
---|
| 1032 | ! gde3w_0(:,:,1) = 0.5 * e3w_0(:,:,1) |
---|
| 1033 | ! DO jk = 2, N |
---|
| 1034 | ! gde3w_0(:,:,jk) = gde3w_0(:,:,jk-1) + e3w_0(:,:,jk) |
---|
| 1035 | ! END DO |
---|
| 1036 | ! END IF |
---|
| 1037 | ! |
---|
| 1038 | DEALLOCATE( zprt, gdept_0, gdepw_0 ) |
---|
| 1039 | ! |
---|
| 1040 | END SUBROUTINE zgr_zps |
---|
| 1041 | |
---|
[10025] | 1042 | !***************************************************** |
---|
[2136] | 1043 | ! function set_child_name(Parentname,Childname) |
---|
| 1044 | !***************************************************** |
---|
| 1045 | ! |
---|
| 1046 | SUBROUTINE set_child_name(Parentname,Childname) |
---|
| 1047 | ! |
---|
| 1048 | CHARACTER(*),INTENT(in) :: Parentname |
---|
| 1049 | CHARACTER(*),INTENT(out) :: Childname |
---|
| 1050 | CHARACTER(2) :: prefix |
---|
| 1051 | INTEGER :: pos |
---|
| 1052 | ! |
---|
| 1053 | pos = INDEX(TRIM(Parentname),'/',back=.TRUE.) |
---|
| 1054 | ! |
---|
| 1055 | prefix=Parentname(pos+1:pos+2) |
---|
| 1056 | IF (prefix == '1_') THEN |
---|
| 1057 | Childname = '2_'//Parentname(pos+3:LEN(Parentname)) |
---|
| 1058 | ELSEIF (prefix == '2_') THEN |
---|
| 1059 | Childname = '3_'//Parentname(pos+3:LEN(Parentname)) |
---|
| 1060 | ELSEIF (prefix == '3_') THEN |
---|
| 1061 | Childname = '4_'//Parentname(pos+3:LEN(Parentname)) |
---|
| 1062 | ELSEIF (prefix == '4_') THEN |
---|
| 1063 | Childname = '5_'//Parentname(pos+3:LEN(Parentname)) |
---|
| 1064 | ELSE |
---|
| 1065 | Childname = '1_'//Parentname(pos+1:LEN(Parentname)) |
---|
| 1066 | ENDIF |
---|
| 1067 | ! |
---|
| 1068 | END SUBROUTINE set_child_name |
---|
| 1069 | ! |
---|
| 1070 | !***************************************************** |
---|
| 1071 | ! subroutine get_interptype(varname,interp_type,conservation) |
---|
| 1072 | !***************************************************** |
---|
| 1073 | ! |
---|
| 1074 | SUBROUTINE get_interptype( varname,interp_type,conservation ) |
---|
| 1075 | ! |
---|
| 1076 | LOGICAL,OPTIONAL :: conservation |
---|
| 1077 | CHARACTER(*) :: interp_type,varname |
---|
| 1078 | INTEGER :: pos,pos1,pos2,pos3,i,k |
---|
| 1079 | LOGICAL :: find |
---|
| 1080 | i=1 |
---|
[10381] | 1081 | DO WHILE ( TRIM(VAR_INTERP(i)) /= '' ) |
---|
[2136] | 1082 | pos = INDEX( TRIM(VAR_INTERP(i)) , TRIM(varname) ) |
---|
| 1083 | IF ( pos .NE. 0 ) THEN |
---|
| 1084 | pos1 = INDEX( TRIM(VAR_INTERP(i)) , 'bicubic' ) |
---|
| 1085 | pos2 = INDEX( TRIM(VAR_INTERP(i)) , 'bilinear' ) |
---|
| 1086 | pos3 = INDEX( TRIM(VAR_INTERP(i)) , 'conservative' ) |
---|
| 1087 | ! initialize interp_type |
---|
| 1088 | IF( pos1 .NE. 0 ) interp_type = 'bicubic' |
---|
| 1089 | IF( pos2 .NE. 0 ) interp_type = 'bilinear' |
---|
| 1090 | IF( pos1 .EQ. 0 .AND. pos2 .EQ. 0) interp_type = 'bicubic' |
---|
| 1091 | ! initialize conservation |
---|
| 1092 | IF( pos3 .NE. 0 .AND. PRESENT(conservation) ) THEN |
---|
| 1093 | conservation = .TRUE. |
---|
| 1094 | RETURN |
---|
| 1095 | ELSE |
---|
| 1096 | conservation = .FALSE. |
---|
| 1097 | ENDIF |
---|
| 1098 | find = .FALSE. |
---|
| 1099 | IF( PRESENT(conservation) ) THEN |
---|
| 1100 | k=0 |
---|
| 1101 | conservation = .FALSE. |
---|
[2455] | 1102 | DO WHILE( k < SIZE(flxtab) .AND. .NOT.find ) |
---|
[2136] | 1103 | k = k+1 |
---|
| 1104 | IF( TRIM(varname) .EQ. TRIM(flxtab(k)) ) THEN |
---|
| 1105 | conservation = .TRUE. |
---|
| 1106 | find = .TRUE. |
---|
| 1107 | ENDIF |
---|
| 1108 | END DO |
---|
| 1109 | ENDIF |
---|
| 1110 | RETURN |
---|
| 1111 | ENDIF |
---|
| 1112 | i = i+1 |
---|
| 1113 | END DO |
---|
| 1114 | !default values interp_type = bicubic // conservation = false |
---|
| 1115 | interp_type = 'bicubic' |
---|
| 1116 | IF( PRESENT(conservation) ) conservation = .FALSE. |
---|
| 1117 | |
---|
| 1118 | RETURN |
---|
| 1119 | ! |
---|
| 1120 | END SUBROUTINE get_interptype |
---|
| 1121 | ! |
---|
| 1122 | !***************************************************** |
---|
| 1123 | ! subroutine Init_mask(name,Grid) |
---|
| 1124 | !***************************************************** |
---|
| 1125 | ! |
---|
| 1126 | SUBROUTINE Init_mask(name,Grid,jpiglo,jpjglo) |
---|
| 1127 | ! |
---|
| 1128 | USE io_netcdf |
---|
| 1129 | ! |
---|
| 1130 | CHARACTER(*) name |
---|
| 1131 | INTEGER :: nx,ny,k,i,j,jpiglo,jpjglo |
---|
| 1132 | TYPE(Coordinates) :: Grid |
---|
| 1133 | REAL*8, POINTER, DIMENSION(:,:) ::zwf => NULL() |
---|
| 1134 | ! |
---|
| 1135 | IF(jpiglo == 1 .AND. jpjglo == 1) THEN |
---|
[10025] | 1136 | CALL read_ncdf_var('Bathy_level',name,Grid%Bathy_level) |
---|
[2136] | 1137 | ELSE |
---|
[10025] | 1138 | CALL read_ncdf_var('Bathy_level',name,Grid%Bathy_level,(/jpizoom,jpjzoom/),(/jpiglo,jpjglo/) ) |
---|
[2136] | 1139 | ENDIF |
---|
| 1140 | |
---|
| 1141 | |
---|
| 1142 | ! |
---|
| 1143 | WRITE(*,*) 'Init masks in T,U,V,F points' |
---|
| 1144 | ! |
---|
| 1145 | nx = SIZE(Grid%Bathy_level,1) |
---|
| 1146 | ny = SIZE(Grid%Bathy_level,2) |
---|
| 1147 | ! |
---|
| 1148 | ! |
---|
| 1149 | ALLOCATE(Grid%tmask(nx,ny,N), & |
---|
| 1150 | Grid%umask(nx,ny,N), & |
---|
| 1151 | Grid%vmask(nx,ny,N), & |
---|
| 1152 | Grid%fmask(nx,ny,N)) |
---|
| 1153 | ! |
---|
| 1154 | DO k = 1,N |
---|
| 1155 | ! |
---|
| 1156 | WHERE(Grid%Bathy_level(:,:) <= k-1 ) |
---|
| 1157 | Grid%tmask(:,:,k) = 0 |
---|
| 1158 | ELSEWHERE |
---|
| 1159 | Grid%tmask(:,:,k) = 1 |
---|
| 1160 | END WHERE |
---|
| 1161 | ! |
---|
| 1162 | END DO |
---|
| 1163 | ! |
---|
| 1164 | Grid%umask(1:nx-1,:,:) = Grid%tmask(1:nx-1,:,:)*Grid%tmask(2:nx,:,:) |
---|
| 1165 | Grid%vmask(:,1:ny-1,:) = Grid%tmask(:,1:ny-1,:)*Grid%tmask(:,2:ny,:) |
---|
| 1166 | ! |
---|
| 1167 | Grid%umask(nx,:,:) = Grid%tmask(nx,:,:) |
---|
| 1168 | Grid%vmask(:,ny,:) = Grid%tmask(:,ny,:) |
---|
| 1169 | ! |
---|
| 1170 | Grid%fmask(1:nx-1,1:ny-1,:) = Grid%tmask(1:nx-1,1:ny-1,:)*Grid%tmask(2:nx,1:ny-1,:)* & |
---|
| 1171 | Grid%tmask(1:nx-1,2:ny,:)*Grid%tmask(2:nx,2:ny,:) |
---|
| 1172 | ! |
---|
| 1173 | Grid%fmask(nx,:,:) = Grid%tmask(nx,:,:) |
---|
| 1174 | Grid%fmask(:,ny,:) = Grid%tmask(:,ny,:) |
---|
| 1175 | ! |
---|
| 1176 | ALLOCATE(zwf(nx,ny)) |
---|
| 1177 | ! |
---|
| 1178 | DO k = 1,N |
---|
| 1179 | ! |
---|
| 1180 | zwf(:,:) = Grid%fmask(:,:,k) |
---|
| 1181 | ! |
---|
| 1182 | DO j = 2, ny-1 |
---|
| 1183 | DO i = 2,nx-1 |
---|
| 1184 | IF( Grid%fmask(i,j,k) == 0. ) THEN |
---|
| 1185 | Grid%fmask(i,j,k) = shlat * MIN(1.,MAX( zwf(i+1,j),zwf(i,j+1),zwf(i-1,j),zwf(i,j-1))) |
---|
| 1186 | END IF |
---|
| 1187 | END DO |
---|
| 1188 | END DO |
---|
| 1189 | ! |
---|
| 1190 | DO j = 2, ny-1 |
---|
| 1191 | IF( Grid%fmask(1,j,k) == 0. ) THEN |
---|
| 1192 | Grid%fmask(1,j,k) = shlat * MIN(1.,MAX(zwf(2,j),zwf(1,j+1),zwf(1,j-1))) |
---|
| 1193 | ENDIF |
---|
| 1194 | |
---|
| 1195 | IF( Grid%fmask(nx,j,k) == 0. ) THEN |
---|
| 1196 | Grid%fmask(nx,j,k) = shlat * MIN(1.,MAX(zwf(nx,j+1),zwf(nx-1,j),zwf(nx,j-1))) |
---|
| 1197 | ENDIF |
---|
| 1198 | END DO |
---|
| 1199 | ! |
---|
| 1200 | DO i = 2, nx-1 |
---|
| 1201 | IF( Grid%fmask(i,1,k) == 0. ) THEN |
---|
| 1202 | Grid%fmask(i, 1 ,k) = shlat*MIN( 1.,MAX(zwf(i+1,1),zwf(i,2),zwf(i-1,1))) |
---|
| 1203 | ENDIF |
---|
| 1204 | ! |
---|
| 1205 | IF( Grid%fmask(i,ny,k) == 0. ) THEN |
---|
| 1206 | Grid%fmask(i,ny,k) = shlat * MIN(1.,MAX(zwf(i+1,ny),zwf(i-1,ny),zwf(i,ny-1))) |
---|
| 1207 | ENDIF |
---|
| 1208 | END DO |
---|
| 1209 | !! |
---|
| 1210 | END DO |
---|
| 1211 | !! |
---|
| 1212 | END SUBROUTINE Init_mask |
---|
| 1213 | ! |
---|
| 1214 | !***************************************************** |
---|
| 1215 | ! subroutine Init_Tmask(name,Grid) |
---|
| 1216 | !***************************************************** |
---|
| 1217 | ! |
---|
| 1218 | SUBROUTINE Init_Tmask(name,Grid,jpiglo,jpjglo) |
---|
| 1219 | ! |
---|
| 1220 | USE io_netcdf |
---|
| 1221 | ! |
---|
| 1222 | CHARACTER(*) name |
---|
| 1223 | INTEGER :: nx,ny,k,i,j,jpiglo,jpjglo |
---|
| 1224 | TYPE(Coordinates) :: Grid |
---|
| 1225 | REAL*8, POINTER, DIMENSION(:,:) ::zwf => NULL() |
---|
| 1226 | ! |
---|
| 1227 | IF(jpiglo == 1 .AND. jpjglo == 1) THEN |
---|
[10025] | 1228 | CALL read_ncdf_var('Bathy_level',name,Grid%Bathy_level) |
---|
[2136] | 1229 | ELSE |
---|
[10025] | 1230 | CALL read_ncdf_var('Bathy_level',name,Grid%Bathy_level,(/jpizoom,jpjzoom/),(/jpiglo,jpjglo/) ) |
---|
[2136] | 1231 | ENDIF |
---|
| 1232 | ! |
---|
| 1233 | nx = SIZE(Grid%Bathy_level,1) |
---|
| 1234 | ny = SIZE(Grid%Bathy_level,2) |
---|
| 1235 | ! |
---|
| 1236 | WRITE(*,*) 'Init masks in T points' |
---|
| 1237 | ! |
---|
[10383] | 1238 | ALLOCATE( Grid%tmask(nx,ny,N) ) |
---|
[2136] | 1239 | ! |
---|
| 1240 | DO k = 1,N |
---|
| 1241 | ! |
---|
| 1242 | WHERE(Grid%Bathy_level(:,:) <= k-1 ) |
---|
| 1243 | Grid%tmask(:,:,k) = 0. |
---|
| 1244 | ELSEWHERE |
---|
| 1245 | Grid%tmask(:,:,k) = 1. |
---|
| 1246 | END WHERE |
---|
| 1247 | ! |
---|
| 1248 | END DO |
---|
| 1249 | ! |
---|
| 1250 | END SUBROUTINE Init_Tmask |
---|
| 1251 | ! |
---|
| 1252 | !***************************************************** |
---|
[10383] | 1253 | ! subroutine Init_ssmask(name,Grid) |
---|
| 1254 | !***************************************************** |
---|
| 1255 | ! |
---|
| 1256 | ! SUBROUTINE Init_ssmask(varname,filename,Grid,jpiglo,jpjglo) |
---|
| 1257 | ! ! |
---|
| 1258 | ! USE io_netcdf |
---|
| 1259 | ! ! |
---|
| 1260 | ! CHARACTER(*) varname,filename |
---|
| 1261 | ! INTEGER :: nx,ny,k,i,j,jpiglo,jpjglo |
---|
| 1262 | ! TYPE(Coordinates) :: Grid |
---|
| 1263 | ! REAL*8, POINTER, DIMENSION(:,:) ::zwf => NULL() |
---|
| 1264 | ! ! |
---|
| 1265 | ! IF(jpiglo == 1 .AND. jpjglo == 1) THEN |
---|
| 1266 | ! CALL read_ncdf_var(varname,filename,Grid%bathy_level) |
---|
| 1267 | ! ELSE |
---|
| 1268 | ! CALL read_ncdf_var(varname,filename,Grid%bathy_level,(/jpizoom,jpjzoom/),(/jpiglo,jpjglo/) ) |
---|
| 1269 | ! ENDIF |
---|
| 1270 | ! ! |
---|
| 1271 | ! nx = SIZE(Grid%bathy_level,1) |
---|
| 1272 | ! ny = SIZE(Grid%bathy_level,2) |
---|
| 1273 | ! ! |
---|
| 1274 | ! WRITE(*,*) 'Init surface masks in T points' |
---|
| 1275 | ! ! |
---|
| 1276 | ! ALLOCATE( Grid%ssmask(nx,ny), Grid%ssumask(nx,ny), Grid%ssvmask(nx,ny) ) |
---|
| 1277 | ! ! |
---|
| 1278 | ! WHERE(Grid%bathy_level(:,:) <= 0. ) |
---|
| 1279 | ! Grid%ssmask(:,:) = 0. |
---|
| 1280 | ! ELSEWHERE |
---|
| 1281 | ! Grid%ssmask(:,:) = 1. |
---|
| 1282 | ! END WHERE |
---|
| 1283 | ! ! |
---|
| 1284 | ! Grid%ssumask(1:nx-1,:) = Grid%ssmask(1:nx-1,:)*Grid%ssmask(2:nx,:) |
---|
| 1285 | ! Grid%ssvmask(:,1:ny-1) = Grid%ssmask(:,1:ny-1)*Grid%ssmask(:,2:ny) |
---|
| 1286 | ! ! |
---|
| 1287 | ! Grid%ssumask(nx,:) = Grid%ssmask(nx,:) |
---|
| 1288 | ! Grid%ssvmask(:,ny) = Grid%ssmask(:,ny) |
---|
| 1289 | ! ! |
---|
| 1290 | ! END SUBROUTINE Init_ssmask |
---|
| 1291 | ! |
---|
| 1292 | !***************************************************** |
---|
[2136] | 1293 | ! subroutine get_mask(name,Grid) |
---|
| 1294 | !***************************************************** |
---|
| 1295 | ! |
---|
| 1296 | SUBROUTINE get_mask(level,posvar,mask,filename) |
---|
| 1297 | ! |
---|
| 1298 | USE io_netcdf |
---|
| 1299 | ! |
---|
| 1300 | CHARACTER(*) filename |
---|
| 1301 | CHARACTER(*) posvar |
---|
| 1302 | INTEGER :: level, nx, ny |
---|
| 1303 | LOGICAL,DIMENSION(:,:),POINTER :: mask |
---|
| 1304 | INTEGER,DIMENSION(:,:),POINTER :: maskT,maskU,maskV |
---|
| 1305 | ! |
---|
| 1306 | TYPE(Coordinates) :: Grid |
---|
| 1307 | ! |
---|
[10025] | 1308 | CALL read_ncdf_var('Bathy_level',filename,Grid%Bathy_level) |
---|
[2136] | 1309 | ! |
---|
| 1310 | nx = SIZE(Grid%Bathy_level,1) |
---|
| 1311 | ny = SIZE(Grid%Bathy_level,2) |
---|
| 1312 | ALLOCATE(maskT(nx,ny),mask(nx,ny)) |
---|
| 1313 | mask = .TRUE. |
---|
| 1314 | ! |
---|
| 1315 | WHERE(Grid%Bathy_level(:,:) <= level-1 ) |
---|
| 1316 | maskT(:,:) = 0 |
---|
| 1317 | ELSEWHERE |
---|
| 1318 | maskT(:,:) = 1 |
---|
| 1319 | END WHERE |
---|
| 1320 | ! |
---|
| 1321 | SELECT CASE(posvar) |
---|
| 1322 | ! |
---|
| 1323 | CASE('T') |
---|
| 1324 | ! |
---|
| 1325 | WHERE(maskT > 0) |
---|
| 1326 | mask = .TRUE. |
---|
| 1327 | ELSEWHERE |
---|
| 1328 | mask = .FALSE. |
---|
| 1329 | END WHERE |
---|
| 1330 | DEALLOCATE(maskT) |
---|
| 1331 | ! |
---|
| 1332 | CASE('U') |
---|
| 1333 | ! |
---|
| 1334 | ALLOCATE(maskU(nx,ny)) |
---|
| 1335 | maskU(1:nx-1,:) = maskT(1:nx-1,:)*maskT(2:nx,:) |
---|
| 1336 | maskU(nx,:) = maskT(nx,:) |
---|
| 1337 | WHERE(maskU > 0) |
---|
| 1338 | mask = .TRUE. |
---|
| 1339 | ELSEWHERE |
---|
| 1340 | mask = .FALSE. |
---|
| 1341 | END WHERE |
---|
| 1342 | DEALLOCATE(maskU,maskT) |
---|
| 1343 | ! |
---|
| 1344 | CASE('V') |
---|
| 1345 | ! |
---|
| 1346 | ALLOCATE(maskV(nx,ny)) |
---|
| 1347 | maskV(:,1:ny-1) = maskT(:,1:ny-1)*maskT(:,2:ny) |
---|
| 1348 | maskV(:,ny) = maskT(:,ny) |
---|
| 1349 | WHERE(maskT > 0) |
---|
| 1350 | mask = .TRUE. |
---|
| 1351 | ELSEWHERE |
---|
| 1352 | mask = .FALSE. |
---|
| 1353 | END WHERE |
---|
| 1354 | DEALLOCATE(maskV,maskT) |
---|
| 1355 | ! |
---|
| 1356 | END SELECT |
---|
| 1357 | ! |
---|
| 1358 | END SUBROUTINE get_mask |
---|
| 1359 | ! |
---|
| 1360 | ! |
---|
| 1361 | !***************************************************** |
---|
| 1362 | ! subroutine read_dimg_var(unit,irec,field) |
---|
| 1363 | !***************************************************** |
---|
| 1364 | ! |
---|
| 1365 | SUBROUTINE read_dimg_var(unit,irec,field,jpk) |
---|
| 1366 | ! |
---|
| 1367 | INTEGER :: unit,irec,jpk |
---|
| 1368 | REAL*8,DIMENSION(:,:,:,:),POINTER :: field |
---|
| 1369 | INTEGER :: k |
---|
| 1370 | ! |
---|
| 1371 | DO k = 1,jpk |
---|
| 1372 | READ(unit,REC=irec) field(:,:,k,1) |
---|
| 1373 | irec = irec + 1 |
---|
| 1374 | ENDDO |
---|
| 1375 | ! |
---|
| 1376 | END SUBROUTINE read_dimg_var |
---|
| 1377 | ! |
---|
| 1378 | ! |
---|
| 1379 | !***************************************************** |
---|
| 1380 | ! subroutine read_dimg_var(unit,irec,field) |
---|
| 1381 | !***************************************************** |
---|
| 1382 | ! |
---|
| 1383 | SUBROUTINE write_dimg_var(unit,irec,field,jpk) |
---|
| 1384 | ! |
---|
| 1385 | INTEGER :: unit,irec,jpk |
---|
| 1386 | REAL*8,DIMENSION(:,:,:,:),POINTER :: field |
---|
| 1387 | INTEGER :: k |
---|
| 1388 | ! |
---|
| 1389 | DO k = 1,jpk |
---|
| 1390 | WRITE(unit,REC=irec) field(:,:,k,1) |
---|
| 1391 | irec = irec + 1 |
---|
| 1392 | ENDDO |
---|
| 1393 | ! |
---|
| 1394 | END SUBROUTINE write_dimg_var |
---|
[10025] | 1395 | |
---|
[2136] | 1396 | END MODULE agrif_readwrite |
---|