  SUBROUTINE cas_bc_current

!&E---------------------------------------------------------------------
!&E                 ***  ROUTINE cas_advcst  ***
!&E
!&E ** Purpose : definition of a test case with constant zonal and baroclinic current
!&E
!&E ** Description :
!&E
!&E ** Called by : cas_init
!&E
!&E ** External calls :
!&E
!&E ** Used ij-arrays :
!&E
!&E ** Modified variables : hx, hy, h0, sal, temp, dy, dhx, dhy, dsigu, dsigw
!&E                         surf, fisac, gwac, dgac, dfiac, fu, fv, sig, sigw
!&E                         tfiusr, tfivsr, goulet, govlet, ig, id, jb, jh
!&E                         ssh,u,v,uz,vz
!&E
!&E ** Reference :
!&E
!&E ** History :
!&E       !  2006     (F. Vandermeirsch) module -Original code
!&E       !  2008-12  (V. Garnier) shift of sigma values : -1 < sig < 0
!&E       !  2009-05  (F. Vandermeirsch) characteristics method with new time scheme
!&E       !  2011-11  (V. Garnier) Spatial extension of ssh(liminm2:limaxp2,ljminm2:ljmaxp2)
!&E
!&E---------------------------------------------------------------------





























































!! * Modules used
   USE parameters
   USE comvars2d
   USE comvarp2d, ONLY : ssh,u,v,sal,temp
   USE comvarp3d, ONLY : uz,vz,ect
   USE comvars3d
   USE output, ONLY : fileresu,tdebsor,tfinsor,nb_out,tinter,tsort,pasor,   &
                      l_out_precdate,l_out_demer,l_out_ssh,l_out_ubt,l_out_u3d, &
                      l_out_sal,l_out_temp,l_out_wind,l_out_buoy,           &
                      l_out_visc,l_out_turb,l_out_off,l_out_atmsflx,nbatm,  &
                      iextmin,iextmax,l_out_nc4par,l_out_pack,              &
                      jextmin,jextmax,kextmin,kextmax,stepi,stepj,stepk,    &
                      igextmin,igextmax,jgextmin,jgextmax,l_filebydate,     &
                      l_out_zaxis_w
   USE obcogcm, ONLY : obc_coefrel
   USE cometdyn
   USE comdiag, ONLY : l_diag,file_diag,tdebdiag,tfindiag
   USE cometemp
   USE comfluxdo

   USE comaree
   USE comcfo
   USE obccombine, ONLY : l_obc_south,l_obc_north,l_obc_west,l_obc_east
   USE obccombine, ONLY : obc_ue,obc_uw,obc_us,obc_un,obc_uze,obc_uzw,obc_uzs,obc_uzn
   USE obccombine, ONLY : obc_ve,obc_vw,obc_vs,obc_vn,obc_vze,obc_vzw,obc_vzs,obc_vzn
   USE obccombine, ONLY : obc_sshe,obc_sshw,obc_sshs,obc_sshn
   USE obcogcm,    ONLY : obc_tw,obc_se,obc_te,obc_ss,obc_ts,obc_sn,obc_tn,obc_sw,obc_tsref
   USE obccombine, ONLY : l_obc_ogcm_rssh,l_obc_ogcm_rt,l_obc_ogcm_rs,l_obc_ogcm_ruvz
   USE obccombine, ONLY : l_obc_mars,l_obc_diri,l_obc_char,l_obc_cycl, &
                          l_obc_tide,l_obc_ogcm
   USE toolmpi, ONLY: ADD_TAG,MASTER,min_all_mpi_rsh,max_all_mpi_rsh

   IMPLICIT NONE

!! * Arguments

!! * Local declarations
   CHARACTER(LEN=19) :: date_start,date_end,date_start3d,date_startagrif
   CHARACTER(LEN=19) :: date_startout,date_endout
   CHARACTER(LEN=19) :: date_startdiag,date_enddiag
   CHARACTER(LEN=19) :: tool_sectodat
   CHARACTER(LEN=lchain) :: logfilename
   INTEGER           :: i,j,k
   REAL(KIND=rsh),DIMENSION(imin:imax,jmin:jmax) :: hx_t,hy_t
   REAL(KIND=rsh)    :: d,lambda,k0,omega0,Ro,amplitude,rt,radi,x,y,dtemp,alpha,beta
   REAL(KIND=rsh)    :: uinit,vinit,hcst,ubar,vbar,usur,ubot
   REAL(KIND=rsh)    :: rhoL,rhoH,rho0,L,eta_cst,dxval,dyval,fuval
   REAL(KIND=rlg)    :: tool_datosec
   REAL(KIND=rsh)    :: value

!!----------------------------------------------------------------------
!! * Executable part

!------------------------------------------------
!------------------------------------------------
!!     courant constant zonal et barocline     !!
!!     surface libre lineaire                  !!
!!     densite constant                        !!
!------------------------------------------------
!------------------------------------------------

! parametres de forme
!---------------------

   fisac=0.0_riolg
   gwac=0.0_riolg
   dfiac=0.01_riolg
   dgac=0.01_riolg
   dxval=20000.0_rsh
   dyval=20000.0_rsh
   dy(:,:)=dyval
   dx(:,:)=dxval
   dx_u(:,:)=dxval
   dx_v(:,:)=dxval
   dx_f(:,:)=dxval
   dy_u(:,:)=dyval
   dy_v(:,:)=dyval
   dy_f(:,:)=dyval
   surf(:,:)=dyval**2
   h0fond=0.01_rsh

   fuval=0.0001_rsh
   fu(:,:)=0.0001_rsh
   fv(:,:)=0.0001_rsh
   radi=57.29578_rsh
!fu(:,:)=0.0001458_rsh*sin(45.0_rsh/radi)
!fv(:,:)=fu(:,:)

   tfiusr(:,:)=0.0_rsh
   tfivsr(:,:)=0.0_rsh
   goulet(:,:)=1.0_rsh
   govlet(:,:)=1.0_rsh
   jb(:)=jmin+1
   jh(:)=jmax
   ig(:)=imin+1
   id(:)=imax
   dhx(:,:)=0.0_rsh
   dhy(:,:)=0.0_rsh
   dhx_g(:,:)=0.0_rsh
   dhy_g(:,:)=0.0_rsh

   DO k=1,kmax
     sig(k)=(REAL(k-kmax,rsh)-0.5_rsh)/REAL(kmax,rsh)
   END DO
   sig(kmax+1)=-sig(kmax)
   sig(0)=ABS(sig(1))-2.0_rsh

   DO k=1,kmax
     dsigw(k)=sig(k+1)-sig(k)
     sigw(k)=0.5_rsh*(sig(k)+sig(k+1))
   END DO
   sigw(kmax+1)=-sigw(kmax)
   sigw(0)=-1.0_rsh

   DO k=1,kmax
     dsigu(k)=sigw(k)-sigw(k-1)
   END DO



! longitude latitude
! ------------------

   lon2d_g(:,0)=(/ (dxval *(i-1), i=0,imax) /)
   DO j=1,jmax
     lon2d_g(:,j)=lon2d_g(:,0)
   END DO
   lon2du_g(:,:)=lon2d_g(:,:)+dxval/2.0_rsh
   lon2dv_g(:,:)=lon2d_g(:,:)
   lat2d_g(0,:)=(/ (dyval *(j-1), j=0,jmax) /)
   DO i=1,imax
     lat2d_g(i,:)=lat2d_g(0,:)
   END DO
   lat2du_g(:,:)=lat2d_g(:,:)
   lat2dv_g(:,:)=lat2d_g(:,:)+dyval/2.0_rsh

   lon2du(liminm1:limax,ljmin:ljmaxp1)=lon2du_g(liminm1:limax,ljmin:ljmaxp1)
   lon2d(limin:limax,ljmin:ljmax)=lon2d_g(limin:limax,ljmin:ljmax)
   lon2dv(limin:limax,ljmin:ljmax)=lon2dv_g(limin:limax,ljmin:ljmax)

   lat2dv(limin:limax,ljminm1:ljmax)=lat2dv_g(limin:limax,ljminm1:ljmax)
   lat2d(limin:limax,ljmin:ljmax)=lat2d_g(limin:limax,ljmin:ljmax)
   lat2du(limin:limax,ljmin:ljmax)=lat2du_g(limin:limax,ljmin:ljmax)

! dates
! -----

   date_start=   '01/01/1900 00:00:00'
   date_end=   '15/01/1900 00:00:00'
!date_end=   '31/01/1900 00:00:00'
!date_end=   '25/02/1900 00:00:00'
!date_end=   '30/06/1900 12:00:00'
!date_end=   '31/12/1900 00:00:00'
   date_start3d=  '01/01/1900 00:00:00'
   date_startagrif='01/01/2002 00:00:00'
   date_ref = '01-Jan-1900 00:00:00'
   tdeb     = tool_datosec(date_start) 
   t=tdeb
   tfin     = tool_datosec(date_end)
   tdeb3d   = tool_datosec(date_start3d)
   tdebagrif = tool_datosec(date_startagrif)
   cdate = tool_sectodat(t)
   date_startshift= 0.0_rlg
   dtmin=200.0_rlg
   dtmax=200.0_rlg
   dt=200.0_rlg
!*******************************************************
   l_modele2d=.TRUE.
   l_modele2d=.FALSE.
!********************************************************
   cflcrt=0.6_rsh
   tobs=12.0_rlg*3600.0_rlg
   tdebobs=t
   cflumaxmax=0.0_rsh
   cflvmaxmax=0.0_rsh

   l_initfromfile=.FALSE.
   l_saverestart_1file=.FALSE.
   l_saverestart_bydate=.FALSE.
   file_saverestart='save.nc'
   saverestart_step=5.0_rlg*3600.0_rlg   ! in hours
   name_in_dt='dt'
   name_bathy_h0='h0'
   name_bathy_hx='hx'
   name_bathy_hy='hy'
   name_in_xev='xev'
   name_in_u='u'
   name_in_v='v'
   name_in_uz='uz'
   name_in_vz='vz'
   name_in_wz='wz'
   name_in_temp='temp'
   name_in_sal='sal'
   name_in_ect='ect'
   name_in_psi='psi'
   name_in_kz='kz'
   name_in_nz='nz'
   name_in_vis2d='name_in_vis2d'
   name_in_vis3d_xe='name_in_vis3d_xe'
   name_in_vis3d_phi='name_in_vis3d_phi'

! outputs
! -------

   nb_out=1
   ALLOCATE( fileresu(nb_out),tdebsor(nb_out),tfinsor(nb_out),tinter(nb_out) )
   ALLOCATE( tsort(nb_out),pasor(nb_out),l_filebydate(nb_out),l_out_precdate(nb_out) )
   ALLOCATE( l_out_off(nb_out),l_out_atmsflx(nb_out),l_out_visc(nb_out),l_out_turb(nb_out) )
   ALLOCATE( l_out_buoy(nb_out), l_out_zaxis_w(nb_out) )
   ALLOCATE( l_out_demer(nb_out),l_out_ubt(nb_out),l_out_u3d(nb_out) )
   ALLOCATE( l_out_sal(nb_out),l_out_temp(nb_out),l_out_ssh(nb_out),l_out_wind(nb_out) )
   ALLOCATE( iextmin(nb_out),iextmax(nb_out),jextmin(nb_out),jextmax(nb_out) )
   ALLOCATE( igextmin(nb_out),igextmax(nb_out),jgextmin(nb_out),jgextmax(nb_out) )
   ALLOCATE( kextmin(nb_out),kextmax(nb_out),stepi(nb_out),stepj(nb_out),stepk(nb_out) )
   fileresu(1)='results.nc'
   suffix='courant'
   date_startout='01/01/1900 00:00:00'
   date_endout='01/01/2010 12:00:00'
   l_maxu=.TRUE.
   tdebsor(1)=tool_datosec(date_startout)
   tfinsor(1)=tool_datosec(date_endout)
!   pasor(1)=450.0_rlg*192.0_rlg
   pasor(1)= 3600.0_rlg ! 24*3600
   l_filebydate(1)=.FALSE.
   IF ( .NOT. l_filebydate(1) ) THEN
     CALL ADD_TAG(fileresu(1),len(fileresu(1)))
   ENDIF
   l_out_off(1)=.FALSE.
   l_out_precdate(1)=.FALSE.
   l_out_demer(1)=.FALSE.
   igextmin(1)=imin
   igextmax(1)=imax
   jgextmin(1)=jmin
   jgextmax(1)=jmax
   iextmin(1)=limin
   iextmax(1)=limax
   jextmin(1)=ljmin
   jextmax(1)=ljmax
   kextmin(1)=1
   kextmax(1)=kmax
   stepi(1)=1
   stepj(1)=1
   stepk(1)=1
   l_out_ssh(1)=.TRUE.
   l_out_ubt(1)=.TRUE.
   l_out_u3d(1)=.TRUE.
   l_out_sal(1)=.TRUE.
   l_out_temp(1)=.TRUE.
   l_out_wind(1)=.FALSE.
   l_out_atmsflx(1)=.FALSE.
   l_out_visc(1)=.FALSE.
   l_out_turb(1)=.FALSE.
   l_out_buoy(1)=.FALSE.
   l_out_zaxis_w(1)=.FALSE.
   l_obc2drank_save=.FALSE.
   l_incbotstress=.FALSE.
   botstressmax=100.0_rsh
   l_stability=.FALSE.
   iscreen=6
   iscreenlog=6
   ierrorlog=17
   iwarnlog=16
   l_out_nc4par = .false.
   l_out_pack = .false.

   IF (MASTER) THEN
   logfilename= 'simu_'//TRIM(suffix)//'.log'
   IF (iscreenlog/=6) OPEN(unit=iscreenlog,file=logfilename,status="replace")
   logfilename='error_'//TRIM(suffix)//'.log'
   IF (ierrorlog/=6) OPEN(unit=ierrorlog,file=logfilename,status="replace")
   logfilename='warning_'//TRIM(suffix)//'.log'
   IF (iwarnlog/=6) OPEN(unit=iwarnlog,file=logfilename,status="replace")
   ENDIF

! diagnostics
!------------

   l_diag=.FALSE.
   date_startdiag = '01/01/1900 00:00:00'
   date_enddiag = '02/01/2001 00:00:00'
   tdebdiag=tool_datosec(date_startdiag)
   tfindiag=tool_datosec(date_enddiag)
   file_diag='diag0.txt'

! meteo
!-------

   l_meteo_stat=.TRUE.
   l_meteo_hom=.TRUE.
   wind_veloc=0.0_rsh
   wind_dir=0.0_rsh
   l_sflx_solarcst=.TRUE.
   sflx_solarcst=0.0_rsh
   nbatm=0

! turbulence
!------------

   fvisc=0.0_rsh
   vis2d(:,:)      =fvisc*0.01_rsh*dyval**1.15_rsh
   vish_xe(:,:,:) = fvisc*0.01_rsh*dyval**1.15_rsh
   vish_phi(:,:,:)= fvisc*0.01_rsh*dyval**1.15_rsh

   turb_nbeq=0
   turb_0eq_option=1
   z0b(:,:)  =1.0e-13
   z0s(:,:)   =0.0_rsh
   hminfrot=0.5_rsh
   hmaxfrot=5.0_rsh
   nzinit=0.0_rsh
   kzinit=0.0_rsh
   frot_lin=0.0_rsh
   nz(:,:,:) = nzinit
   kz(:,:,:) = kzinit
   ect(:,:,:)=0.0_rsh
   l_smagor=.FALSE.

   vismin=0.0_rsh
   vismax=0.0_rsh
   cosmag=0.2_rsh
   icon=2
   kx=0.0_rsh
   ky=0.0_rsh
   qmax=12.0_rsh
   qmaxz=10.0_rsh
   rhoref=1024.0_rsh
   saliref=31.0_rsh
   tetaref=10.0_rsh
   l_stateeq_lin=.TRUE.

   WRITE(*,*)'######################################################'
   PRINT*,'Vous devez compiler avec la clef -Dkey_equeta_temp' 
   WRITE(*,*)'######################################################'
   STOP

   saliref_lin=saliref
   tetaref_lin=tetaref
   rhoair=1.25_rsh
   chp=3986.0_rsh
   coext=0.1_rsh
   file_river='e'

!-----------------------------------------------------

   WRITE(*,*)'####################################'
   WRITE(*,*)'   CAS TEST n 12 courant barocline'
   WRITE(*,*)'####################################'

! bathymetrie
!-------------

   hcst = 4000.0_rsh
   hx_g(:,:)=hcst
   hy_g(:,:)=hcst
   h0_g(:,:)=hcst
   h0(:,:)=hcst

!----------------------------
! vitesses et surface libre:
! solution exterieure +
! initialisation
!----------------------------

   ubot = 0.0_rsh ! vitesse au fond
   usur = 1.0_rsh ! vitesse en surface

   IF (limin==imin .AND. ljmin==jmin) WRITE(*,*)'      k ,       uz(k,1,1),      dsigu(k),       sig(k),       sigw(k)'  
   DO k=1,kmax
! expression analytique de u(z) :
!  u(z)= (usur-ubot)/H * (z+H)
     uzext(k,:,:)=((ubot-usur)*k + (usur-kmax*ubot))/REAL(1-kmax)
     IF (limin==imin .AND. ljmin==jmin) WRITE(*,*)k,uzext(k,1,1),dsigu(k),sig(k),sigw(k)
   END DO
   vzext(:,:,:)= 0.0_rsh

! calcul de la vitesse barotrope
   ubar = 0.0_rsh
   vbar = 0.0_rsh
   DO k=1,kmax
     ubar = ubar + uzext(k,1,1)*dsigu(k)
   END DO

   WRITE(*,*)'ubar modele     ',ubar 
! ubar est la moyenne de usur et ubot
   ubar = (usur+ubot)/2  
   WRITE(*,*)'ubar analytique',ubar 

   WRITE(*,*)'icon = ',icon 

! initialisation
   uinit = ubar  !-0.1_rsh
   vinit = 0.0_rsh
   eta_cst = 0.0_rsh
   xext(:,:)= eta_cst
!------------------------------------------------------------------------------
! Condition aux limites
! devrait etre obc_ogcm_ssh ou obc_tide_ssh
   obc_sshe(:)=xext(imax,:)
   obc_sshw(:)=xext(imin+1,:)
   obc_sshn(:)=xext(:,jmax)
   obc_sshs(:)=xext(:,jmin+1)

! meme remarque
   obc_ue(:,:) = ubar
   obc_uw(:,:) = ubar
   obc_us(:,:) = ubar
   obc_un(:,:) = ubar

   DO k=1,kmax
     obc_uze(k,imax-obc_width+1:imax,jmin:jmax)    = uzext(k,imax-obc_width+1:imax,jmin:jmax) 
     obc_uzw(k,imin:imin+obc_width,jmin:jmax)      = uzext(k,imin:imin+obc_width,jmin:jmax)
     obc_uzs(k,imin:imax,jmin:jmin+obc_width)      = uzext(k,imin:imax,jmin:jmin+obc_width)
     obc_uzn(k,imin:imax,jmax-obc_width+1:jmax)    = uzext(k,imin:imax,jmax-obc_width+1:jmax) 
   END DO

   obc_ve(:,:) = vbar
   obc_vw(:,:) = vbar
   obc_vs(:,:) = vbar
   obc_vn(:,:) = vbar
   obc_vze(:,:,:) = vinit
   obc_vzw(:,:,:) = vinit
   obc_vzs(:,:,:) = vinit
   obc_vzn(:,:,:) = vinit

! initialisation
!----------------------

   ssh(liminm2:limaxp2,liminm2:limaxp2) = xext(liminm2:limaxp2,liminm2:limaxp2)

   u(:,:)    = uinit
   v(:,:)    = vinit 

   uz(:,liminm2:limaxp2,ljminm2:ljmaxp2) = uzext(:,liminm2:limaxp2,ljminm2:ljmaxp2)
   vz(:,liminm2:limaxp2,ljminm2:ljmaxp2) = vzext(:,liminm2:limaxp2,ljminm2:ljmaxp2)

   wz(:,:,:) = 0.0_rsh
   nz(:,:,:) = nzinit
   kz(:,:,:) = kzinit   

   rhoref=1040.0_rsh
   rho0=rhoref
   L = (jmax-1)*dyval 
   rhoL = rhoref*L*fuval/((hcst+eta_cst)*grav)*(usur-ubot)/(sig(kmax)-sig(1))
   rhoH = 6.0_rsh
   PRINT*,'rhoL = ',rhoL 

   DO j=ljminm2,ljmaxp2
! expression analytique de pa(y):
!  pa(:,:) = -f0.rhoref.ubar.y
     pa(:,j) = (-fuval*rhoref*((sig(kmax)+1.0_rsh)*ubot-(sig(1)+1.0_rsh)*usur)/(sig(kmax)-sig(1)) &
                -(hcst+eta_cst)*grav*rhoL/L)*(j-jmin)*dyval
!      pa(:,j) = -fuval*rhoref*usur*(j-jmin)*dyval
   END DO

! open boundaries
! ---------------

! geometry (where are the open boundaries ?)

! cas 1 : 4 frontieres ouvertes
   l_obc_south=.TRUE.
   l_obc_north=.TRUE.
   l_obc_west=.TRUE.
   l_obc_east=.TRUE.

! cas 2 : 2 frontieres ouvertes
!   l_obc_south=.false.
!   l_obc_north=.false.
!   l_obc_west=.true.
!   l_obc_east=.true.

! used scheme to manage obc

   l_obc_mars = .FALSE.
   l_obc_diri = .TRUE.
   l_obc_char = .FALSE.
   l_obc_cycl = .FALSE.
   l_obc_cycl_x=.FALSE.
   l_obc_cycl_y=.FALSE.

! physics introduced into the simulated area (tide or ogcm solution)

   l_obc_tide=.FALSE.
   l_tide_harmcompo=.FALSE.
   l_obc_ogcm=.FALSE.

   obc_coefrel=0.0_rsh

   IF (l_obc_south .AND. l_obc_north .AND. l_obc_west .AND. l_obc_east) THEN
     WRITE(*,*)' *******************************************'
     WRITE(*,*) 'ATTENTION configuration sans mur'
     WRITE(*,*)' *******************************************'
   END IF

   IF (.NOT. l_obc_south) THEN
     WRITE(*,*)' ***************************************'
     WRITE(*,*) 'ATTENTION configuration avec un mur sud'
     WRITE(*,*)' ***************************************'
     hy_g(imin+2:imax,jmin+1)=-valmanq-h0fond
   END IF

   IF (.NOT. l_obc_north) THEN
     WRITE(*,*)' ****************************************'
     WRITE(*,*) 'ATTENTION configuration avec un mur nord'
     WRITE(*,*)' ******************************************'
     hy_g(imin+2:imax,jmax-1)=-valmanq-h0fond
   END IF

   IF (.NOT. l_obc_west) THEN
     WRITE(*,*)' *****************************************'
     WRITE(*,*) 'ATTENTION configuration avec un mur ouest'
     WRITE(*,*)' *****************************************'
     hx_g(imin+1,:)=-valmanq-h0fond
   END IF

   IF (.NOT. l_obc_east) THEN
     WRITE(*,*)' ***************************************'
     WRITE(*,*) 'ATTENTION configuration avec un mur est'
     WRITE(*,*)' ***************************************'
     hx_g(imax-1,:)=-valmanq-h0fond
   END IF
   hx(:,:) = hx_g(liminm2:limaxp2,ljminm2:ljmaxp2)
   hy(:,:) = hy_g(liminm2:limaxp2,ljminm2:ljmaxp2)

!-----------------------
! traceurs :
! solution exterieure +
! initialisation
!----------------------

   CALL obc_tsref

   tetaref=rhoref    

   DO k=1,kmax
     DO j=jmin,jmax    
! expression analytique de rho(z):
! rho(y,z)= rhoL/L.Y - rhoH/H.Z +  rho0
       obc_tw(k,:,j) = -(rhoL/REAL(jmax-1)*(j-jmin-1) &
                         +rho0-rhoH/(sig(kmax)+1.0_rsh)*(sig(k)+1.0_rsh)-rhoref)*grav/rhoref
       obc_te(k,:,j) = -(rhoL/REAL(jmax-1)*(j-jmin-1) &
                         +rho0-rhoH/(sig(kmax)+1.0_rsh)*(sig(k)+1.0_rsh)-rhoref)*grav/rhoref
       temp(k,:,j) = -(rhoL/REAL(jmax-1)*(j-jmin-1) &
                       +rho0-rhoH/(sig(kmax)+1.0_rsh)*(sig(k)+1.0_rsh)-rhoref)*grav/rhoref
     END DO
   END DO

   sal(:,:,:)  = saliref

   DO k=1,kmax
     DO j=jmax-obc_width+1,jmax
       obc_tn(k,:,j) = -(rhoL/REAL(jmax-1)*(j-jmin-1) &
                         +rho0-rhoH/(sig(kmax)+1.0_rsh)*(sig(k)+1.0_rsh)-rhoref)*grav/rhoref
     END DO
     DO j=jmin,jmin+obc_width
       obc_ts(k,:,j) = -(rhoL/REAL(jmax-1)*(j-jmin-1) &
                         +rho0-rhoH/(sig(kmax)+1.0_rsh)*(sig(k)+1.0_rsh)-rhoref)*grav/rhoref
     END DO
   END DO

   obc_se(:,:,:) = saliref
   obc_sw(:,:,:) = saliref
   obc_ss(:,:,:) = saliref
   obc_sn(:,:,:) = saliref

! Percent of wet surface on each cell
! -----------------------------------

   hx_t(:,:)=100000.0_rsh
   hy_t(:,:)=100000.0_rsh
   WHERE (hx_g(:,:) /= -valmanq -h0fond) hx_t(:,:)=hx_g(:,:)
   WHERE (hy_g(:,:) /= -valmanq -h0fond) hy_t(:,:)=hy_g(:,:)
   dhxy(:,:)=1.0_rsh
   DO j=ljminp1,ljmax
      DO i=liminp1,limax
         dhxy(i,j)=MAX(hx(i,j),hx(i-1,j),hy(i,j),hy(i,j-1))- &
                   MIN(hx_t(i,j),hx_t(i-1,j),hy_t(i,j),hy_t(i,j-1))+epsilon
      END DO
   END DO

   fwet(:,:) = 1.0_rsh
   fwetp(:,:) = 1.0_rsh
   fwetmin=1.0_rsh
   IF (fwetmin/=1.0_rsh) THEN
     DO j=ljmin,ljmax
       DO i=limin,limax
         IF (dhxy(i,j) == 0.0_rsh) THEN
            fwet(i,j)=1.0_rsh
           fwetp(i,j)=1.0_rsh
         ELSE
            fwet(i,j)=MIN(MAX(fwetmin,(ssh(i,j)+h0(i,j))/dhxy(i,j)),1.0_rsh)
           fwetp(i,j)=MIN(MAX(fwetmin,(ssh(i,j)+h0(i,j))/dhxy(i,j)),1.0_rsh)
         ENDIF
       END DO
     END DO
   END IF

! initialisation des noms de variables du fichier sortie
!-------------------------------------------------------

   name_out_h0='h0'
   name_out_hx='hx'
   name_out_hy='hy'
   name_out_u='u'
   name_out_v='v'
   name_out_uz='uz'
   name_out_vz='vz'
   name_out_temp='temp'
   name_out_sal='sal'
   name_out_xe='xe'
   name_out_vis2d='VIS2D'
   name_out_vishxe='VISH_XE'
   name_out_vishphi='VISH_PHI'
   name_out_nz='NZ'
   name_out_kz='KZ'
   name_out_bz='BZ'

! initialization of attributes (max and min valid values for each variable)
! -------------------------------------------------------------------------

   riog_valid_min_h0=-50.0_riosh
   riog_valid_max_h0=10000.0_riosh
   riog_valid_min_u=-100.0_riosh
   riog_valid_max_u=100.0_riosh
   riog_valid_min_temp=0.0_riosh
   riog_valid_max_temp=50.0_riosh
   riog_valid_min_sal=0.0_riosh
   riog_valid_max_sal=50.0_riosh
   riog_valid_min_xe=-50.0_riosh
   riog_valid_max_xe=50.0_riosh
   riog_valid_min_wstress=0.0_riosh
   riog_valid_max_wstress=2000.0_riosh
   riog_valid_min_visc=0.0_riosh
   riog_valid_max_visc=100000.0_riosh
   riog_valid_min_nz=0.0_riosh
   riog_valid_max_nz=1000.0_riosh
   riog_valid_min_kz=0.0_riosh
   riog_valid_max_kz=1000.0_riosh
   riog_valid_min_bz=-50.0_riosh
   riog_valid_max_bz=50.0_riosh
   value=MINVAL(dx(:,:))
   CALL MIN_ALL_MPI_RSH(value)
   riog_valid_min_dx=value-1.0_riosh
   value=MAXVAL(dx(:,:))
   CALL MAX_ALL_MPI_RSH(value)
   riog_valid_max_dx=value+1.0_riosh
   value=MINVAL(dy(:,:))
   CALL MIN_ALL_MPI_RSH(value)
   riog_valid_min_dy=value-1.0_riosh
   value=MAXVAL(dy(:,:))
   CALL MAX_ALL_MPI_RSH(value)
   riog_valid_max_dy=value+1.0_riosh

! initialization of other attributes
! ----------------------------------

   cg_conv_title="name of area"
   cg_conv_comment=" "
   cg_conv_source="MARS CONF"
   cg_conv_area="area"
   cg_conv_institution="IFREMER"
   cg_conv_institution_references="http://www.ifremer.fr"
   cg_conv_data_centre="DYNECO/PHYSED"
   cg_conv_data_centre_references="http://www.ifremer.physed.fr"
   cg_conv_forecast_range=" "
   cg_conv_forecast_type=" "
   cg_conv_field_type="n/a"
   cg_conv_operational_status=" "

  END SUBROUTINE cas_bc_current
