Ignore:
Timestamp:
01/20/22 11:39:11 (2 years ago)
Author:
aquiquet
Message:

Add proper scalar outputs to Laure16 geometry

Location:
trunk/SOURCES/Laure16_files
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • trunk/SOURCES/Laure16_files/lect-laurentide_mod.f90

    r318 r355  
    135135  enddo 
    136136  close(2004) 
    137                
     137 
     138  where(xlong(:,:).lt.0.) xlong(:,:)=xlong(:,:)+360 
     139 
    138140  xmin=xcc(1,1)/1000. 
    139141  ymin=ycc(1,1)/1000. 

  • trunk/SOURCES/Laure16_files/output_laure16_mod.f90

    r318 r355  
    11module  output_laure16_mod 
    22 
    3        USE module3D_phy 
    4        use bilan_eau_mod 
     3  use module3D_phy 
     4  use bilan_eau_mod 
     5  use bilan_flux_mod 
     6  use netcdf 
     7  use io_netcdf_grisli 
    58 
    69implicit none 
    710 
    8 real ::  bmean                        ! 
    9 real ::  accmean                      ! accumulation moyenne 
    10 real ::  ablmean                      ! ablation moyenne 
    11 real ::  calvmean                     ! moyenne calving 
    12 real ::  ablbordmean                  ! 
    13 real ::  ablatotmean 
    14 real ::  bmeltmean                    ! moyenne fusion basale 
    15 real ::  tbmean                       ! temperature basale moyenne 
    16 real ::  tbdotmean                    ! moyenne variation / temps temperature basale 
    17 real ::  vsmean                       ! vitesse de surface moyenne 
    18 !real ::  vsdotmean                    ! moyenne variation / temps vitesse de surface 
    19 real ::  uzsmean   !!!! utilise ?     ! vitesse verticale de surface moyenne 
    20 real ::  uzsdotmean                   ! moyenne variation / temps vitesse verticale de surface 
    21 real ::  uzkmean                      ! moyenne vitesse verticale de surface 
    22 real ::  hdotmean                     ! moyenne derivee / temps de H 
    23 real ::  bdotmean                     ! moyenne bedrock derive / temps 
    24 !real ::  pf0mean                      ! moyenne de PF0 
    25 !real ::  pf1mean                      ! moyenne de PF1 
    26 !real ::  evapmean                     ! moyenne de EVAP 
    27 !real ::  pwmean                       ! moyenne PW 
    28 !real ::  pdfmean                      ! moyenne PDF 
    29 logical,dimension(nx,ny,7) :: mask_cal ! masque regions calotte  
    30 integer, dimension(nx,ny) :: write_mask 
    31 CONTAINS 
     11logical,dimension(nx,ny,7) :: mask_cal !< masque regions calotte  
     12real,dimension(nx,ny) :: corrsurf     !< facteur de correction de la surface 
     13 
     14 
     15! variables netcdf 
     16integer,parameter :: ncshortout=1 ! 1 sorties netcdf short initMIP 
     17integer,parameter :: nvar=10 ! nombre de variables dans le fichier de sortie temporel netcdf  
     18integer :: ncid 
     19integer :: status 
     20integer :: timeDimID 
     21integer :: timeVarID 
     22integer,dimension(nvar) :: varID 
     23integer :: nbtimeout  ! index time output 
     24 
     25real,dimension(nvar) :: var_shortoutput 
     26 
     27contains 
    3228!_________________________________________________________________________ 
    3329subroutine init_outshort 
     30 
     31 
     32double precision,dimension(:,:),pointer      :: tab               !< tableau 2d real pointer 
     33character(len=100),dimension(nvar) :: namevar ! name, standard_name, long_name,units 
     34character(len=100),dimension(nvar) :: standard_name 
     35character(len=100),dimension(nvar) :: long_name 
     36character(len=100),dimension(nvar) :: units 
     37 
    3438 
    3539!ndisp sorite courte tous les ndisp 
    … …  
    4145   do i=1,nx 
    4246! -- Laurentide 
    43       IF( (((xlong(i,j).ge.190).AND.(xlong(i,j).lt.210)).AND. & 
    44            ((ylat(i,j).ge.50).AND.(ylat(i,j).le.70))).OR.      & 
    45            (((xlong(i,j).ge.210).AND.(xlong(i,j).lt.220)).AND. & 
    46            ((ylat(i,j).ge.55).AND.(ylat(i,j).le.75))).OR.      & 
    47            (((xlong(i,j).ge.220).AND.(xlong(i,j).lt.250)).AND. &  
    48            ((ylat(i,j).ge.40).AND.(ylat(i,j).le.85))).OR.      & 
    49            (((xlong(i,j).ge.250).AND.(xlong(i,j).lt.290)).AND. & 
    50            ((ylat(i,j).ge.35).AND.(ylat(i,j).le.85))).OR.      & 
    51            (((xlong(i,j).ge.290).AND.(xlong(i,j).lt.300)).AND. & 
    52            ((ylat(i,j).ge.35).AND.(ylat(i,j).le.75))).OR.      & 
    53            (((xlong(i,j).ge.300).AND.(xlong(i,j).le.310)).AND. & 
    54            ((ylat(i,j).ge.35).AND.(ylat(i,j).lt.60))) ) mask_cal(i,j,2)=.true. 
     47      if( (((xlong(i,j).ge.190).and.(xlong(i,j).lt.210)).and. & 
     48           ((ylat(i,j).ge.50).and.(ylat(i,j).le.70))).or.      & 
     49           (((xlong(i,j).ge.210).and.(xlong(i,j).lt.220)).and. & 
     50           ((ylat(i,j).ge.55).and.(ylat(i,j).le.75))).or.      & 
     51           (((xlong(i,j).ge.220).and.(xlong(i,j).lt.250)).and. &  
     52           ((ylat(i,j).ge.40).and.(ylat(i,j).le.85))).or.      & 
     53           (((xlong(i,j).ge.250).and.(xlong(i,j).lt.290)).and. & 
     54           ((ylat(i,j).ge.35).and.(ylat(i,j).le.85))).or.      & 
     55           (((xlong(i,j).ge.290).and.(xlong(i,j).lt.300)).and. & 
     56           ((ylat(i,j).ge.35).and.(ylat(i,j).le.75))).or.      & 
     57           (((xlong(i,j).ge.300).and.(xlong(i,j).le.310)).and. & 
     58           ((ylat(i,j).ge.35).and.(ylat(i,j).lt.60))) ) mask_cal(i,j,2)=.true. 
    5559! -- Groenland 
    56  
     60      if( (((xlong(i,j).ge.290).and.(xlong(i,j).le.310)).and. & 
     61           ((ylat(i,j).ge.75).and.(ylat(i,j).le.85))).or.      & 
     62           (((xlong(i,j).ge.300).and.(xlong(i,j).lt.310)).and. & 
     63           ((ylat(i,j).ge.60).and.(ylat(i,j).le.75))).or.      & 
     64           (((xlong(i,j).ge.310).and.(xlong(i,j).le.345)).and. & 
     65           ((ylat(i,j).ge.54).and.(ylat(i,j).le.85))) ) mask_cal(i,j,3)=.true. 
    5766! -- Labrador Sector  
    58         IF((((xlong(i,j).ge.285).AND.(xlong(i,j).le.300)).AND.  & 
    59             ((ylat(i,j).ge.35).AND.(ylat(i,j).le.70))).OR.      & 
    60             (((xlong(i,j).ge.300).AND.(xlong(i,j).lt.310)).AND. & 
    61             ((ylat(i,j).ge.35).AND.(ylat(i,j).lt.60)))) mask_cal(i,j,4)=.true. 
     67      if((((xlong(i,j).ge.285).and.(xlong(i,j).le.300)).and.  & 
     68           ((ylat(i,j).ge.35).and.(ylat(i,j).le.70))).or.      & 
     69           (((xlong(i,j).ge.300).and.(xlong(i,j).lt.310)).and. & 
     70           ((ylat(i,j).ge.35).and.(ylat(i,j).lt.60)))) mask_cal(i,j,4)=.true. 
    6271! -- Keewatin Sector  
    63         IF(  ((xlong(i,j).gt.240).AND.(xlong(i,j).lt.285)).AND. & 
    64              ((ylat(i,j).ge.35).AND.(ylat(i,j).le.70)) ) mask_cal(i,j,5)=.true. 
     72      if(  ((xlong(i,j).gt.240).and.(xlong(i,j).lt.285)).and. & 
     73           ((ylat(i,j).ge.35).and.(ylat(i,j).le.70)) ) mask_cal(i,j,5)=.true. 
    6574! -- Innuitian Ice Sheet 
    66         IF( ((xlong(i,j).gt.230).AND.(xlong(i,j).lt.290)).AND. & 
    67             ((ylat(i,j).gt.70).AND.(ylat(i,j).le.85))) mask_cal(i,j,6)=.true. 
     75      if( ((xlong(i,j).gt.230).and.(xlong(i,j).lt.290)).and. & 
     76           ((ylat(i,j).gt.70).and.(ylat(i,j).le.85))) mask_cal(i,j,6)=.true. 
    6877! -- Cordilleran Ice Sheet 
    69         IF( (((xlong(i,j).ge.190).AND.(xlong(i,j).lt.210)).AND. & 
    70             ((ylat(i,j).ge.50).AND.(ylat(i,j).le.70))).OR.      & 
    71             (((xlong(i,j).ge.210).AND.(xlong(i,j).lt.220)).AND. & 
    72             ((ylat(i,j).ge.55).AND.(ylat(i,j).le.75))).OR.      & 
    73             (((xlong(i,j).ge.210).AND.(xlong(i,j).lt.220)).AND. & 
    74              ((ylat(i,j).ge.55).AND.(ylat(i,j).le.75))).OR.     & 
    75             (((xlong(i,j).ge.220).AND.(xlong(i,j).le.240)).AND. & 
    76              ((ylat(i,j).ge.40).AND.(ylat(i,j).le.70))) ) mask_cal(i,j,7)=.true. 
    77 ! -- Groenland  
    78         IF( (((xlong(i,j).ge.290).AND.(xlong(i,j).le.310)).AND. & 
    79             ((ylat(i,j).ge.75).AND.(ylat(i,j).le.85))).OR.      & 
    80             (((xlong(i,j).ge.300).AND.(xlong(i,j).lt.310)).AND. & 
    81             ((ylat(i,j).ge.60).AND.(ylat(i,j).le.75))).OR.      & 
    82 !           (((xlong(i,j).ge.310).AND.(xlong(i,j).le.350)).AND. & 
    83 !           ((ylat(i,j).ge.40).AND.(ylat(i,j).le.85))) ) THEN 
    84             (((xlong(i,j).ge.310).AND.(xlong(i,j).le.345)).AND. & 
    85             ((ylat(i,j).ge.54).AND.(ylat(i,j).le.85))) ) mask_cal(i,j,3)=.true. 
    86        enddo 
    87     enddo 
    88 !        open (UNIT=7777,file='mask_hemin40-02.dat') 
    89 !   do j=1,ny 
    90 !        do i=1,nx 
    91 !        if (mask_cal(i,j,2)) then 
    92 !                write_mask(i,j)=2 
    93 !        else 
    94 !                write_mask(i,j)=0 
    95 !        endif     
    96 !        write(7777,*) i,j,write_mask(i,j) 
    97 !        enddo      
    98 !   enddo 
    99 !   close(UNIT=7777) 
    100 !   stop 
    101  
     78      if( (((xlong(i,j).ge.190).and.(xlong(i,j).lt.210)).and. & 
     79           ((ylat(i,j).ge.50).and.(ylat(i,j).le.70))).or.      & 
     80           (((xlong(i,j).ge.210).and.(xlong(i,j).lt.220)).and. & 
     81           ((ylat(i,j).ge.55).and.(ylat(i,j).le.75))).or.      & 
     82           (((xlong(i,j).ge.210).and.(xlong(i,j).lt.220)).and. & 
     83           ((ylat(i,j).ge.55).and.(ylat(i,j).le.75))).or.     & 
     84           (((xlong(i,j).ge.220).and.(xlong(i,j).le.240)).and. & 
     85           ((ylat(i,j).ge.40).and.(ylat(i,j).le.70))) ) mask_cal(i,j,7)=.true. 
     86   enddo 
     87enddo 
     88 
     89! hemin40    : time,isvol,deltavol,inp,isvolf,inf 
    10290! ecriture entete avec les colonnes du fichier short : 
    103 WRITE(num_ritz, '(a91)')"! hemin40    : time,isvol,deltavol,inp,isvolf,inf,isacc,isabl,ISABLBORD,ABLATOT,ISCALV,ISBM" 
    104 WRITE(num_ritz, '(a90)')"! colonne    :  1    2      3       4    5     6    7     8      9        10      11    12" 
    105 WRITE(num_ritz,'(a157)')"! calotte    : Laurentide Groenland LabradorSector KeewatinSector InnuitianIS CordilleranIS" 
    106 WRITE(num_ritz,'(a153)')"! isvol      :    13        26          39            52            65" 
    107 WRITE(num_ritz,'(a153)')"! inp        :    14        27          40            53            66" 
    108 WRITE(num_ritz,'(a153)')"! isvolf     :    15        28          41            54            67" 
    109 WRITE(num_ritz,'(a153)')"! inf        :    16        29          42            55            68" 
    110 WRITE(num_ritz,'(a153)')"! isacc      :    17        30          43            56            69" 
    111 WRITE(num_ritz,'(a153)')"! isabl      :    18        31          44            57            70" 
    112 WRITE(num_ritz,'(a153)')"! isablbord  :    19        32          45            58            71" 
    113 WRITE(num_ritz,'(a153)')"! ablatot    :    20        33          46            59            72" 
    114 WRITE(num_ritz,'(a153)')"! iscalv     :    21        34          47            60            73" 
    115 WRITE(num_ritz,'(a153)')"! isbm       :    22        35          48            61            74" 
    116 WRITE(num_ritz,'(a153)')"! itjjamean_ :    23        36          49            62            75" 
    117 WRITE(num_ritz,'(a153)')"! hmean_     :    24        37          50            63            76" 
    118 WRITE(num_ritz,'(a153)')"! Hmax_      :    25        38          51            64            77" 
    119 WRITE(num_ritz,*) 
     91write(num_ritz, '(a45)')"! laure16    : time,isvol,inp,isvolf,inf,volf" 
     92write(num_ritz, '(a82)')"! colonne    :    1         2           3             4             5            6" 
     93write(num_ritz,'(a91)')"! calotte    : Laurentide Groenland LabradorSector KeewatinSector InnuitianIS CordilleranIS" 
     94write(num_ritz,'(a83)')"! isvol      :    7        12          17            22            27            32" 
     95write(num_ritz,'(a83)')"! inp        :    8        13          18            23            28            33" 
     96write(num_ritz,'(a83)')"! isvolf     :    9        14          19            24            29            34" 
     97write(num_ritz,'(a83)')"! inf        :    10       15          20            25            30            35" 
     98write(num_ritz,'(a83)')"! volaf      :    11       16          21            26            31            36" 
     99write(num_ritz,'(a61)')"! water bilan: diff_H,water_bilan, calving, ablbord, Bm, bmelt" 
     100write(num_ritz,'(a83)')"! colonne    :    37       38          39            40            41            42" 
     101write(num_ritz,*) 
     102 
     103 
     104if (ncshortout.eq.1) then  ! ecriture netcdf 
     105 
     106! lecture du fichier avec les corrections de surface 
     107   if (geoplace.eq.'ant16km') then 
     108      call Read_Ncdf_var('z',trim(DIRNAMEINP)//'/corrsurf-initMIP-16km.grd',tab) 
     109      corrsurf(:,:)=tab(:,:) 
     110   else 
     111      corrsurf(:,:)= 1. 
     112   end if 
     113 
     114   open(568,file=trim(dirsource)//'/Fichiers-parametres/short-ISMIPnc.dat',status='old') 
     115! lecture en-tete 
     116   read(568,*) 
     117   read(568,*) 
     118   read(568,*) 
     119   read(568,*) 
     120   read(568,*) 
     121! lecture des infos sur les variables : 
     122   do k=1,nvar 
     123      read(568,'(a100)') namevar(k) 
     124      read(568,'(a100)') standard_name(k) 
     125      read(568,'(a100)') long_name(k) 
     126      read(568,'(a100)') units(k) 
     127      read(568,*) 
     128   enddo 
     129   close(568) 
     130! Fichier Netcdf initMIP 
     131! creation du fichier Netcdf : 
     132   status=nf90_create(path = trim(dirnameout)//'short'//runname//'.nc', cmode = nf90_clobber, ncid = ncid) 
     133   if (status /= nf90_noerr) call handle_err(status) 
     134 
     135! definition des dimension : 
     136   status=nf90_def_dim(ncid, name="time", len=NF90_UNLIMITED, dimid=timeDimID) 
     137   if (status /= nf90_noerr) call handle_err(status) 
     138   status=nf90_def_var(ncid, name="time", xtype=nf90_float, dimids=(/ timeDimID/), varid=timeVarID)  
     139   if (status /= nf90_noerr) call handle_err(status) 
     140   status=nf90_put_att(ncid, timeVarID, "standard_name", "time") 
     141   if (status /= nf90_noerr) call handle_err(status) 
     142   status=nf90_put_att(ncid, timeVarID,"units", "years since 1995-01-01 00:00:00") 
     143   if (status /= nf90_noerr) call handle_err(status) 
     144   status=nf90_put_att(ncid, timeVarID,"calendar", "360_day") 
     145   if (status /= nf90_noerr) call handle_err(status) 
     146   
     147! definition des variables de sortie : 
     148   do k=1,nvar  ! boucle sur le nbr de variable a definir 
     149      status=nf90_def_var(ncid, name=trim(namevar(k)), xtype=nf90_float, dimids= & 
     150           (/ timeDimID /), varid=varID(k)) 
     151      if (status /= nf90_noerr) call handle_err(status) 
     152      status=nf90_put_att(ncid, varID(k), "standard_name", trim(standard_name(k))) 
     153      if (status /= nf90_noerr) call handle_err(status) 
     154      status=nf90_put_att(ncid, varID(k), "long_name", trim(long_name(k))) 
     155      if (status /= nf90_noerr) call handle_err(status) 
     156      status=nf90_put_att(ncid, varID(k), "units", trim(units(k))) 
     157      if (status /= nf90_noerr) call handle_err(status) 
     158   enddo 
     159 
     160! fin de la definition du fchier : 
     161   status=nf90_enddef(ncid) 
     162   if (status /= nf90_noerr) call handle_err(status) 
     163   nbtimeout = 0 ! initialisation compteur sorties axe time 
     164   
     165else ! pas de sortie netcdf et sans correction de surface 
     166   corrsurf(:,:)=1. 
     167endif 
     168 
    120169end subroutine init_outshort 
    121170 
    … …  
    125174! 1_initialization 
    126175!------------------ 
    127       integer KK 
    128 !     integer inpl, INPG, INPF 
    129 !     integer inplab,inpkeew,inpinn,inpcord 
    130 !     integer inpscaN,inpbar,inpkara,inpsib,inpbrit 
    131      integer inp(7) ! Surface posee (nb de noeuds) 
    132      integer inf(7) ! surface flottante 
    133  
    134 !Variables pour sommer       
    135       real isvol(7),isvolf(7) ! volume posé et flottants 
    136       REAL ISCALV(7),ISACC(7),ISBM(7),ISABL(7) 
    137       REAL ISABLBORD(7),ABLATOT(7),TACC(7),TBM(7) 
    138       REAL ITJJA(7) 
    139 !      REAL ABLAMEAN 
    140 !moyennes utilisées en output 
    141       REAL HMAX_(7) , HMEAN_(7)  
    142       REAL BMEAN_(7) , ACCMEAN_(7) , ABLMEAN_(7) , ABLBORDMEAN_(7) 
    143       REAL ABLATOTMEAN_(7) , CALVMEAN_(7) , ITJJAMEAN_(7) 
    144  
    145  
    146  
    147 !     REAL ITJJAMEAN_L, ITJJAMEAN_G, ITJJAMEAN_F ... 
    148 !     REAL ITJJAMEAN_LAB,ITJJAMEAN_KEEW,ITJJAMEAN_INN,ITJJAMEAN_CORD 
    149 !     REAL ITJJAMEAN_SCAN,ITJJAMEAN_BAR,ITJJAMEAN_KARA,ITJJAMEAN_SIB,ITJJAMEAN_BRIT 
    150  
    151       REAL HDOTMEAN_G 
    152 !      REAL ABLA(NX,NY) 
    153       REAL DELTAVOL 
    154  
    155  
    156 !     open(unit=4145,file='reg_output_nord.dat') 
    157  
    158       BMELTMEAN=0. 
    159       TBMEAN=0. 
    160       TBDOTMEAN=0. 
    161       VSMEAN=0. 
    162       UZSMEAN=0. 
    163       UZSDOTMEAN=0. 
    164       UZKMEAN=0. 
    165       HDOTMEAN=0. 
    166       HDOTMEAN_G = 0. 
    167       BDOTMEAN=0. 
    168  
    169       DO kk = 1,13 
    170         INP(kk) = 0 
    171         INF(kk) = 0 
    172         ISVOL(kk) = 0. 
    173         ISVOLF(kk) = 0. 
    174         ISBM(kk) = 0. 
    175         ISACC(kk) = 0. 
    176         ISCALV(kk) = 0. 
    177         ISABL(kk) = 0. 
    178         ISABLBORD(kk) = 0. 
    179         ABLATOT(kk)= 0. 
    180         TACC(kk) = 0. 
    181         TBM(kk) = 0. 
    182         ITJJA(kk) = 0. 
    183          
    184 ! nouveau tof mai 2009 
    185 !        where (mask_cal(:,:,kk).and.H(:,:).gt.2..and.flot(:,:)) ISVOLF(kk) = ISVOLF(kk) + H(I,J) 
    186         isvolf(kk) = sum(H(:,:),mask=(mask_cal(:,:,kk).and.(H(:,:).gt.2.).and.flot(:,:))) 
    187         INF(kk) = count(mask_cal(:,:,kk).and.(H(:,:).gt.2.).and.flot(:,:)) 
    188  
    189         INP(kk) = count(mask_cal(:,:,kk).and.(H(:,:).gt.2.).and..not.flot(:,:)) 
    190         isvol(kk) = sum(H(:,:),mask=(mask_cal(:,:,kk).and.(H(:,:).gt.2.).and..not.flot(:,:))) 
    191         isacc(kk) = sum(Acc(:,:),mask=(mask_cal(:,:,kk).and.(H(:,:).gt.2.).and..not.flot(:,:))) 
    192         isabl(kk) = sum(Abl(:,:),mask=(mask_cal(:,:,kk).and.(H(:,:).gt.2.).and..not.flot(:,:))) 
    193         itjja(kk) = sum(Tjuly(:,:),mask=(mask_cal(:,:,kk).and.(H(:,:).gt.2.).and..not.flot(:,:))) 
    194         Hmax_(kk) = maxval(H(:,:),mask=(mask_cal(:,:,kk).and.(H(:,:).gt.2.).and..not.flot(:,:))) 
    195         Tacc(kk) = sum(Acc(:,:),mask=(mask_cal(:,:,kk).and.(H(:,:).gt.2.))) 
    196         isablbord(kk) = sum(ablbord(:,:),mask=(mask_cal(:,:,kk))) 
    197         iscalv(kk) =  sum(calv(:,:),mask=(mask_cal(:,:,kk).and.(H(:,:).gt.2.)))/dt 
    198  
    199         ablatot(kk) = isabl(kk) + isablbord(kk) 
    200         isbm(kk) = isacc(kk)+isabl(kk)+iscalv(kk)+isablbord(kk) 
    201         Tbm(kk) = Tacc(kk)+isabl(kk)+iscalv(kk)+isablbord(kk) 
    202  
    203  
    204      enddo 
    205  
    206  
    207 !!$ 
    208 !!$ 
    209 !!$            IF(H(i,j).gt.2.)  THEN  
    210 !!$             if (mk(i,j).eq.1) then 
    211 !!$               INF(2) = INF(2) + 1  
    212 !!$               ISVOLF(2) = ISVOLF(2) + H(I,J) 
    213 !!$             else 
    214 !!$               INP(2) = INP(2) + 1 
    215 !!$               ISVOL(2) = ISVOL(2) + H(I,J) 
    216 !!$               ISACC(2) = ISACC(2) + ACC(I,J) 
    217 !!$               ISABL(2) = ISABL(2) + ABL(I,J) 
    218 !!$               ITJJA(2) = ITJJA(2) + TJULY(I,J) 
    219 !!$               if (H(I,J).gt.HMAX_(2)) HMAX_(2)=H(I,J)  
    220 !!$             endif 
    221 !!$            ENDIF 
    222 !!$            TACC(2) = TACC(2) + ACC(I,J) 
    223 !!$             ISABLBORD(2) = ISABLBORD(2) + ABLBORD(I,J) 
    224 !!$             ABLATOT(2) = ISABL(2) + ISABLBORD(2) 
    225 !!$             ISCALV(2) = ISCALV(2) + CALV(I,J) 
    226 !!$             ISBM(2) = ISACC(2)+ISABL(2)+ISCALV(2)+ISABLBORD(2) 
    227 !!$             TBM(2) = TACC(2)+ISABL(2)+ISCALV(2)+ISABLBORD(2) 
    228 !!$        ENDIF 
    229 !!$ 
    230 !!$ 
    231 !!$      END DO   
    232  
    233  
    234    
    235 do K=1,7 
    236  
    237 ! == Les moyennes       
    238       IF(INP(K).ne.0)THEN 
    239           HMEAN_(K) = ISVOL(K)   /INP(K) ! /(DX*DY*INP(K)) 
    240           BMEAN_(K) = ISBM(K)    /INP(K) 
    241           ACCMEAN_(K) = ISACC(K) /INP(K) 
    242           ABLMEAN_(K) = ISABL(K) /INP(K) 
    243           ABLBORDMEAN_(K) = ISABLBORD(K) /INP(K) 
    244           CALVMEAN_(K)    = ISCALV(K)    /INP(K) 
    245           ABLATOTMEAN_(K) = ABLATOT(K)   /INP(K) 
    246           ITJJAMEAN_(K)   = ITJJA(K)     /INP(K) 
    247        ENDIF 
    248 ! == Les volmes intergrées (3D) 
    249       ISVOL(K) = ISVOL(K)*DX*DY      
    250               isacc(k)=isacc(k)*dx*dy 
    251               isabl(k)=isabl(k)*dx*dy 
    252               isablbord(k)=isablbord(k)*dx*dy 
    253               ablatot(k)=ablatot(k)*dx*dy 
    254               iscalv(k)=iscalv(k)*dx*dy 
    255               isbm(k)=isbm(k)*dx*dy 
    256               tacc(k)=tacc(k)*dx*dy 
    257               tbm(k)=tbm(k)*dx*dy 
    258 enddo 
    259 !nom_table='abl' 
    260 !call printtable_r(abl,nom_table) 
    261  
    262 !==================================================== 
    263 ! -- Hemin40                : : 1    
    264  
    265       if (NP.ne.0) then  
    266 !        HMEAN=VOL/NP  
    267 !        VOL=VOL*DX*DY  
    268 !        ISVOL(1) = VOL 
    269 !        BMEAN=BMEAN/NP  
    270 !        ACCMEAN=ACCMEAN/NP  
    271 !        ITJJAMEAN = ITJJA(1)/NP 
    272 !        ABLMEAN=BMEAN-ACCMEAN  
    273 !!!     ABLMEAN = ABLMEAN/NP 
    274 !        CALVMEAN=CALVMEAN/NP  
    275 !!!     BMELTMEAN=BMELTMEAN/NP 
    276 !!!     ABLBORDMEAN=ABLBORDMEAN/NP 
    277 !!!     ABLAMEAN = ABLAMEAN/NP 
    278         TBMEAN=TBMEAN/NP 
    279         TBDOTMEAN=TBDOTMEAN/NP 
    280         VSMEAN=VSMEAN/NP 
    281 !       VSDOTMEAN=VSDOTMEAN/NP 
    282         UZSMEAN=UZSMEAN/NP 
    283         UZSDOTMEAN=UZSDOTMEAN/NP 
    284         UZKMEAN=UZKMEAN/NP 
    285         DELTAVOL = HDOTMEAN*DX*DY 
    286         HDOTMEAN=HDOTMEAN/NP  
    287     ELSE 
    288         TBMEAN=0. 
    289         TBDOTMEAN=0. 
    290         VSMEAN=0. 
    291         UZSMEAN=0. 
    292         UZSDOTMEAN=0. 
    293         UZKMEAN=0. 
    294         DELTAVOL = 0. 
    295         HDOTMEAN=0. 
    296       endif  
    297  
    298       BDOTMEAN=BDOTMEAN/NX/NY   
     176 
     177  ! ISMIP type outputs: (nc file) 
     178  real :: lim                           !< Total ice mass 
     179  real :: volf                          !< volume above flotation (kg) 
     180  real :: iareag                        !< surface posee 
     181  real :: iareaf                        !< surface flottante 
     182  real :: tendacabf                     !< Total SMB flux 
     183  real :: tendlibmassbf                 !< Total Basal mass balance flux 
     184  real :: tendlibmassbffl               !< Total Basal mass balance flux on floating ice 
     185  real :: tendlicalvf                   !< Total calving flux 
     186  real :: tendlifmassbf                 !< Total calving and ice front melting flux 
     187  real :: tendligroundf                 !< Total grounding line flux 
     188 
     189  ! old-version outputs: (ritz file) 
     190  real inp(7) ! Surface posee (nb de noeuds) 
     191  real inf(7) ! surface flottante 
     192  real isvol(7),isvolf(7) ! volume posé et flottants      
     193  real volfm(7)           ! volume above floatation (m3) 
     194 
     195 
     196      do k = 1,7 
     197          
     198         inp(k) = sum(corrsurf(:,:),mask=(mask_cal(:,:,k).and.(ice(:,:).eq.1).and..not.flot(:,:))) 
     199         inf(k) = sum(corrsurf(:,:),mask=(mask_cal(:,:,k).and.(ice(:,:).eq.1).and.flot(:,:))) 
     200         isvol(k) = sum(h(:,:)*corrsurf(:,:),mask=(mask_cal(:,:,k).and.(ice(:,:).eq.1).and..not.flot(:,:))) 
     201         isvolf(k) = sum(h(:,:)*corrsurf(:,:),mask=(mask_cal(:,:,k).and.(ice(:,:).eq.1).and.flot(:,:))) 
     202 
     203         volfm(k) = sum(h(:,:)*corrsurf(:,:), mask=(mask_cal(:,:,k).and.(ice(:,:).eq.1).and.(sealevel_2d(:,:)-B(:,:).le.0.))) 
     204         volfm(k) = volfm(k)+sum ((h(:,:)-row/ro*(sealevel_2d(:,:)-b(:,:)))*corrsurf(:,:), mask=(mask_cal(:,:,k).and.(ice(:,:).eq.1).and.(sealevel_2d(:,:)-B(:,:).gt.0.))) 
     205          
     206      enddo 
     207 
     208      inp(:)    = inp(:)*dx*dy 
     209      inf(:)    = inf(:)*dx*dy 
     210      isvol(:)  = isvol(:)*dx*dy 
     211      isvolf(:) = isvol(:)*dx*dy 
     212      volfm(:)  = volfm(:)*dx*dy 
     213       
     214 
     215      !lim=sum(h(:,:)*corrsurf(:,:),mask=(ice(:,:).eq.1))                                              ! grounded volume 
     216      !iareag=sum(corrsurf(:,:),mask=(ice(:,:).eq.1).and.(.not.flot(:,:)))                             ! grounded area 
     217      !iareaf=sum(corrsurf(:,:),mask=(ice(:,:).eq.1).and.flot(:,:))                                    ! floating area 
     218      tendacabf=sum(corrsurf(:,:)*bm_dtt(:,:),mask=(ice(:,:).eq.1))/dtt_flux                           ! smb 
     219      tendlibmassbf=sum(-corrsurf(:,:)*bmelt_dtt(:,:),mask=(ice(:,:).eq.1))/dtt_flux                   ! basal melt 
     220      tendlibmassbffl=sum(-corrsurf(:,:)*bmelt_dtt(:,:),mask=((ice(:,:).eq.1).and.flot(:,:)))/dtt_flux ! subshelf melt 
     221      tendlicalvf=sum(corrsurf(:,:)*calv_dtt(:,:))/dtt_flux                                            ! calving flux 
     222      tendlifmassbf=sum(corrsurf(:,:)*(calv_dtt(:,:)-ablbord_dtt))/dtt_flux                            ! calving+ablbord flux 
     223      tendligroundf=sum(corrsurf(:,:)*grline_dtt(:,:))/dtt_flux                                        ! grounding line flux 
     224      
     225      lim=isvol(1)*ice_density !lim*dx*dy*ice_density 
     226      iareag=inp(1)            !iareag*dx*dy  
     227      iareaf=inf(1)            !iareaf*dx*dy 
     228      volf=volfm(1)*ice_density 
     229      tendacabf=tendacabf*dx*dy*ice_density/secyear 
     230      tendlibmassbf=tendlibmassbf*dx*dy*ice_density/secyear 
     231      tendlibmassbffl=tendlibmassbffl*dx*dy*ice_density/secyear 
     232      tendlicalvf=tendlicalvf*dx*dy*ice_density/secyear 
     233      tendlifmassbf=tendlifmassbf*dx*dy*ice_density/secyear 
     234      tendligroundf=tendligroundf*ice_density/secyear 
     235 
    299236 
    300237! 2_writing outputs 
    301238!------------------ 
    302239 
    303 !         WRITE(35,904)TIME,VOL,DELTAVOL,NP,ISACC(1),ISABL(1),      & 
    304 !          ISABLBORD(1),ABLATOT(1),ISCALV(1),ISBM(1),               & 
    305 !          ISVOL(2),INPL,ISACC(2),                                  & 
    306 !          ISABL(2),ISABLBORD(2),ABLATOT(2),ISCALV(2),ISBM(2),      & 
    307 !          ISVOL(3),INPG,ISACC(3),ISABL(3),ISABLBORD(3),            & 
    308 !          ABLATOT(3),ISCALV(3),ISBM(3),ISVOL(4),INPF,ISACC(4),     & 
    309 !          ISABL(4),ISABLBORD(4),ABLATOT(4),ISCALV(4),ISBM(4),      & 
    310 !          nint(HMEAN),nint(HMAX),BMEAN,TBMEAN,nint(VSMEAN),        & 
    311 !          TBDOTMEAN,HDOTMEAN,BDOTMEAN,BMELTMEAN,NPAB,    & 
    312 !!          TBDOTMEAN,VSDOTMEAN,HDOTMEAN,BDOTMEAN,BMELTMEAN,NPAB,    & 
    313 !          NPCAL,HDOTMEAN_G,TACC(1),TBM(1),TACC(2),TBM(2),          & 
    314 !          TACC(3),TBM(3),TACC(4),TBM(4),                           & 
    315 !           ITJJAMEAN,ITJJAMEAN_L,ITJJAMEAN_G,ITJJAMEAN_F,          & 
    316 !           ISVOL(5),INPLAB,ISACC(5),ABLATOT(5),ITJJAMEAN_LAB,      & 
    317 !           ISVOL(6),INPKEEW,ISACC(6),ABLATOT(6),ITJJAMEAN_KEEW,    & 
    318 !           ISVOL(7),INPINN,ISACC(7),ABLATOT(7),ITJJAMEAN_INN,      &  
    319 !           ISVOL(8),INPCORD,ISACC(8),ABLATOT(8),ITJJAMEAN_CORD,    &  
    320 !           ISVOL(9),INPSCAN,ISACC(9),ABLATOT(9),ITJJAMEAN_SCAN,    &  
    321 !           ISVOL(10),INPBAR,ISACC(10),ABLATOT(10),ITJJAMEAN_BAR,   &  
    322 !           ISVOL(11),INPKARA,ISACC(11),ABLATOT(11),ITJJAMEAN_KARA, &      
    323 !           ISVOL(12),INPSIB,ISACC(12),ABLATOT(12),ITJJAMEAN_SIB,   &      
    324 !           ISVOL(13),INPBRIT,ISACC(13),ABLATOT(13),ITJJAMEAN_BRIT, & 
    325 !           HMEAN,HMEAN_L,HMEAN_G,HMEAN_F,HMEAN_LAB,HMEAN_KEEW,     & 
    326 !           HMEAN_INN,HMEAN_CORD,HMEAN_SCAN,HMEAN_BAR,HMEAN_KARA,   & 
    327 !           HMEAN_SIB,HMEAN_BRIT,                                   & 
    328 !           HMAX,HMAX_L,HMAX_G,HMAX_F,HMAX_LAB,HMAX_KEEW,HMAX_INN,  & 
    329 !           HMAX_CORD, HMAX_SCAN,HMAX_BAR,HMAX_KARA,HMAX_SIB,HMAX_BRIT      
    330  
    331  
    332  
    333  
    334 !print*,'=========================' 
    335 !print*,'ecrit short, time = ',time,dt 
    336 !print*,'hdot ISABLBORD',hdot(30:40,100),ISABLBORD(2) 
    337  
    338           write(num_ritz,905)   time,isvol(1),deltavol,inp(1),isvolf(1),inf(1),isacc(1),isabl(1),    &   ! 8 
    339               ISABLBORD(1),ABLATOT(1),ISCALV(1),ISBM(1),     &                                     ! 4 
    340   isvol(2),inp(2),isvolf(2),inf(2),isacc(2),isabl(2),isablbord(2),ablatot(2),iscalv(2),isbm(2),itjjamean_(2),hmean_(2),Hmax_(2), & ! 13 
    341   isvol(3),inp(3),isvolf(3),inf(3),isacc(3),isabl(3),isablbord(3),ablatot(3),iscalv(3),isbm(3),itjjamean_(3),hmean_(3),hmax_(3), & ! 13 
    342   isvol(4),inp(4),isvolf(4),inf(4),isacc(4),isabl(4),isablbord(4),ablatot(4),iscalv(4),isbm(4),itjjamean_(4),hmean_(4),hmax_(4), & ! 13 
    343   isvol(5),inp(5),isvolf(5),inf(5),isacc(5),isabl(5),isablbord(5),ablatot(5),iscalv(5),isbm(5),itjjamean_(5),hmean_(5),hmax_(5), & ! 13 
    344   isvol(6),inp(6),isvolf(6),inf(6),isacc(6),isabl(6),isablbord(6),ablatot(6),iscalv(6),isbm(6),itjjamean_(6),hmean_(6),hmax_(6), & ! 13 
    345   isvol(7),inp(7),isvolf(7),inf(7),isacc(7),isabl(7),isablbord(7),ablatot(7),iscalv(7),isbm(7),itjjamean_(7),hmean_(7),hmax_(7), & ! 13 
    346   diff_H, water_bilan, sum(calv_dtt(:,:))/dtt, sum(ablbord_dtt(:,:))/dtt, sum(Bm_dtt(:,:))/dtt, sum(bmelt_dtt(:,:))/dtt !6 
    347  
    348   
    349 ! !! format 900 faux, a reprendre 
    350  
    351 905 format(f9.1,1x, e11.4,1x,e11.5, 1x,i5, 1x,e10.4,1x,i5 , 6(1x,e12.5), & 
    352          6(2(1x,e10.4,1x,i5), 9(1x,e11.4) ) , 6(1x,e15.8))  
    353   
    354      
    355 940   format('%%%% ',a,'   time=',f8.0,' %%%%') 
     240      write(num_ritz,905)   time,isvol(1),inp(1),isvolf(1),inf(1),volfm(1),    &   ! 6 
     241           isvol(2),inp(2),isvolf(2),inf(2),volfm(2), & ! 5 
     242           isvol(3),inp(3),isvolf(3),inf(3),volfm(3), & ! 5 
     243           isvol(4),inp(4),isvolf(4),inf(4),volfm(4), & ! 5 
     244           isvol(5),inp(5),isvolf(5),inf(5),volfm(5), & ! 5 
     245           isvol(6),inp(6),isvolf(6),inf(6),volfm(6), & ! 5 
     246           isvol(7),inp(7),isvolf(7),inf(7),volfm(7), & ! 5 
     247           diff_H, water_bilan, sum(calv_dtt(:,:))/dtt, sum(ablbord_dtt(:,:))/dtt, sum(Bm_dtt(:,:))/dtt, sum(bmelt_dtt(:,:))/dtt !6 
     248 
     249 
     250 
     251!!905 format(f0.2,7(1x,2(e10.5,i6,1x),e10.5), 6(1x,e15.8)) 
     252905 format(f0.2,7(1x,e15.8,1x,e15.8,1x,e15.8,1x,e15.8,1x,e15.8), 6(1x,e15.8)) 
     253 
     254      if (ncshortout.eq.1) then  ! ecriture netcdf 
     255         ! Total ice mass 
     256         var_shortoutput(1)=lim 
     257         ! Mass above floatation 
     258         var_shortoutput(2)=volf 
     259         ! Grounded ice area 
     260         var_shortoutput(3)=iareag 
     261         ! Floating ice area 
     262         var_shortoutput(4)=iareaf 
     263         ! Total SMB flux  
     264         var_shortoutput(5)=tendacabf 
     265         ! Total Basal mass balance flux 
     266         var_shortoutput(6)=tendlibmassbf 
     267         ! Total Basal mass balance flux beneath floating ice 
     268         var_shortoutput(7)=tendlibmassbffl 
     269         ! Total calving flux  
     270         var_shortoutput(8)=tendlicalvf 
     271         ! Total calving and ice front melting flux 
     272         var_shortoutput(9)=tendlifmassbf 
     273         ! Total grounding line flux 
     274         var_shortoutput(10)=tendligroundf 
     275 
     276         nbtimeout=nbtimeout+1 
     277          
     278         status=nf90_put_var(ncid, timeVarID, time, start=(/nbtimeout/)) 
     279         if (status /= nf90_noerr) call handle_err(status) 
     280          
     281         do k=1,nvar  ! boucle sur le nbr de variable a ecrire 
     282            status=nf90_put_var(ncid, varID(k), var_shortoutput(k),start=(/nbtimeout/)) 
     283            if (status /= nf90_noerr) call handle_err(status) 
     284         enddo 
     285         status=nf90_sync(ncid) 
     286         if (status /= nf90_noerr) call handle_err(status) 
     287      endif 
     288 
     289 
    356290end subroutine shortoutput 
    357291 
     292subroutine handle_err(status) 
     293  integer, intent(in) :: status 
     294  if (status /= nf90_noerr) then 
     295     print*,trim(nf90_strerror(status)) 
     296     stop "stopped" 
     297  end if 
     298end subroutine handle_err 
    358299 
    359300end module output_laure16_mod 
Note: See TracChangeset for help on using the changeset viewer.