[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 | ! |
---|
| 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) |
---|
| 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 | ! |
---|
| 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) |
---|
| 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' /) |
---|
| 148 | CALL Write_Ncdf_dim(dimnames(1),name,nxfin) |
---|
| 149 | CALL Write_Ncdf_dim(dimnames(2),name,nyfin) |
---|
[2136] | 150 | ! |
---|
[2455] | 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 | ! |
---|
[2455] | 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 | ! |
---|
[2455] | 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)) |
---|
[2136] | 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) |
---|
| 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 | ! |
---|
[2455] | 211 | CALL Read_Ncdf_var('mbathy',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 | ! |
---|
[8656] | 234 | status = nf90_create(name,NF90_NOCLOBBER,ncid) |
---|
[2136] | 235 | status = nf90_close(ncid) |
---|
| 236 | ! |
---|
[2455] | 237 | dimnames = (/ 'x','y' /) |
---|
| 238 | CALL Write_Ncdf_dim(dimnames(1),name,nxfin) |
---|
| 239 | CALL Write_Ncdf_dim(dimnames(2),name,nyfin) |
---|
[2136] | 240 | ! |
---|
[2455] | 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') |
---|
| 243 | CALL Write_Ncdf_var('mbathy' ,dimnames,name,Grid%bathy_level,'float') |
---|
[2136] | 244 | ! |
---|
[2455] | 245 | CALL Copy_Ncdf_att('nav_lon',TRIM(parent_meshmask_file),name,MINVAL(Grid%nav_lon),MAXVAL(Grid%nav_lon)) |
---|
| 246 | CALL Copy_Ncdf_att('nav_lat',TRIM(parent_meshmask_file),name,MINVAL(Grid%nav_lat),MAXVAL(Grid%nav_lat)) |
---|
| 247 | CALL Copy_Ncdf_att('mbathy' ,TRIM(parent_meshmask_file),name) |
---|
[2136] | 248 | ! |
---|
| 249 | WRITE(*,*) ' ' |
---|
| 250 | WRITE(*,*) 'Writing bathymetry file: ',name |
---|
| 251 | WRITE(*,*) ' ' |
---|
| 252 | ! |
---|
| 253 | Write_Bathy_level = 1 |
---|
| 254 | ! |
---|
| 255 | END FUNCTION Write_Bathy_level |
---|
| 256 | ! |
---|
| 257 | !***************************************************** |
---|
| 258 | ! function Read_Bathy_meter(name,CoarseGrid,ChildGrid) |
---|
| 259 | !***************************************************** |
---|
| 260 | ! |
---|
| 261 | INTEGER FUNCTION Read_Bathy_meter(name,CoarseGrid,ChildGrid,Pacifique) |
---|
| 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(*,*) '****' |
---|
| 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(*,*) '****' |
---|
| 284 | CALL Read_Ncdf_var('nav_lon',name,CoarseGrid%nav_lon) |
---|
| 285 | CALL Read_Ncdf_var('nav_lat',name,CoarseGrid%nav_lat) |
---|
[2455] | 286 | CALL Read_Ncdf_var(parent_batmet_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 | ! |
---|
| 296 | Read_Bathy_meter = 1 |
---|
| 297 | RETURN |
---|
| 298 | ELSE |
---|
| 299 | WRITE(*,*) '****' |
---|
| 300 | WRITE(*,*) '*** ERROR Bad format for external high resolution database' |
---|
| 301 | WRITE(*,*) '****' |
---|
| 302 | STOP |
---|
| 303 | ENDIF |
---|
| 304 | ! |
---|
| 305 | IF( MAXVAL(ChildGrid%glamt) > 180 ) THEN |
---|
| 306 | ! |
---|
| 307 | WHERE( topo_lon < 0 ) |
---|
| 308 | topo_lon = topo_lon + 360. |
---|
| 309 | END WHERE |
---|
| 310 | ! |
---|
[9149] | 311 | i_min = MAXLOC(topo_lon,mask = topo_lon < MINVAL(ChildGrid%nav_lon)-zdel) |
---|
| 312 | i_max = MINLOC(topo_lon,mask = topo_lon > MAXVAL(ChildGrid%nav_lon)+zdel) |
---|
| 313 | j_min = MAXLOC(topo_lat,mask = topo_lat < MINVAL(ChildGrid%nav_lat)-zdel) |
---|
| 314 | j_max = MINLOC(topo_lat,mask = topo_lat > MAXVAL(ChildGrid%nav_lat)+zdel) |
---|
[2136] | 315 | ! |
---|
| 316 | tabdim1 = ( SIZE(topo_lon) - i_min(1) + 1 ) + i_max(1) |
---|
| 317 | ! |
---|
| 318 | IF(j_min(1)-2 >= 1 .AND. j_max(1)+3 <= SIZE(topo_lat,1) ) THEN |
---|
[2455] | 319 | j_min(1) = j_min(1)-2 |
---|
| 320 | j_max(1) = j_max(1)+3 |
---|
[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 | ! |
---|
[9149] | 354 | i_min = MAXLOC(topo_lon,mask = topo_lon < MINVAL(ChildGrid%nav_lon)-zdel) |
---|
| 355 | i_max = MINLOC(topo_lon,mask = topo_lon > MAXVAL(ChildGrid%nav_lon)+zdel) |
---|
| 356 | j_min = MAXLOC(topo_lat,mask = topo_lat < MINVAL(ChildGrid%nav_lat)-zdel) |
---|
| 357 | j_max = MINLOC(topo_lat,mask = topo_lat > MAXVAL(ChildGrid%nav_lat)+zdel) |
---|
[2136] | 358 | ! |
---|
| 359 | IF(i_min(1)-2 >= 1 .AND. i_max(1)+3 <= SIZE(topo_lon,1) ) THEN |
---|
[2455] | 360 | i_min(1) = i_min(1)-2 |
---|
| 361 | i_max(1) = i_max(1)+3 |
---|
[2136] | 362 | ENDIF |
---|
[2455] | 363 | tabdim1 = i_max(1) - i_min(1) + 1 |
---|
[2136] | 364 | ! |
---|
| 365 | IF(j_min(1)-2 >= 1 .AND. j_max(1)+3 <= SIZE(topo_lat,1) ) THEN |
---|
[2455] | 366 | j_min(1) = j_min(1)-2 |
---|
| 367 | j_max(1) = j_max(1)+3 |
---|
[2136] | 368 | ENDIF |
---|
[2455] | 369 | tabdim2 = j_max(1) - j_min(1) + 1 |
---|
[2136] | 370 | ! |
---|
| 371 | WRITE(*,*) ' ' |
---|
| 372 | WRITE(*,*) 'Reading bathymetry file: ',name |
---|
| 373 | WRITE(*,*) ' ' |
---|
| 374 | ! |
---|
| 375 | ALLOCATE(CoarseGrid%nav_lon(tabdim1,tabdim2)) |
---|
| 376 | ALLOCATE(CoarseGrid%nav_lat(tabdim1,tabdim2)) |
---|
| 377 | ALLOCATE(CoarseGrid%Bathy_meter(tabdim1,tabdim2)) |
---|
| 378 | ! |
---|
[2455] | 379 | DO j = 1,tabdim2 |
---|
| 380 | CoarseGrid%nav_lon(:,j) = topo_lon(i_min(1):i_max(1)) |
---|
[2136] | 381 | END DO |
---|
| 382 | ! |
---|
[2455] | 383 | DO i = 1,tabdim1 |
---|
| 384 | CoarseGrid%nav_lat(i,:) = topo_lat(j_min(1):j_max(1)) |
---|
[2136] | 385 | END DO |
---|
| 386 | ! |
---|
| 387 | status = nf90_open(name,NF90_NOWRITE,ncid) |
---|
[2455] | 388 | status = nf90_inq_varid(ncid,elevation_name,varid) |
---|
[9149] | 389 | IF (status/=nf90_noerr) THEN |
---|
| 390 | WRITE(*,*)"Can't find variable: ", elevation_name |
---|
| 391 | STOP |
---|
| 392 | ENDIF |
---|
| 393 | |
---|
[2455] | 394 | status = nf90_get_var(ncid,varid,CoarseGrid%Bathy_meter, & |
---|
| 395 | & start=(/i_min(1),j_min(1)/),count=(/tabdim1,tabdim2/)) |
---|
[2136] | 396 | ! |
---|
| 397 | ENDIF |
---|
| 398 | ! |
---|
| 399 | status = nf90_close(ncid) |
---|
| 400 | ! |
---|
| 401 | WHERE(CoarseGrid%Bathy_meter.GE.0) |
---|
| 402 | CoarseGrid%Bathy_meter = 0.0 |
---|
| 403 | END WHERE |
---|
| 404 | ! |
---|
[2455] | 405 | CoarseGrid%Bathy_meter(:,:) = -1.0 * CoarseGrid%Bathy_meter(:,:) |
---|
[2136] | 406 | ! |
---|
| 407 | Read_Bathy_meter = 1 |
---|
| 408 | RETURN |
---|
| 409 | ! |
---|
| 410 | END FUNCTION Read_Bathy_meter |
---|
| 411 | ! |
---|
| 412 | ! |
---|
| 413 | !***************************************************** |
---|
| 414 | ! function Read_Bathy_meter(name,CoarseGrid,ChildGrid) |
---|
| 415 | !***************************************************** |
---|
| 416 | ! |
---|
| 417 | INTEGER FUNCTION Read_Bathymeter(name,Grid) |
---|
| 418 | ! |
---|
| 419 | ! |
---|
| 420 | USE io_netcdf |
---|
| 421 | ! |
---|
| 422 | CHARACTER(*) name |
---|
| 423 | TYPE(Coordinates) :: Grid |
---|
| 424 | ! |
---|
[2455] | 425 | CALL Read_Ncdf_var(parent_batmet_name,name,Grid%Bathy_meter) |
---|
[2136] | 426 | ! |
---|
| 427 | WRITE(*,*) ' ' |
---|
| 428 | WRITE(*,*) 'Reading bathymetry file: ',name |
---|
| 429 | WRITE(*,*) ' ' |
---|
| 430 | ! |
---|
| 431 | Read_Bathymeter = 1 |
---|
| 432 | ! |
---|
| 433 | END FUNCTION Read_Bathymeter |
---|
| 434 | ! |
---|
| 435 | !***************************************************** |
---|
| 436 | ! function Write_Bathy_meter(name,Grid) |
---|
| 437 | !***************************************************** |
---|
| 438 | ! |
---|
| 439 | INTEGER FUNCTION Write_Bathy_meter(name,Grid) |
---|
| 440 | ! |
---|
| 441 | USE io_netcdf |
---|
| 442 | ! |
---|
| 443 | CHARACTER(*) name |
---|
| 444 | TYPE(Coordinates) :: Grid |
---|
[2455] | 445 | INTEGER :: status,ncid |
---|
| 446 | CHARACTER(len=1),DIMENSION(2) :: dimnames |
---|
[2136] | 447 | INTEGER :: nx,ny |
---|
| 448 | ! |
---|
| 449 | status = nf90_create(name,NF90_WRITE,ncid) |
---|
| 450 | status = nf90_close(ncid) |
---|
| 451 | ! |
---|
| 452 | nx = SIZE(Grid%bathy_meter,1) |
---|
| 453 | ny = SIZE(Grid%bathy_meter,2) |
---|
[2455] | 454 | dimnames = (/ 'x','y' /) |
---|
[2136] | 455 | |
---|
[2455] | 456 | CALL Write_Ncdf_dim(dimnames(1),name,nx) |
---|
| 457 | CALL Write_Ncdf_dim(dimnames(2),name,ny) |
---|
[2136] | 458 | ! |
---|
[2455] | 459 | CALL Write_Ncdf_var('nav_lon' ,dimnames,name,Grid%nav_lon ,'float') |
---|
| 460 | CALL Write_Ncdf_var('nav_lat' ,dimnames,name,Grid%nav_lat ,'float') |
---|
| 461 | CALL Write_Ncdf_var(parent_batmet_name,dimnames,name,Grid%bathy_meter,'float') |
---|
[9628] | 462 | CALL Write_Ncdf_var('weight' ,dimnames,name,Grid%wgt ,'float') |
---|
[2136] | 463 | ! |
---|
[2455] | 464 | CALL Copy_Ncdf_att('nav_lon' ,TRIM(parent_bathy_meter),name,MINVAL(Grid%nav_lon),MAXVAL(Grid%nav_lon)) |
---|
| 465 | CALL Copy_Ncdf_att('nav_lat' ,TRIM(parent_bathy_meter),name,MINVAL(Grid%nav_lat),MAXVAL(Grid%nav_lat)) |
---|
| 466 | CALL Copy_Ncdf_att(parent_batmet_name,TRIM(parent_bathy_meter),name) |
---|
[2136] | 467 | ! |
---|
| 468 | WRITE(*,*) ' ' |
---|
| 469 | WRITE(*,*) 'Writing bathymetry file: ',name |
---|
| 470 | WRITE(*,*) ' ' |
---|
| 471 | ! |
---|
| 472 | Write_Bathy_meter = 1 |
---|
| 473 | ! |
---|
| 474 | END FUNCTION Write_Bathy_meter |
---|
| 475 | ! |
---|
| 476 | !***************************************************** |
---|
| 477 | ! function set_child_name(Parentname,Childname) |
---|
| 478 | !***************************************************** |
---|
| 479 | ! |
---|
| 480 | SUBROUTINE set_child_name(Parentname,Childname) |
---|
| 481 | ! |
---|
| 482 | CHARACTER(*),INTENT(in) :: Parentname |
---|
| 483 | CHARACTER(*),INTENT(out) :: Childname |
---|
| 484 | CHARACTER(2) :: prefix |
---|
| 485 | INTEGER :: pos |
---|
| 486 | ! |
---|
| 487 | pos = INDEX(TRIM(Parentname),'/',back=.TRUE.) |
---|
| 488 | ! |
---|
| 489 | prefix=Parentname(pos+1:pos+2) |
---|
| 490 | IF (prefix == '1_') THEN |
---|
| 491 | Childname = '2_'//Parentname(pos+3:LEN(Parentname)) |
---|
| 492 | ELSEIF (prefix == '2_') THEN |
---|
| 493 | Childname = '3_'//Parentname(pos+3:LEN(Parentname)) |
---|
| 494 | ELSEIF (prefix == '3_') THEN |
---|
| 495 | Childname = '4_'//Parentname(pos+3:LEN(Parentname)) |
---|
| 496 | ELSEIF (prefix == '4_') THEN |
---|
| 497 | Childname = '5_'//Parentname(pos+3:LEN(Parentname)) |
---|
| 498 | ELSE |
---|
| 499 | Childname = '1_'//Parentname(pos+1:LEN(Parentname)) |
---|
| 500 | ENDIF |
---|
| 501 | ! |
---|
| 502 | END SUBROUTINE set_child_name |
---|
| 503 | ! |
---|
| 504 | !***************************************************** |
---|
| 505 | ! function set_child_name(Parentname,Childname) |
---|
| 506 | !***************************************************** |
---|
| 507 | ! |
---|
| 508 | !***************************************************** |
---|
| 509 | ! subroutine get_interptype(varname,interp_type,conservation) |
---|
| 510 | !***************************************************** |
---|
| 511 | ! |
---|
| 512 | SUBROUTINE get_interptype( varname,interp_type,conservation ) |
---|
| 513 | ! |
---|
| 514 | LOGICAL,OPTIONAL :: conservation |
---|
| 515 | CHARACTER(*) :: interp_type,varname |
---|
| 516 | INTEGER :: pos,pos1,pos2,pos3,i,k |
---|
| 517 | LOGICAL :: find |
---|
| 518 | i=1 |
---|
| 519 | DO WHILE ( TRIM(VAR_INTERP(i)) .NE. 'NULL' ) |
---|
| 520 | pos = INDEX( TRIM(VAR_INTERP(i)) , TRIM(varname) ) |
---|
| 521 | IF ( pos .NE. 0 ) THEN |
---|
| 522 | pos1 = INDEX( TRIM(VAR_INTERP(i)) , 'bicubic' ) |
---|
| 523 | pos2 = INDEX( TRIM(VAR_INTERP(i)) , 'bilinear' ) |
---|
| 524 | pos3 = INDEX( TRIM(VAR_INTERP(i)) , 'conservative' ) |
---|
| 525 | ! initialize interp_type |
---|
| 526 | IF( pos1 .NE. 0 ) interp_type = 'bicubic' |
---|
| 527 | IF( pos2 .NE. 0 ) interp_type = 'bilinear' |
---|
| 528 | IF( pos1 .EQ. 0 .AND. pos2 .EQ. 0) interp_type = 'bicubic' |
---|
| 529 | ! initialize conservation |
---|
| 530 | IF( pos3 .NE. 0 .AND. PRESENT(conservation) ) THEN |
---|
| 531 | conservation = .TRUE. |
---|
| 532 | RETURN |
---|
| 533 | ELSE |
---|
| 534 | conservation = .FALSE. |
---|
| 535 | ENDIF |
---|
| 536 | find = .FALSE. |
---|
| 537 | IF( PRESENT(conservation) ) THEN |
---|
| 538 | k=0 |
---|
| 539 | conservation = .FALSE. |
---|
[2455] | 540 | DO WHILE( k < SIZE(flxtab) .AND. .NOT.find ) |
---|
[2136] | 541 | k = k+1 |
---|
| 542 | IF( TRIM(varname) .EQ. TRIM(flxtab(k)) ) THEN |
---|
| 543 | conservation = .TRUE. |
---|
| 544 | find = .TRUE. |
---|
| 545 | ENDIF |
---|
| 546 | END DO |
---|
| 547 | ENDIF |
---|
| 548 | RETURN |
---|
| 549 | ENDIF |
---|
| 550 | i = i+1 |
---|
| 551 | END DO |
---|
| 552 | !default values interp_type = bicubic // conservation = false |
---|
| 553 | interp_type = 'bicubic' |
---|
| 554 | IF( PRESENT(conservation) ) conservation = .FALSE. |
---|
| 555 | |
---|
| 556 | RETURN |
---|
| 557 | ! |
---|
| 558 | END SUBROUTINE get_interptype |
---|
| 559 | ! |
---|
| 560 | !***************************************************** |
---|
| 561 | ! end subroutine get_interptype |
---|
| 562 | !***************************************************** |
---|
| 563 | ! |
---|
| 564 | !***************************************************** |
---|
| 565 | ! subroutine Init_mask(name,Grid) |
---|
| 566 | !***************************************************** |
---|
| 567 | ! |
---|
| 568 | SUBROUTINE Init_mask(name,Grid,jpiglo,jpjglo) |
---|
| 569 | ! |
---|
| 570 | USE io_netcdf |
---|
| 571 | ! |
---|
| 572 | CHARACTER(*) name |
---|
| 573 | INTEGER :: nx,ny,k,i,j,jpiglo,jpjglo |
---|
| 574 | TYPE(Coordinates) :: Grid |
---|
| 575 | REAL*8, POINTER, DIMENSION(:,:) ::zwf => NULL() |
---|
| 576 | ! |
---|
| 577 | IF(jpiglo == 1 .AND. jpjglo == 1) THEN |
---|
| 578 | CALL Read_Ncdf_var('Bathy_level',name,Grid%Bathy_level) |
---|
| 579 | ELSE |
---|
| 580 | CALL Read_Ncdf_var('Bathy_level',name,Grid%Bathy_level,(/jpizoom,jpjzoom/),(/jpiglo,jpjglo/) ) |
---|
| 581 | ENDIF |
---|
| 582 | |
---|
| 583 | |
---|
| 584 | ! |
---|
| 585 | WRITE(*,*) 'Init masks in T,U,V,F points' |
---|
| 586 | ! |
---|
| 587 | nx = SIZE(Grid%Bathy_level,1) |
---|
| 588 | ny = SIZE(Grid%Bathy_level,2) |
---|
| 589 | ! |
---|
| 590 | ! |
---|
| 591 | ALLOCATE(Grid%tmask(nx,ny,N), & |
---|
| 592 | Grid%umask(nx,ny,N), & |
---|
| 593 | Grid%vmask(nx,ny,N), & |
---|
| 594 | Grid%fmask(nx,ny,N)) |
---|
| 595 | ! |
---|
| 596 | DO k = 1,N |
---|
| 597 | ! |
---|
| 598 | WHERE(Grid%Bathy_level(:,:) <= k-1 ) |
---|
| 599 | Grid%tmask(:,:,k) = 0 |
---|
| 600 | ELSEWHERE |
---|
| 601 | Grid%tmask(:,:,k) = 1 |
---|
| 602 | END WHERE |
---|
| 603 | ! |
---|
| 604 | END DO |
---|
| 605 | ! |
---|
| 606 | Grid%umask(1:nx-1,:,:) = Grid%tmask(1:nx-1,:,:)*Grid%tmask(2:nx,:,:) |
---|
| 607 | Grid%vmask(:,1:ny-1,:) = Grid%tmask(:,1:ny-1,:)*Grid%tmask(:,2:ny,:) |
---|
| 608 | ! |
---|
| 609 | Grid%umask(nx,:,:) = Grid%tmask(nx,:,:) |
---|
| 610 | Grid%vmask(:,ny,:) = Grid%tmask(:,ny,:) |
---|
| 611 | ! |
---|
| 612 | Grid%fmask(1:nx-1,1:ny-1,:) = Grid%tmask(1:nx-1,1:ny-1,:)*Grid%tmask(2:nx,1:ny-1,:)* & |
---|
| 613 | Grid%tmask(1:nx-1,2:ny,:)*Grid%tmask(2:nx,2:ny,:) |
---|
| 614 | ! |
---|
| 615 | Grid%fmask(nx,:,:) = Grid%tmask(nx,:,:) |
---|
| 616 | Grid%fmask(:,ny,:) = Grid%tmask(:,ny,:) |
---|
| 617 | ! |
---|
| 618 | ALLOCATE(zwf(nx,ny)) |
---|
| 619 | ! |
---|
| 620 | DO k = 1,N |
---|
| 621 | ! |
---|
| 622 | zwf(:,:) = Grid%fmask(:,:,k) |
---|
| 623 | ! |
---|
| 624 | DO j = 2, ny-1 |
---|
| 625 | DO i = 2,nx-1 |
---|
| 626 | IF( Grid%fmask(i,j,k) == 0. ) THEN |
---|
| 627 | 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))) |
---|
| 628 | END IF |
---|
| 629 | END DO |
---|
| 630 | END DO |
---|
| 631 | ! |
---|
| 632 | DO j = 2, ny-1 |
---|
| 633 | IF( Grid%fmask(1,j,k) == 0. ) THEN |
---|
| 634 | Grid%fmask(1,j,k) = shlat * MIN(1.,MAX(zwf(2,j),zwf(1,j+1),zwf(1,j-1))) |
---|
| 635 | ENDIF |
---|
| 636 | |
---|
| 637 | IF( Grid%fmask(nx,j,k) == 0. ) THEN |
---|
| 638 | Grid%fmask(nx,j,k) = shlat * MIN(1.,MAX(zwf(nx,j+1),zwf(nx-1,j),zwf(nx,j-1))) |
---|
| 639 | ENDIF |
---|
| 640 | END DO |
---|
| 641 | ! |
---|
| 642 | DO i = 2, nx-1 |
---|
| 643 | IF( Grid%fmask(i,1,k) == 0. ) THEN |
---|
| 644 | Grid%fmask(i, 1 ,k) = shlat*MIN( 1.,MAX(zwf(i+1,1),zwf(i,2),zwf(i-1,1))) |
---|
| 645 | ENDIF |
---|
| 646 | ! |
---|
| 647 | IF( Grid%fmask(i,ny,k) == 0. ) THEN |
---|
| 648 | Grid%fmask(i,ny,k) = shlat * MIN(1.,MAX(zwf(i+1,ny),zwf(i-1,ny),zwf(i,ny-1))) |
---|
| 649 | ENDIF |
---|
| 650 | END DO |
---|
| 651 | !! |
---|
| 652 | END DO |
---|
| 653 | !! |
---|
| 654 | END SUBROUTINE Init_mask |
---|
| 655 | ! |
---|
| 656 | !***************************************************** |
---|
| 657 | ! end subroutine Init_mask |
---|
| 658 | !***************************************************** |
---|
| 659 | ! |
---|
| 660 | !***************************************************** |
---|
| 661 | ! subroutine Init_Tmask(name,Grid) |
---|
| 662 | !***************************************************** |
---|
| 663 | ! |
---|
| 664 | SUBROUTINE Init_Tmask(name,Grid,jpiglo,jpjglo) |
---|
| 665 | ! |
---|
| 666 | USE io_netcdf |
---|
| 667 | ! |
---|
| 668 | CHARACTER(*) name |
---|
| 669 | INTEGER :: nx,ny,k,i,j,jpiglo,jpjglo |
---|
| 670 | TYPE(Coordinates) :: Grid |
---|
| 671 | REAL*8, POINTER, DIMENSION(:,:) ::zwf => NULL() |
---|
| 672 | ! |
---|
| 673 | IF(jpiglo == 1 .AND. jpjglo == 1) THEN |
---|
| 674 | CALL Read_Ncdf_var('Bathy_level',name,Grid%Bathy_level) |
---|
| 675 | ELSE |
---|
| 676 | CALL Read_Ncdf_var('Bathy_level',name,Grid%Bathy_level,(/jpizoom,jpjzoom/),(/jpiglo,jpjglo/) ) |
---|
| 677 | ENDIF |
---|
| 678 | ! |
---|
| 679 | nx = SIZE(Grid%Bathy_level,1) |
---|
| 680 | ny = SIZE(Grid%Bathy_level,2) |
---|
| 681 | ! |
---|
| 682 | WRITE(*,*) 'Init masks in T points' |
---|
| 683 | ! |
---|
| 684 | ALLOCATE(Grid%tmask(nx,ny,N)) |
---|
| 685 | ! |
---|
| 686 | DO k = 1,N |
---|
| 687 | ! |
---|
| 688 | WHERE(Grid%Bathy_level(:,:) <= k-1 ) |
---|
| 689 | Grid%tmask(:,:,k) = 0. |
---|
| 690 | ELSEWHERE |
---|
| 691 | Grid%tmask(:,:,k) = 1. |
---|
| 692 | END WHERE |
---|
| 693 | ! |
---|
| 694 | END DO |
---|
| 695 | ! |
---|
| 696 | END SUBROUTINE Init_Tmask |
---|
| 697 | ! |
---|
| 698 | !***************************************************** |
---|
| 699 | ! subroutine get_mask(name,Grid) |
---|
| 700 | !***************************************************** |
---|
| 701 | ! |
---|
| 702 | SUBROUTINE get_mask(level,posvar,mask,filename) |
---|
| 703 | ! |
---|
| 704 | USE io_netcdf |
---|
| 705 | ! |
---|
| 706 | CHARACTER(*) filename |
---|
| 707 | CHARACTER(*) posvar |
---|
| 708 | INTEGER :: level, nx, ny |
---|
| 709 | LOGICAL,DIMENSION(:,:),POINTER :: mask |
---|
| 710 | INTEGER,DIMENSION(:,:),POINTER :: maskT,maskU,maskV |
---|
| 711 | ! |
---|
| 712 | TYPE(Coordinates) :: Grid |
---|
| 713 | ! |
---|
| 714 | CALL Read_Ncdf_var('Bathy_level',filename,Grid%Bathy_level) |
---|
| 715 | ! |
---|
| 716 | nx = SIZE(Grid%Bathy_level,1) |
---|
| 717 | ny = SIZE(Grid%Bathy_level,2) |
---|
| 718 | ALLOCATE(maskT(nx,ny),mask(nx,ny)) |
---|
| 719 | mask = .TRUE. |
---|
| 720 | ! |
---|
| 721 | WHERE(Grid%Bathy_level(:,:) <= level-1 ) |
---|
| 722 | maskT(:,:) = 0 |
---|
| 723 | ELSEWHERE |
---|
| 724 | maskT(:,:) = 1 |
---|
| 725 | END WHERE |
---|
| 726 | ! |
---|
| 727 | SELECT CASE(posvar) |
---|
| 728 | ! |
---|
| 729 | CASE('T') |
---|
| 730 | ! |
---|
| 731 | WHERE(maskT > 0) |
---|
| 732 | mask = .TRUE. |
---|
| 733 | ELSEWHERE |
---|
| 734 | mask = .FALSE. |
---|
| 735 | END WHERE |
---|
| 736 | DEALLOCATE(maskT) |
---|
| 737 | ! |
---|
| 738 | CASE('U') |
---|
| 739 | ! |
---|
| 740 | ALLOCATE(maskU(nx,ny)) |
---|
| 741 | maskU(1:nx-1,:) = maskT(1:nx-1,:)*maskT(2:nx,:) |
---|
| 742 | maskU(nx,:) = maskT(nx,:) |
---|
| 743 | WHERE(maskU > 0) |
---|
| 744 | mask = .TRUE. |
---|
| 745 | ELSEWHERE |
---|
| 746 | mask = .FALSE. |
---|
| 747 | END WHERE |
---|
| 748 | DEALLOCATE(maskU,maskT) |
---|
| 749 | ! |
---|
| 750 | CASE('V') |
---|
| 751 | ! |
---|
| 752 | ALLOCATE(maskV(nx,ny)) |
---|
| 753 | maskV(:,1:ny-1) = maskT(:,1:ny-1)*maskT(:,2:ny) |
---|
| 754 | maskV(:,ny) = maskT(:,ny) |
---|
| 755 | WHERE(maskT > 0) |
---|
| 756 | mask = .TRUE. |
---|
| 757 | ELSEWHERE |
---|
| 758 | mask = .FALSE. |
---|
| 759 | END WHERE |
---|
| 760 | DEALLOCATE(maskV,maskT) |
---|
| 761 | ! |
---|
| 762 | END SELECT |
---|
| 763 | ! |
---|
| 764 | END SUBROUTINE get_mask |
---|
| 765 | ! |
---|
| 766 | !***************************************************** |
---|
| 767 | ! end subroutine get_mask |
---|
| 768 | !***************************************************** |
---|
| 769 | ! |
---|
| 770 | ! |
---|
| 771 | !***************************************************** |
---|
| 772 | ! subroutine read_dimg_var(unit,irec,field) |
---|
| 773 | !***************************************************** |
---|
| 774 | ! |
---|
| 775 | SUBROUTINE read_dimg_var(unit,irec,field,jpk) |
---|
| 776 | ! |
---|
| 777 | INTEGER :: unit,irec,jpk |
---|
| 778 | REAL*8,DIMENSION(:,:,:,:),POINTER :: field |
---|
| 779 | INTEGER :: k |
---|
| 780 | ! |
---|
| 781 | DO k = 1,jpk |
---|
| 782 | READ(unit,REC=irec) field(:,:,k,1) |
---|
| 783 | irec = irec + 1 |
---|
| 784 | ENDDO |
---|
| 785 | ! |
---|
| 786 | END SUBROUTINE read_dimg_var |
---|
| 787 | ! |
---|
| 788 | ! |
---|
| 789 | !***************************************************** |
---|
| 790 | ! subroutine read_dimg_var(unit,irec,field) |
---|
| 791 | !***************************************************** |
---|
| 792 | ! |
---|
| 793 | SUBROUTINE write_dimg_var(unit,irec,field,jpk) |
---|
| 794 | ! |
---|
| 795 | INTEGER :: unit,irec,jpk |
---|
| 796 | REAL*8,DIMENSION(:,:,:,:),POINTER :: field |
---|
| 797 | INTEGER :: k |
---|
| 798 | ! |
---|
| 799 | DO k = 1,jpk |
---|
| 800 | WRITE(unit,REC=irec) field(:,:,k,1) |
---|
| 801 | irec = irec + 1 |
---|
| 802 | ENDDO |
---|
| 803 | ! |
---|
| 804 | END SUBROUTINE write_dimg_var |
---|
| 805 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
| 806 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
| 807 | END MODULE agrif_readwrite |
---|