![(please configure the [header_logo] section in trac.ini)](/grisli/chrome/site/Laurentide-16km-3D-Umag-water-bordure-GRISLI-small.png)
- Timestamp:
- 12/14/23 15:27:11 (6 months ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/SOURCES/Ant40_files/output_anta40_mod-0.4.f90
r275 r457 13 13 module output_antarcti_mod 14 14 15 USE module3D_phy 16 use bilan_eau_mod 17 use netcdf 18 use io_netcdf_grisli 19 use bilan_flux_mod 20 21 implicit none 22 23 !real :: vol ; !integer :: np 24 integer :: nflot !< nbr de point flottant 25 real :: bmean !< 26 real :: accmean !< accumulation moyenne 27 real :: ablmean !< ablation moyenne 28 real :: calvmean !< moyenne calving 29 real :: ablbordmean !< 30 real :: bmeltmean !< moyenne fusion basale 31 real :: tbmean !< temperature basale moyenne 32 real :: tbdotmean !< moyenne variation / temps temperature basale 33 real :: vsmean !< vitesse de surface moyenne 34 !real :: vsdotmean !< moyenne variation / temps vitesse de surface 35 !real :: uzsmean !!!! utilise ? !< vitesse verticale de surface moyenne 36 real :: uzsdotmean !< moyenne variation / temps vitesse verticale de surface 37 real :: uzkmean !< moyenne vitesse verticale de surface 38 real :: hdotmean !< moyenne derivee / temps de h 39 real :: bdotmean !< moyenne bedrock derive / temps 40 real :: volf !< volume au dessus de la flottaison 41 42 real :: lim !< Total ice mass 43 real :: iareag !< surface posee 44 real :: iareaf !< surface flottante 45 real :: tendacabf !< Total SMB flux 46 real :: tendlibmassbf !< Total Basal mass balance flux 47 real :: tendlibmassbffl !< Total Basal mass balance flux on floating ice 48 real :: tendlicalvf !< Total calving flux 49 real :: tendlifmassbf !< Total calving and ice front melting flux 50 real :: tendligroundf !< Total grounding line flux 51 52 real,dimension(nx,ny) :: corrsurf !< facteur de correction de la surface 53 !real,parameter :: ice_density=910. !< densite de la glace pour conversion en masse 54 55 ! variables netcdf 56 integer,parameter :: ncshortout=1 ! 1 sorties netcdf short initMIP 57 integer,parameter :: nvar=10 ! nombre de variables dans le fichier de sortie temporel netcdf 58 integer :: ncid 59 integer :: status 60 integer :: timeDimID 61 integer :: timeVarID 62 integer,dimension(nvar) :: varID 63 integer :: nbtimeout ! index time output 64 65 real,dimension(nvar) :: var_shortoutput 66 15 USE module3D_phy 16 use bilan_eau_mod 17 use netcdf 18 use io_netcdf_grisli 19 use bilan_flux_mod 20 21 implicit none 22 23 !real :: vol ; !integer :: np 24 integer :: nflot !< nbr de point flottant 25 real :: bmean !< 26 real :: accmean !< accumulation moyenne 27 real :: ablmean !< ablation moyenne 28 real :: calvmean !< moyenne calving 29 real :: ablbordmean !< 30 real :: bmeltmean !< moyenne fusion basale 31 real :: tbmean !< temperature basale moyenne 32 real :: tbdotmean !< moyenne variation / temps temperature basale 33 real :: vsmean !< vitesse de surface moyenne 34 !real :: vsdotmean !< moyenne variation / temps vitesse de surface 35 !real :: uzsmean !!!! utilise ? !< vitesse verticale de surface moyenne 36 real :: uzsdotmean !< moyenne variation / temps vitesse verticale de surface 37 real :: uzkmean !< moyenne vitesse verticale de surface 38 real :: hdotmean !< moyenne derivee / temps de h 39 real :: bdotmean !< moyenne bedrock derive / temps 40 real :: volf !< volume au dessus de la flottaison 41 42 real :: lim !< Total ice mass 43 real :: iareag !< surface posee 44 real :: iareaf !< surface flottante 45 real :: tendacabf !< Total SMB flux 46 real :: tendlibmassbf !< Total Basal mass balance flux 47 real :: tendlibmassbffl !< Total Basal mass balance flux on floating ice 48 real :: tendlicalvf !< Total calving flux 49 real :: tendlifmassbf !< Total calving and ice front melting flux 50 real :: tendligroundf !< Total grounding line flux 51 52 real,dimension(nx,ny) :: corrsurf !< facteur de correction de la surface 53 !real,parameter :: ice_density=910. !< densite de la glace pour conversion en masse 54 55 ! variables netcdf 56 integer,parameter :: ncshortout=1 ! 1 sorties netcdf short initMIP 57 integer,parameter :: nvar=10 ! nombre de variables dans le fichier de sortie temporel netcdf 58 integer :: ncid 59 integer :: status 60 integer :: timeDimID 61 integer :: timeVarID 62 integer,dimension(nvar) :: varID 63 integer :: nbtimeout ! index time output 64 real,dimension(nvar) :: var_shortoutput 65 66 ! variables netcdf icecore (ic) 67 integer,parameter :: ncicecoreout=1 ! 1 sorties netcdf icecore / 0 no output 68 integer :: nbic ! number of ice core site in output 69 integer,parameter :: nvaric=10 ! nombre de variables dans le fichier ice core 70 real,parameter :: dticecoreout=50. ! sorties ice core tous les dticecoreout 71 integer :: ncidic 72 integer :: timeDimIDic 73 integer :: timeVarIDic 74 integer :: stationDimIDic 75 integer :: stationVarIDic 76 integer :: station_strlenDimIDic 77 integer :: lonVarIDic 78 integer :: latVarIDic 79 integer,dimension(nvaric) :: varIDic 80 integer :: nbtimeoutic ! index time output 81 real,allocatable,dimension(:,:) :: var_icoutput 82 character(len=100), allocatable, dimension(:) :: site_nameic 83 integer, allocatable, dimension(:) :: iic,jic ! coordonnees i,j des sites ice core 67 84 68 85 CONTAINS 69 86 70 subroutine init_outshort 71 72 double precision,dimension(:,:),pointer :: tab !< tableau 2d real pointer 73 character(len=100),dimension(nvar) :: namevar ! name, standard_name, long_name, units 74 character(len=100),dimension(nvar) :: standard_name 75 character(len=100),dimension(nvar) :: long_name 76 character(len=100),dimension(nvar) :: units 77 78 !ndisp sorite courte tous les ndisp 79 NDISP=100 80 81 82 if (ncshortout.eq.1) then ! ecriture netcdf 83 84 ! lecture du fichier avec les corrections de surface 85 if (geoplace.eq.'ant16km') then 86 call Read_Ncdf_var('z',trim(DIRNAMEINP)//'/corrsurf-initMIP-16km.grd',tab) 87 corrsurf(:,:)=tab(:,:) 88 else 89 corrsurf(:,:)= 1. 90 end if 91 92 open(568,file=trim(dirsource)//'/Fichiers-parametres/short-ISMIPnc.dat',status='old') 93 ! lecture en-tete 94 read(568,*) 95 read(568,*) 96 read(568,*) 97 read(568,*) 98 read(568,*) 99 ! lecture des infos sur les variables : 100 do k=1,nvar 101 read(568,'(a100)') namevar(k) 102 read(568,'(a100)') standard_name(k) 103 read(568,'(a100)') long_name(k) 104 read(568,'(a100)') units(k) 105 read(568,*) 106 enddo 107 close(568) 108 ! Fichier Netcdf initMIP 109 ! creation du fichier Netcdf : 110 status=nf90_create(path = trim(dirnameout)//'short'//runname//'.nc', cmode = nf90_clobber, ncid = ncid) 111 if (status /= nf90_noerr) call handle_err(status) 112 113 ! definition des dimension : 114 status=nf90_def_dim(ncid, name="time", len=NF90_UNLIMITED, dimid=timeDimID) 115 if (status /= nf90_noerr) call handle_err(status) 116 status=nf90_def_var(ncid, name="time", xtype=nf90_float, dimids=(/ timeDimID/), varid=timeVarID) 117 if (status /= nf90_noerr) call handle_err(status) 118 status=nf90_put_att(ncid, timeVarID, "standard_name", "time") 119 if (status /= nf90_noerr) call handle_err(status) 120 status=nf90_put_att(ncid, timeVarID,"units", "years since 1995-01-01 00:00:00") 121 if (status /= nf90_noerr) call handle_err(status) 122 status=nf90_put_att(ncid, timeVarID,"calendar", "360_day") 123 if (status /= nf90_noerr) call handle_err(status) 124 125 ! definition des variables de sortie : 126 do k=1,nvar ! boucle sur le nbr de variable a definir 127 status=nf90_def_var(ncid, name=trim(namevar(k)), xtype=nf90_float, dimids= & 128 (/ timeDimID /), varid=varID(k)) 129 if (status /= nf90_noerr) call handle_err(status) 130 status=nf90_put_att(ncid, varID(k), "standard_name", trim(standard_name(k))) 131 if (status /= nf90_noerr) call handle_err(status) 132 status=nf90_put_att(ncid, varID(k), "long_name", trim(long_name(k))) 133 if (status /= nf90_noerr) call handle_err(status) 134 status=nf90_put_att(ncid, varID(k), "units", trim(units(k))) 135 if (status /= nf90_noerr) call handle_err(status) 136 enddo 137 138 ! fin de la definition du fchier : 139 status=nf90_enddef(ncid) 140 if (status /= nf90_noerr) call handle_err(status) 141 nbtimeout = 0 ! initialisation compteur sorties axe time 142 else ! pas de sortie netcdf et sans correction de surface 143 corrsurf(:,:)=1. 144 endif 145 146 end subroutine init_outshort 147 148 149 150 !_________________________________________________________________________ 151 subroutine shortoutput 152 153 ! 1_initialization 154 !------------------ 155 real :: smax 156 157 vol=0. 158 np=0 159 nflot=0 160 hmax=0. 161 smax=0. 162 bmean=0. 163 accmean=0. 164 ablmean=0. 165 calvmean=0. 166 ablbordmean=0. 167 bmeltmean=0. 168 tbmean=0. 169 tbdotmean=0. 170 vsmean=0. 171 ! vsdotmean=0. 172 ! uzsmean=0. 173 uzsdotmean=0. 174 uzkmean=0. 175 hdotmean=0. 176 bdotmean=0. 177 volf=0. 178 lim=0. 179 tendacabf=0. 180 iareag=0. 181 iareaf=0. 182 tendlicalvf=0. 183 tendligroundf=0. 184 tendlibmassbf=0. 185 tendlibmassbffl=0. 186 tendlifmassbf=0. 187 ! 2_preparing outputs 188 !-------------------- 87 subroutine init_outshort 88 89 double precision,dimension(:,:),pointer :: tab !< tableau 2d real pointer 90 character(len=100),dimension(nvar) :: namevar ! name, standard_name, long_name, units 91 character(len=100),dimension(nvar) :: standard_name 92 character(len=100),dimension(nvar) :: long_name 93 character(len=100),dimension(nvar) :: units 94 95 ! ice core variables 96 character(len=100),allocatable,dimension(:) :: ic_name ! ice core name 97 real,allocatable,dimension(:) :: ic_lon ! ice core lon 98 real,allocatable,dimension(:) :: ic_lat ! ice core lat 99 character(len=100),dimension(nvaric) :: namevaric ! name, standard_name, long_name, units 100 character(len=100),dimension(nvaric) :: standard_nameic 101 character(len=100),dimension(nvaric) :: long_nameic 102 character(len=100),dimension(nvaric) :: unitsic 103 integer :: err ! recuperation erreur 104 integer :: ii,jj 105 integer :: nloop ! nbr de tour dans la boucle while (evite boucle folle) 106 real :: dist_lonlat, dist_lonlat_new, distance ! distance to lon lat ice core site 107 108 109 110 !ndisp sorite courte tous les ndisp 111 NDISP=100 112 113 114 if (ncshortout.eq.1) then ! ecriture netcdf 115 116 ! lecture du fichier avec les corrections de surface 117 if (geoplace.eq.'ant16km') then 118 call Read_Ncdf_var('z',trim(DIRNAMEINP)//'/corrsurf-initMIP-16km.grd',tab) 119 corrsurf(:,:)=tab(:,:) 120 else 121 corrsurf(:,:)= 1. 122 end if 123 124 open(568,file=trim(dirsource)//'/Fichiers-parametres/short-ISMIPnc.dat',status='old') 125 ! lecture en-tete 126 read(568,*) 127 read(568,*) 128 read(568,*) 129 read(568,*) 130 read(568,*) 131 ! lecture des infos sur les variables : 132 do k=1,nvar 133 read(568,'(a100)') namevar(k) 134 read(568,'(a100)') standard_name(k) 135 read(568,'(a100)') long_name(k) 136 read(568,'(a100)') units(k) 137 read(568,*) 138 enddo 139 close(568) 140 141 ! Fichier Netcdf initMIP 142 ! creation du fichier Netcdf : 143 status=nf90_create(path = trim(dirnameout)//'short'//runname//'.nc', cmode = nf90_clobber, ncid = ncid) 144 if (status /= nf90_noerr) call handle_err(status) 145 146 ! definition des dimension : 147 status=nf90_def_dim(ncid, name="time", len=NF90_UNLIMITED, dimid=timeDimID) 148 if (status /= nf90_noerr) call handle_err(status) 149 status=nf90_def_var(ncid, name="time", xtype=nf90_float, dimids=(/ timeDimID/), varid=timeVarID) 150 if (status /= nf90_noerr) call handle_err(status) 151 status=nf90_put_att(ncid, timeVarID, "standard_name", "time") 152 if (status /= nf90_noerr) call handle_err(status) 153 status=nf90_put_att(ncid, timeVarID,"units", "years since 1995-01-01 00:00:00") 154 if (status /= nf90_noerr) call handle_err(status) 155 status=nf90_put_att(ncid, timeVarID,"calendar", "360_day") 156 if (status /= nf90_noerr) call handle_err(status) 157 158 ! definition des variables de sortie : 159 do k=1,nvar ! boucle sur le nbr de variable a definir 160 status=nf90_def_var(ncid, name=trim(namevar(k)), xtype=nf90_float, dimids= & 161 (/ timeDimID /), varid=varID(k)) 162 if (status /= nf90_noerr) call handle_err(status) 163 status=nf90_put_att(ncid, varID(k), "standard_name", trim(standard_name(k))) 164 if (status /= nf90_noerr) call handle_err(status) 165 status=nf90_put_att(ncid, varID(k), "long_name", trim(long_name(k))) 166 if (status /= nf90_noerr) call handle_err(status) 167 status=nf90_put_att(ncid, varID(k), "units", trim(units(k))) 168 if (status /= nf90_noerr) call handle_err(status) 169 enddo 170 171 ! fin de la definition du fichier : 172 status=nf90_enddef(ncid) 173 if (status /= nf90_noerr) call handle_err(status) 174 nbtimeout = 0 ! initialisation compteur sorties axe time 175 else ! pas de sortie netcdf et sans correction de surface 176 corrsurf(:,:)=1. 177 endif 178 179 180 if (ncicecoreout.eq.1) then ! ecriture netcdf Ice Core 181 ! Fichier Netcdf ice core Antarctique 182 ! Lecture du fichier avec la liste des sites de forages 183 ! infos sur les variables : 184 open(568,file=trim(dirsource)//'/Fichiers-parametres/icecore-Ant.dat',status='old') 185 ! lecture en-tete 186 read(568,*) 187 read(568,*) 188 read(568,*) 189 read(568,*) 190 read(568,*) 191 ! lectuire du nombre de sites de forages: 192 read(568,*) nbic 193 read(568,*) 194 ! allocation des variables dependant du nombre de site de forage (nbic) 195 if (.not.allocated(iic)) then 196 allocate(iic(nbic),stat=err) 197 if (err/=0) then 198 print *,"erreur a l'allocation du tableau iic dans init_outshort",err 199 stop 4 200 end if 201 end if 202 if (.not.allocated(jic)) then 203 allocate(jic(nbic),stat=err) 204 if (err/=0) then 205 print *,"erreur a l'allocation du tableau jic dans init_outshort",err 206 stop 4 207 end if 208 end if 209 if (.not.allocated(ic_name)) then 210 allocate(ic_name(nbic),stat=err) 211 if (err/=0) then 212 print *,"erreur a l'allocation du tableau ic_name dans init_outshort",err 213 stop 4 214 end if 215 end if 216 if (.not.allocated(ic_lon)) then 217 allocate(ic_lon(nbic),stat=err) 218 if (err/=0) then 219 print *,"erreur a l'allocation du tableau ic_lon dans init_outshort",err 220 stop 4 221 end if 222 end if 223 if (.not.allocated(ic_lat)) then 224 allocate(ic_lat(nbic),stat=err) 225 if (err/=0) then 226 print *,"erreur a l'allocation du tableau ic_lat dans init_outshort",err 227 stop 4 228 end if 229 end if 230 if (.not.allocated(var_icoutput)) then 231 allocate(var_icoutput(nbic,nvaric),stat=err) 232 if (err/=0) then 233 print *,"erreur a l'allocation du tableau var_icoutput dans init_outshort",err 234 stop 4 235 end if 236 end if 237 ! lecture des infos sur les variables : 238 do k=1,nbic 239 read(568,'(a100)') ic_name(k) 240 !ic_name(k)=trim(ic_name(k)) 241 read(568,*) ic_lon(k) 242 read(568,*) ic_lat(k) 243 read(568,*) 244 enddo 245 close(568) 246 ! recherche du point i,j correspondant aux coordonnées lon,lat de chaque site de forage 247 do k=1,nbic 248 dist_lonlat = 1.e7 ! initialisation 249 dist_lonlat_new = 1.e6 ! initialisation 250 ! on part du centre de la grille 251 ! iic(k)=int(nx/2) 252 ! jic(k)=int(ny/2) 253 iic(k)=6 254 jic(k)=6 255 nloop=0 256 do while (dist_lonlat_new .LT. dist_lonlat .AND. nloop .LT. 10000) 257 nloop=nloop+1 258 dist_lonlat = dist_lonlat_new 259 ii = iic(k) 260 jj = jic(k) 261 do j=max(jj-5,1),min(jj+5,ny) 262 do i=max(ii-5,1),min(ii+5,nx) 263 distance = abs(xlong(i,j) - ic_lon(k)) + abs(ylat(i,j) - ic_lat(k)) 264 if (distance .LT. dist_lonlat) then 265 dist_lonlat_new = distance 266 iic(k) = i 267 jic(k) = j 268 endif 269 enddo 270 enddo 271 enddo 272 !print*,'init_outshort : ice core ', trim(ic_name(k)), nloop 273 !print*,'init_outshort lon,lat', ic_lon(k),ic_lat(k) 274 !print*,'init_outshort i,j',iic(k),jic(k) 275 !print*,'init_outshort lon lat GRISLI',xlong(iic(k),jic(k)),ylat(iic(k),jic(k)) 276 enddo 277 278 279 280 ! creation du fichier Netcdf : 281 status=nf90_create(path = trim(dirnameout)//'icecore'//runname//'.nc', cmode = nf90_clobber, ncid = ncidic) 282 if (status /= nf90_noerr) call handle_err(status) 283 284 ! definition des dimension : 285 status=nf90_def_dim(ncidic, name="time", len=NF90_UNLIMITED, dimid=timeDimIDic) 286 if (status /= nf90_noerr) call handle_err(status) 287 status=nf90_def_var(ncidic, name="time", xtype=nf90_float, dimids=(/ timeDimIDic/), varid=timeVarIDic) 288 if (status /= nf90_noerr) call handle_err(status) 289 status=nf90_put_att(ncidic, timeVarIDic, "standard_name", "time") 290 if (status /= nf90_noerr) call handle_err(status) 291 status=nf90_put_att(ncidic, timeVarIDic,"units", "years") 292 if (status /= nf90_noerr) call handle_err(status) 293 status=nf90_put_att(ncidic, timeVarIDic,"calendar", "360_day") 294 if (status /= nf90_noerr) call handle_err(status) 295 296 ! dimension "station" Sites de forages 297 status=nf90_def_dim(ncidic, name="station", len=nbic, dimid=stationDimIDic) 298 if (status /= nf90_noerr) call handle_err(status) 299 status=nf90_def_dim(ncidic, name="name_strlen", len=len(ic_name(1)), dimid=station_strlenDimIDic) 300 if (status /= nf90_noerr) call handle_err(status) 301 status=nf90_def_var(ncidic, name="station_name", xtype=nf90_char, dimids=(/ station_strlenDimIDic, stationDimIDic /), varid=stationVarIDic) 302 if (status /= nf90_noerr) call handle_err(status) 303 status=nf90_put_att(ncidic, stationVarIDic, "standard_name", "station_name") 304 if (status /= nf90_noerr) call handle_err(status) 305 status=nf90_put_att(ncidic, stationVarIDic, "long_name", "station name") 306 if (status /= nf90_noerr) call handle_err(status) 307 status=nf90_def_var(ncidic, name="lon", xtype=nf90_float, dimids=(/ stationDimIDic/), varid=lonVarIDic) 308 if (status /= nf90_noerr) call handle_err(status) 309 status=nf90_put_att(ncidic, lonVarIDic, "standard_name", "longitude") 310 if (status /= nf90_noerr) call handle_err(status) 311 status=nf90_put_att(ncidic, lonVarIDic, "long_name", "station longitude") 312 if (status /= nf90_noerr) call handle_err(status) 313 status=nf90_def_var(ncidic, name="lat", xtype=nf90_float, dimids=(/ stationDimIDic/), varid=latVarIDic) 314 if (status /= nf90_noerr) call handle_err(status) 315 status=nf90_put_att(ncidic, latVarIDic, "standard_name", "latitude") 316 if (status /= nf90_noerr) call handle_err(status) 317 status=nf90_put_att(ncidic, latVarIDic, "long_name", "station latitude") 318 if (status /= nf90_noerr) call handle_err(status) 319 320 321 ! infos sur les variables : 322 open(568,file=trim(dirsource)//'/Fichiers-parametres/icecore-var.dat',status='old') 323 ! lecture en-tete 324 read(568,*) 325 read(568,*) 326 read(568,*) 327 read(568,*) 328 read(568,*) 329 ! lecture des infos sur les variables : 330 do l=1,nvaric 331 read(568,'(a100)') namevaric(l) 332 read(568,'(a100)') standard_nameic(l) 333 read(568,'(a100)') long_nameic(l) 334 read(568,'(a100)') unitsic(l) 335 read(568,*) 336 enddo 337 close(568) 338 339 ! definition des variables de sortie : 340 do l=1,nvaric ! boucle sur le nbr de variable a definir 341 status=nf90_def_var(ncidic, name=trim(namevaric(l)), xtype=nf90_float, dimids= & 342 (/ stationDimIDic, timeDimIDic /), varid=varIDic(l)) 343 if (status /= nf90_noerr) call handle_err(status) 344 status=nf90_put_att(ncidic, varIDic(l), "standard_name", trim(standard_nameic(l))) 345 if (status /= nf90_noerr) call handle_err(status) 346 status=nf90_put_att(ncidic, varIDic(l), "long_name", trim(long_nameic(l))) 347 if (status /= nf90_noerr) call handle_err(status) 348 status=nf90_put_att(ncidic, varIDic(l), "units", trim(unitsic(l))) 349 if (status /= nf90_noerr) call handle_err(status) 350 enddo 351 352 ! fin de la definition du fichier : 353 status=nf90_enddef(ncidic) 354 if (status /= nf90_noerr) call handle_err(status) 355 356 nbtimeoutic = 0 ! initialisation compteur sorties axe time 357 ! ecriture des noms des sites de forage dans la variable site : 358 ! do k=1,nbic 359 ! status=nf90_put_var(ncidic, stationVarIDic, trim(ic_name(k)),start=(/k/),count=(/1/)) 360 ! if (status /= nf90_noerr) call handle_err(status) 361 ! enddo 362 print*,'ic_name',ic_name 363 status=nf90_put_var(ncidic, stationVarIDic, ic_name) !,start=(/1/),count=(/nbic/)) 364 if (status /= nf90_noerr) call handle_err(status) 365 366 status=nf90_put_var(ncidic, lonVarIDic, ic_lon,start=(/1/),count=(/nbic/)) 367 if (status /= nf90_noerr) call handle_err(status) 368 369 status=nf90_put_var(ncidic, latVarIDic, ic_lat,start=(/1/),count=(/nbic/)) 370 if (status /= nf90_noerr) call handle_err(status) 371 ! print*, 'debug 2 init_outshort' 372 endif 373 374 375 376 end subroutine init_outshort 377 378 379 380 !_________________________________________________________________________ 381 subroutine shortoutput 382 383 ! 1_initialization 384 !------------------ 385 real :: smax 386 387 vol=0. 388 np=0 389 nflot=0 390 hmax=0. 391 smax=0. 392 bmean=0. 393 accmean=0. 394 ablmean=0. 395 calvmean=0. 396 ablbordmean=0. 397 bmeltmean=0. 398 tbmean=0. 399 tbdotmean=0. 400 vsmean=0. 401 ! vsdotmean=0. 402 ! uzsmean=0. 403 uzsdotmean=0. 404 uzkmean=0. 405 hdotmean=0. 406 bdotmean=0. 407 volf=0. 408 lim=0. 409 tendacabf=0. 410 iareag=0. 411 iareaf=0. 412 tendlicalvf=0. 413 tendligroundf=0. 414 tendlibmassbf=0. 415 tendlibmassbffl=0. 416 tendlifmassbf=0. 417 ! 2_preparing outputs 418 !-------------------- 189 419 do j=1,ny 190 do i=1,nx 191 if (ice(i,j).eq.1) then ! point englace 192 if (.not.flot(i,j)) then ! point pose 193 np=np+1 194 vol=vol+h(i,j) 195 iareag=iareag+1.*corrsurf(i,j) ! surface englacee posee 196 197 ! calcul de la hauteur au dessus de la flottaison 198 if (sealevel_2d(i,j)-B(i,j).le.0.) then ! socle au dessus du niveau des mers 199 volf=volf+h(i,j)*corrsurf(i,j) ! volume au-dessus de la flottaison 200 else 201 volf=volf+(h(i,j)-row/ro*(sealevel_2d(i,j)-b(i,j)))*corrsurf(i,j) ! volume au-dessus de la flottaison 202 endif 203 204 if (h(i,j).gt.hmax) hmax=h(i,j) 205 if (s(i,j).gt.smax) smax=s(i,j) 206 bmean=bm(i,j)+bmean 207 accmean=acc(i,j)+accmean 208 tbmean=tbmean+t(i,j,nz) 209 tbdotmean=tbdotmean+tbdot(i,j) 210 vsmean=vsmean+sqrt(ux(i,j,1)**2+uy(i,j,1)**2) 211 ! vsdotmean=vsdotmean+vsdot(i,j) 212 ! uzsmean=uzsmean+uz(i,j,1) 213 uzsdotmean=uzsdotmean+uzsdot(i,j) 214 uzkmean=uzkmean+uzk(i,j) 215 hdotmean=hdotmean+abs(hdot(i,j)) 216 bdotmean=bdotmean+abs(bdot(i,j)) 217 bmeltmean=bmeltmean+bmelt(i,j) 218 else ! point flottant 219 iareaf=iareaf+1.*corrsurf(i,j) ! surface flottante 220 tendlibmassbffl=tendlibmassbffl-bmelt_dtt(i,j)*corrsurf(i,j)/dtt_flux ! fonte basale sur zone flottante 420 do i=1,nx 421 if (ice(i,j).eq.1) then ! point englace 422 if (.not.flot(i,j)) then ! point pose 423 np=np+1 424 vol=vol+h(i,j) 425 iareag=iareag+1.*corrsurf(i,j) ! surface englacee posee 426 427 ! calcul de la hauteur au dessus de la flottaison 428 if (sealevel_2d(i,j)-B(i,j).le.0.) then ! socle au dessus du niveau des mers 429 volf=volf+h(i,j)*corrsurf(i,j) ! volume au-dessus de la flottaison 430 else 431 volf=volf+(h(i,j)-row/ro*(sealevel_2d(i,j)-b(i,j)))*corrsurf(i,j) ! volume au-dessus de la flottaison 432 endif 433 434 if (h(i,j).gt.hmax) hmax=h(i,j) 435 if (s(i,j).gt.smax) smax=s(i,j) 436 bmean=bm(i,j)+bmean 437 accmean=acc(i,j)+accmean 438 tbmean=tbmean+t(i,j,nz) 439 tbdotmean=tbdotmean+tbdot(i,j) 440 vsmean=vsmean+sqrt(ux(i,j,1)**2+uy(i,j,1)**2) 441 ! vsdotmean=vsdotmean+vsdot(i,j) 442 ! uzsmean=uzsmean+uz(i,j,1) 443 uzsdotmean=uzsdotmean+uzsdot(i,j) 444 uzkmean=uzkmean+uzk(i,j) 445 hdotmean=hdotmean+abs(hdot(i,j)) 446 bdotmean=bdotmean+abs(bdot(i,j)) 447 bmeltmean=bmeltmean+bmelt(i,j) 448 else ! point flottant 449 iareaf=iareaf+1.*corrsurf(i,j) ! surface flottante 450 tendlibmassbffl=tendlibmassbffl-bmelt_dtt(i,j)*corrsurf(i,j)/dtt_flux ! fonte basale sur zone flottante 451 endif 452 lim=lim+h(i,j)*corrsurf(i,j) ! volume total de glace 453 tendacabf=tendacabf+bm_dtt(i,j)*corrsurf(i,j)/dtt_flux ! smb surface 454 tendlibmassbf=tendlibmassbf-bmelt_dtt(i,j)*corrsurf(i,j)/dtt_flux ! fonte basale 221 455 endif 222 lim=lim+h(i,j)*corrsurf(i,j) ! volume total de glace 223 tendacabf=tendacabf+bm_dtt(i,j)*corrsurf(i,j)/dtt_flux ! smb surface 224 tendlibmassbf=tendlibmassbf-bmelt_dtt(i,j)*corrsurf(i,j)/dtt_flux ! fonte basale 225 endif 226 tendlicalvf=tendlicalvf+calv_dtt(i,j)*corrsurf(i,j)/dtt_flux ! calving 227 tendlifmassbf=tendlifmassbf+(calv_dtt(i,j)-ablbord_dtt(i,j))*corrsurf(i,j)/dtt_flux ! calving + ablbord (ablbord est positif) 228 tendligroundf=tendligroundf+grline_dtt(i,j)*corrsurf(i,j)/dtt_flux ! flux a la grounding line 229 end do 456 tendlicalvf=tendlicalvf+calv_dtt(i,j)*corrsurf(i,j)/dtt_flux ! calving 457 tendlifmassbf=tendlifmassbf+(calv_dtt(i,j)-ablbord_dtt(i,j))*corrsurf(i,j)/dtt_flux ! calving + ablbord (ablbord est positif) 458 tendligroundf=tendligroundf+grline_dtt(i,j)*corrsurf(i,j)/dtt_flux ! flux a la grounding line 459 end do 230 460 end do 231 461 232 462 233 if (np.ne.0) then 234 hmean=vol/np 235 vol=vol*dx*dy 236 volf=volf*dx*dy*ice_density 237 bmean=bmean/np 238 accmean=accmean/np 239 ablmean=bmean-accmean 240 calvmean=calvmean/np 241 bmeltmean=bmeltmean/np 242 ablbordmean=ablbordmean/np 243 tbmean=tbmean/np 244 tbdotmean=tbdotmean/np 245 vsmean=vsmean/np 246 ! vsdotmean=vsdotmean/np 247 ! uzsmean=uzsmean/np 248 uzsdotmean=uzsdotmean/np 249 uzkmean=uzkmean/np 250 hdotmean=hdotmean/np 251 endif 252 lim=lim*dx*dy*ice_density 253 iareag=iareag*dx*dy 254 iareaf=iareaf*dx*dy 255 tendacabf=tendacabf*dx*dy*ice_density/secyear 256 tendlibmassbf=tendlibmassbf*dx*dy*ice_density/secyear 257 tendlibmassbffl=tendlibmassbffl*dx*dy*ice_density/secyear 258 tendlicalvf=tendlicalvf*dx*dy*ice_density/secyear 259 tendlifmassbf=tendlifmassbf*dx*dy*ice_density/secyear 260 tendligroundf=tendligroundf*ice_density/secyear 261 262 263 bdotmean=bdotmean/nx/ny 264 265 266 ! 2_writing outputs 267 !------------------ 268 ! **** short display **** 269 270 write(num_ritz,904) nt,time,tafor,sealevel,vol,volf,np, & 271 nint(hmean),nint(smax), & 272 bmean,tbmean,nint(vsmean), & 273 tbdotmean,hdotmean,dt,accmean, & 274 diff_H,water_bilan,sum(calv_dtt(:,:))/dtt_flux, & 275 sum(ablbord_dtt(:,:))/dtt_flux, & 276 sum(Bm_dtt(:,:))/dtt_flux,sum(bmelt_dtt(:,:))/dtt_flux 277 278 279 903 format(i8,1x,f0.2,1x,f0.4,1x,f0.2,1x,2(e10.5,1x),i6,1x,i4,1x,i5,1x, & 280 f0.4,1x,f0.3,1x,i3,4(1x,e8.2),1x,f0.4,1x,f0.4) 281 904 format(i8,1x,f0.2,4(1x,e15.8),3(1x,i8),2(1x,e15.8),1x,i8,10(1x,e15.8)) 282 !940 format('%%%% ',a,' time=',f8.0,' %%%%') 283 284 285 if (ncshortout.eq.1) then ! ecriture netcdf 286 ! Total ice mass 287 var_shortoutput(1)=lim 288 ! Mass above floatation 289 var_shortoutput(2)=volf 290 ! Grounded ice area 291 var_shortoutput(3)=iareag 292 ! Floating ice area 293 var_shortoutput(4)=iareaf 294 ! Total SMB flux 295 var_shortoutput(5)=tendacabf 296 ! Total Basal mass balance flux 297 var_shortoutput(6)=tendlibmassbf 298 ! Total Basal mass balance flux beneath floating ice 299 var_shortoutput(7)=tendlibmassbffl 300 ! Total calving flux 301 var_shortoutput(8)=tendlicalvf 302 ! Total calving and ice front melting flux 303 var_shortoutput(9)=tendlifmassbf 304 ! Total grounding line flux 305 var_shortoutput(10)=tendligroundf 306 307 nbtimeout=nbtimeout+1 308 309 status=nf90_put_var(ncid, timeVarID, time, start=(/nbtimeout/)) 310 if (status /= nf90_noerr) call handle_err(status) 311 312 do k=1,nvar ! boucle sur le nbr de variable a ecrire 313 status=nf90_put_var(ncid, varID(k), var_shortoutput(k),start=(/nbtimeout/)) 314 if (status /= nf90_noerr) call handle_err(status) 315 enddo 316 status=nf90_sync(ncid) 317 if (status /= nf90_noerr) call handle_err(status) 318 endif 319 320 321 end subroutine shortoutput 322 323 subroutine handle_err(status) 324 integer, intent(in) :: status 325 if (status /= nf90_noerr) then 326 print*,trim(nf90_strerror(status)) 327 stop "stopped" 328 end if 329 end subroutine handle_err 463 if (np.ne.0) then 464 hmean=vol/np 465 vol=vol*dx*dy 466 volf=volf*dx*dy*ice_density 467 bmean=bmean/np 468 accmean=accmean/np 469 ablmean=bmean-accmean 470 calvmean=calvmean/np 471 bmeltmean=bmeltmean/np 472 ablbordmean=ablbordmean/np 473 tbmean=tbmean/np 474 tbdotmean=tbdotmean/np 475 vsmean=vsmean/np 476 ! vsdotmean=vsdotmean/np 477 ! uzsmean=uzsmean/np 478 uzsdotmean=uzsdotmean/np 479 uzkmean=uzkmean/np 480 hdotmean=hdotmean/np 481 endif 482 lim=lim*dx*dy*ice_density 483 iareag=iareag*dx*dy 484 iareaf=iareaf*dx*dy 485 tendacabf=tendacabf*dx*dy*ice_density/secyear 486 tendlibmassbf=tendlibmassbf*dx*dy*ice_density/secyear 487 tendlibmassbffl=tendlibmassbffl*dx*dy*ice_density/secyear 488 tendlicalvf=tendlicalvf*dx*dy*ice_density/secyear 489 tendlifmassbf=tendlifmassbf*dx*dy*ice_density/secyear 490 tendligroundf=tendligroundf*ice_density/secyear 491 492 493 bdotmean=bdotmean/nx/ny 494 495 496 ! 2_writing outputs 497 !------------------ 498 ! **** short display **** 499 500 write(num_ritz,904) nt,time,tafor,sealevel,vol,volf,np, & 501 nint(hmean),nint(smax), & 502 bmean,tbmean,nint(vsmean), & 503 tbdotmean,hdotmean,dt,accmean, & 504 diff_H,water_bilan,sum(calv_dtt(:,:))/dtt_flux, & 505 sum(ablbord_dtt(:,:))/dtt_flux, & 506 sum(Bm_dtt(:,:))/dtt_flux,sum(bmelt_dtt(:,:))/dtt_flux 507 508 509 903 format(i8,1x,f0.2,1x,f0.4,1x,f0.2,1x,2(e10.5,1x),i6,1x,i4,1x,i5,1x, & 510 f0.4,1x,f0.3,1x,i3,4(1x,e8.2),1x,f0.4,1x,f0.4) 511 904 format(i8,1x,f0.2,4(1x,e15.8),3(1x,i8),2(1x,e15.8),1x,i8,10(1x,e15.8)) 512 !940 format('%%%% ',a,' time=',f8.0,' %%%%') 513 514 515 if (ncshortout.eq.1) then ! ecriture netcdf 516 ! Total ice mass 517 var_shortoutput(1)=lim 518 ! Mass above floatation 519 var_shortoutput(2)=volf 520 ! Grounded ice area 521 var_shortoutput(3)=iareag 522 ! Floating ice area 523 var_shortoutput(4)=iareaf 524 ! Total SMB flux 525 var_shortoutput(5)=tendacabf 526 ! Total Basal mass balance flux 527 var_shortoutput(6)=tendlibmassbf 528 ! Total Basal mass balance flux beneath floating ice 529 var_shortoutput(7)=tendlibmassbffl 530 ! Total calving flux 531 var_shortoutput(8)=tendlicalvf 532 ! Total calving and ice front melting flux 533 var_shortoutput(9)=tendlifmassbf 534 ! Total grounding line flux 535 var_shortoutput(10)=tendligroundf 536 537 nbtimeout=nbtimeout+1 538 539 status=nf90_put_var(ncid, timeVarID, time, start=(/nbtimeout/)) 540 if (status /= nf90_noerr) call handle_err(status) 541 542 do k=1,nvar ! boucle sur le nbr de variable a ecrire 543 status=nf90_put_var(ncid, varID(k), var_shortoutput(k),start=(/nbtimeout/)) 544 if (status /= nf90_noerr) call handle_err(status) 545 enddo 546 status=nf90_sync(ncid) 547 if (status /= nf90_noerr) call handle_err(status) 548 endif 549 550 551 ! ecriture netcdf ice core output (tous les 50 ans) 552 if (ncicecoreout.eq.1 .and. (mod(abs(time),dticecoreout).lt.dtmin)) then ! ecriture netcdf 553 nbtimeoutic=nbtimeoutic+1 554 do k=1,nbic ! boucle sur les stations 555 ! Surface elevation 556 var_icoutput(k,1)=S(iic(k),jic(k)) 557 ! Ice thickness 558 var_icoutput(k,2)=H(iic(k),jic(k)) 559 ! Bedrock elevation 560 var_icoutput(k,3)=Bsoc(iic(k),jic(k)) 561 ! Surface Mass Balance 562 var_icoutput(k,4)=Bm(iic(k),jic(k)) 563 ! Basal melt 564 var_icoutput(k,5)=Bmelt(iic(k),jic(k)) 565 ! Surface velocity in x 566 var_icoutput(k,6)=Uxbar(iic(k),jic(k)) 567 ! Surface velocity in y 568 var_icoutput(k,7)=Uybar(iic(k),jic(k)) 569 ! Velocity magnitude 570 var_icoutput(k,8)=sqrt(Uxbar(iic(k),jic(k))**2 + Uybar(iic(k),jic(k))**2) 571 ! Surface Temperature 572 var_icoutput(k,9)=Ts(iic(k),jic(k)) 573 ! Accumulation 574 var_icoutput(k,10)=Acc(iic(k),jic(k)) 575 enddo 576 ! ecriture de time 577 status=nf90_put_var(ncidic, timeVarIDic, time, start=(/nbtimeoutic/)) 578 if (status /= nf90_noerr) call handle_err(status) 579 580 do l=1,nvaric ! boucle sur le nbr de variables a ecrire 581 status=nf90_put_var(ncidic, varIDic(l), var_icoutput(:,l),start=(/1,nbtimeoutic/)) 582 if (status /= nf90_noerr) call handle_err(status) 583 enddo 584 status=nf90_sync(ncidic) 585 if (status /= nf90_noerr) call handle_err(status) 586 endif 587 588 589 590 end subroutine shortoutput 591 592 subroutine handle_err(status) 593 integer, intent(in) :: status 594 if (status /= nf90_noerr) then 595 print*,trim(nf90_strerror(status)) 596 stop "stopped" 597 end if 598 end subroutine handle_err 330 599 331 600 end module output_antarcti_mod
Note: See TracChangeset
for help on using the changeset viewer.
