!> \file GrIce2sea_files/output_Grice2sea-regions_mod.f90 !> \namespace output_grice2sea_mod !! Module de sortie temporelle Grice2sea avec regions !! \author Cat !! \date ... !! @note Used module !! @note - use module3D_phy !! @note - use param_phy_mod !< module output_grice2sea_mod USE module3D_phy USE param_phy_mod use io_netcdf_grisli implicit none !real :: vol ; !integer :: np real :: bmean !< real :: accmean !< accumulation moyenne real :: ablmean !< ablation moyenne real :: calvmean !< moyenne calving real :: ablbordmean !< real :: bmeltmean !< moyenne fusion basale real :: tbmean !< temperature basale moyenne real :: tbdotmean !< moyenne variation / temps temperature basale real :: vsmean !< vitesse de surface moyenne real :: uzsdotmean !< moyenne variation / temps vitesse verticale de surface real :: uzkmean !< moyenne vitesse verticale de surface real :: hdotmean !< moyenne derivee / temps de h real :: bdotmean !< moyenne bedrock derive / temps real :: volf !< volume au dessus de la flottaison ! dimensionnement pour les bassins ou regions integer,parameter :: nbreg_max=20 !< nombre maxi de regions integer :: nbregions !< nombre de region integer :: i_region !< index des regions integer :: scal_np !< nombre de points grounded par region integer, dimension(nbreg_max) :: np_reg !< nombre de points grounded par region integer, dimension(nx,ny) :: mask_region !< carte des region integer, dimension(nbreg_max) :: nbpoint_reg !< nb de points par region real, dimension(nbreg_max) :: vol_reg !< volumes par region (points grounded) real, dimension(nbreg_max) :: volbuoy_reg !< volumes au dessus flottaison par region real, dimension(nbreg_max) :: meanhdot_reg !< hdot moyens par region real, dimension(nbreg_max) :: sigma_hdot_reg !< sigma hdot par region real, dimension(nbreg_max) :: smb_reg !< bilan de masse en surface real :: scal_vol !< volumes par region (points grounded) real :: scal_volbuoy !< volumes au dessus flottaison par region real :: scal_meanhdot !< hdot moyens par region real :: scal_sigma_hdot !< sigma hdot par region real :: scal_smb ! bilan de mass en surface par region integer :: nvol = 250 !< numero du fichier volume integer :: nvolbuoy = 251 !< numero du fichier volume au dessus flottaison integer :: nhdot = 252 !< numero du fichier hdot mean integer :: nsigma_hdot = 253 !< numero du fichier hdot sigma integer :: npoints = 254 !< numero du fichier nb de points integer :: nsmb = 255 !< numero du fichier smb contains !_______________________________________________________________________________ !< initialise les sorties subroutine init_outshort implicit none real*8, dimension(:,:), pointer :: work_tab !< tableau 2d real ecrit dans le fichier character(len=100) :: file_ncdf !< fichier netcdf issue des fichiers .dat character(len=100) :: region_file ! Fichier dans lequel sont les regions character(len=100) :: filout ! Pour fichier ecriture namelist/output_regions/nbregions,region_file if (itracebug.eq.1) call tracebug('output_grice2sea_mod subroutine init_out_short') ndisp = 1 ! ndisp sortie courte tous les ndisp ! lecture des regions ----------------------------------------------------- rewind(num_param) ! pour revenir au debut du fichier param_list.dat 428 format(A) read(num_param,output_regions) write(num_rep_42,428) '!___________________________________________________________' write(num_rep_42,428) '! sorties temporelles par regions' write(num_rep_42,output_regions) write(num_rep_42,428) '! nbregions = nombre de regions' write(num_rep_42,428) '! region_file = carte des regions ' write(num_rep_42,428) '!___________________________________________________________' ! lecture du masque region file_ncdf = trim(dirnameinp)//'test_map_region'//'.nc' region_file = trim(dirnameinp)//region_file write(6,*) file_ncdf write(6,*) region_file ! call lect_input(1,'mask_region',1,work_tab,region_file,file_ncdf) call Read_Ncdf_var('z',region_file,work_tab) mask_region(:,:) = work_tab(:,:) if (mask_region((nx-1)/2,(ny-1)/2).eq.0) then write(6,*) 'mauvaise lecture du fichier region' STOP end if ! ouverture des fichiers par type de variable, toutes les regions dans le meme fichier filout = runname//'_vol_regions.dat' filout = TRIM(DIRNAMEOUT)//TRIM(filout) open(nvol,file=filout) filout = runname// '_volbuoy_regions.dat' filout = TRIM(DIRNAMEOUT)//TRIM(filout) open(nvolbuoy,file=filout) filout = runname//'_meanhdot_regions.dat' filout = TRIM(DIRNAMEOUT)//TRIM(filout) open(nhdot,file=filout) filout = runname//'_sigmahdot_regions.dat' filout = TRIM(DIRNAMEOUT)//TRIM(filout) open(nsigma_hdot,file=filout) filout = runname//'_nbpoints_regions.dat' filout = TRIM(DIRNAMEOUT)//TRIM(filout) open(npoints,file=filout) filout = runname//'_smb_regions.dat' filout = TRIM(DIRNAMEOUT)//TRIM(filout) open(nsmb,file=filout) end subroutine init_outshort !_________________________________________________________________________ subroutine shortoutput implicit none real :: smax if (itracebug.eq.1) call tracebug('output_grice2sea_mod subroutine shortoutput') ! Sorties pour toute la calotte (ancienne version) ! 1_initialization !------------------ vol=0. np=0 hmax=0. smax=0. bmean=0. accmean=0. ablmean=0. calvmean=0. ablbordmean=0. bmeltmean=0. tbmean=0. tbdotmean=0. vsmean=0. uzsdotmean=0. uzkmean=0. hdotmean=0. bdotmean=0. volf=0. ! 2_preparing outputs !-------------------- do i=1,nx do j=1,ny if (.not.flot(i,j)) then ! if (h(i,j).gt.1.) then np=np+1 vol=vol+h(i,j) ! calcul de la hauteur au dessus de la flottaison if (sealevel_2d(i,j)-B(i,j).le.0.) then ! socle au dessus du niveau des mers volf=volf+h(i,j) else volf=volf+h(i,j)-row/ro*(sealevel_2d(i,j)-b(i,j)) endif if (h(i,j).gt.hmax) hmax=h(i,j) if (s(i,j).gt.smax) smax=s(i,j) bmean=bm(i,j)+bmean accmean=acc(i,j)+accmean tbmean=tbmean+t(i,j,nz) tbdotmean=tbdotmean+tbdot(i,j) vsmean=vsmean+sqrt(ux(i,j,1)**2+uy(i,j,1)**2) uzsdotmean=uzsdotmean+uzsdot(i,j) uzkmean=uzkmean+uzk(i,j) hdotmean=hdotmean+abs(hdot(i,j)) bdotmean=bdotmean+abs(bdot(i,j)) bmeltmean=bmeltmean+bmelt(i,j) endif calvmean=calvmean+calv(i,j) ablbordmean=ablbordmean+ablbord(i,j) end do end do if (np.ne.0) then hmean=vol/np vol=vol*dx*dy volf=volf*dx*dy bmean=bmean/np accmean=accmean/np ablmean=bmean-accmean calvmean=calvmean/np bmeltmean=bmeltmean/np ablbordmean=ablbordmean/np tbmean=tbmean/np tbdotmean=tbdotmean/np vsmean=vsmean/np ! vsdotmean=vsdotmean/np ! uzsmean=uzsmean/np uzsdotmean=uzsdotmean/np uzkmean=uzkmean/np hdotmean=hdotmean/np endif bdotmean=bdotmean/nx/ny ! 2_writing outputs !------------------ ! **** short display **** write(num_ritz,903) nt,time,tafor,sealevel,vol,volf,np, & nint(hmean),nint(smax), & bmean,tbmean,nint(vsmean), & ! tbdotmean,vsdotmean,hdotmean,bdotmean, & ! tbdotmean,hdotmean,dt,bmeltmean,accmean tbdotmean,hdotmean,dt,bmelt(3,3),accmean 903 format(i8,1x,f0.2,1x,f0.4,1x,f0.2,1x,2(e15.8,1x),i0,1x,i0,1x,i0,1x, & f0.4,1x,f0.3,1x,i0,4(1x,e11.4),1x,f0.4,1x,f0.4) !940 format('%%%% ',a,' time=',f8.0,' %%%%') call output_regions end subroutine shortoutput !---------------------------------------------------------------------------- subroutine output_regions implicit none character(len=100) :: fmt1,fmt2,fmt3 character(len=2) :: string_region if (itracebug.eq.1) call tracebug('output_grice2sea_mod subroutine output_regions') region: do i_region = 1, nbregions ! initialisation scal_vol = 0. scal_volbuoy = 0. scal_meanhdot = 0. scal_sigma_hdot = 0. scal_np = 0 scal_smb = 0. do i=1,nx do j=1,ny ! write(165,*) i_region,i,j,mask_region(i,j),flot(i,j) ! Attention, pour le Groenland le test porte aussi sur l'epaisseur ! pour le smb, on perd sans doute la demi maille des points ou H=0 if ((mask_region(i,j).eq.i_region).and.(.not.flot(i,j)).and.H(i,j).gt.1.) then scal_np = scal_np + 1 scal_vol = scal_vol + H(i,j) scal_volbuoy = scal_volbuoy + H(i,j) - row/ro * max(0.,sealevel_2d(i,j)-b(i,j)) scal_meanhdot = scal_meanhdot + Hdot(i,j) scal_sigma_hdot = scal_sigma_hdot + Hdot(i,j)**2 scal_smb = scal_smb + bm(i,j) end if end do end do np_reg (i_region) = scal_np vol_reg(i_region) = scal_vol*dx*dx volbuoy_reg(i_region) = scal_volbuoy*dx*dx smb_reg(i_region) = scal_smb*dx*dx if (scal_np.ne.0) then meanhdot_reg(i_region) = scal_meanhdot/scal_np sigma_hdot_reg(i_region) = (scal_sigma_hdot/scal_np)**0.5 else write(6,*) 'region',i_region,' n a aucun point grounded' end if end do region ! ecriture dans les fichiers write(string_region,"(i0)") nbregions+1 fmt1 = '(f10.1,4x,'//string_region//'(es15.8,1x))' fmt2 = '(f10.1,4x,'//string_region//'(i8,1x))' write(string_region,"(i0)") nbregions fmt3 = '(f10.1,4x,'//string_region//'(es15.8,1x))' ! write(6,*) fmt1 ! write(6,*) fmt2 write(nvol,fmt1) time, sum(vol_reg(1:nbregions)), & (vol_reg(i_region) , i_region = 1,nbregions) write(nvolbuoy,fmt1) time, sum(volbuoy_reg(1:nbregions)), & (volbuoy_reg(i_region) , i_region = 1,nbregions) write(nsmb,fmt1) time, sum(smb_reg(1:nbregions)), & (smb_reg(i_region) , i_region = 1,nbregions) write(nhdot,fmt1) time, sum(meanhdot_reg(1:nbregions)), & (meanhdot_reg(i_region) , i_region = 1,nbregions) write(npoints,fmt2) time, sum(np_reg(1:nbregions)), & (np_reg(i_region) , i_region = 1,nbregions) write(nsigma_hdot,fmt3) time, (sigma_hdot_reg(i_region), i_region = 1,nbregions) end subroutine output_regions end module output_grice2sea_mod