  SUBROUTINE cas_rip

!&E-----------------------------------------------------------------------------
!&E                 ***  ROUTINE cas_rip  ***
!&E
!&E ** Purpose :
!&E
!&E ** Description :
!&E
!&E ** Called by : cas_init
!&E
!&E ** Reference :
!&E
!&E ** History :
!&E        !  2010-03  (A.-C. Bennis) Creation
!&E        !  2014-03  (V. Garnier) New formalism for the test case initialization
!&E        !  2014-05  (B.Thouvenin) compatibility with substance (first corrections)
!&E
!&E-----------------------------------------------------------------------------





























































!! * Modules used
   USE parameters,ONLY : rsh,rlg,riosh,riolg,imin,imax,jmin,jmax,kmax, &
                         jshift_obc_north,jshift_obc_south,valmanq,pi, &
                         ljmin,ljminp1,ljmax,limin,limax,liminm1,limaxp1,ljminm1,ljmaxp1
   USE comvars2d, ONLY : fisac,gwac,fu,fv,dtmin,dtmax,dt,halfdt,suffix,     &
                         vish_xe,vish_phi,l_maxu,icon,h0fond,hm,     &
                         hx_g,hy_g,h0_g,hx,hy,dhx,dhy,dhx_g,dhy_g,   &
                         iscreen,l_obc_cycl_y,frot_lin,tobs,jb,jh,   &
                         l_spongw,ig,l_modele2d,vismul,sponge_nbcell,&
                         fvisc,vis2d,l_smagor,iscreenlog
   USE comvars2d, ONLY : l_saverestart_bydate,saverestart_step
   USE comvars3d, ONLY : kzinit,nzinit,l_stateeq_lin,bgdiff,         &
                         turb_nbeq,turb_2eq_option
   USE comvarp3d, ONLY : ect
   USE obccombine,ONLY : l_obc_cycl,l_obc_east,l_obc_mars,l_obc_tide,l_obc_diri

   USE output, ONLY : pasor,l_out_ssh,l_out_ubt,l_out_u3d,l_out_sal,l_out_temp, &
                      l_out_wave

   USE toolmpi, ONLY : ADD_TAG,ex_i_rsh,ex_j_rsh,MASTER
   USE casinit,   ONLY : cas_init_param,cas_init_grid,cas_init_bathy, &
                         cas_init_io,cas_init_obc,cas_init_meteo,     &
                         cas_init_visc,cas_init_turb,cas_init_tspar,  &
                         cas_init_progvar,cas_init_misc,cas_init_river

   USE casinit,   ONLY : cas_init_wave


   USE comsubstance, ONLY : nvpc,nvp,nv_cp,nv_dis,nv_fix,nv_dri,nv_int,nv_sorb,     &
                            name_var,long_name_var,standard_name_var,nv_ncp,        &
                            name_filesubs,cv_wat,                        &
                            valid_min_var,valid_max_var,cv_wat,nb_var,date_startsub
   USE casinit,   ONLY : cas_init_subs

   USE comoasis,  ONLY :    init_comoasis


   USE wave,      ONLY : l_wave_rdir,l_wave_rfoc,l_wave_rubot,l_wave_ruubot,l_wave_rvubot, &
                         l_wave_rtp, l_wave_rhs,l_wave_rt0m1,l_wave_rtauwx,l_wave_rtauwy,  &
                         l_wave_rtauox,l_wave_rtauoy,l_wave_rwlm,l_wave_rbhd,file_wave,    &
                         l_wave_rtus,l_wave_russ,l_wave_rfbb,l_wave_rtbb
   USE wave,      ONLY : name_out_wave_tauox,name_out_wave_tauoy


   IMPLICIT NONE

!! * Arguments

!! * Local declarations
   LOGICAL           :: l_nc4par
   CHARACTER(LEN=19) :: date_start,date_end,date_start3d,date_startagrif
   INTEGER           :: i,j,k
   INTEGER           :: ni_wall_w,nb_outfile
   REAL(KIND=rsh)    :: z0bottom,z0surface,depthcst,limit,nzinit0
   REAL(kind=rsh)    :: vismuloc
   REAL(KIND=rsh),DIMENSION(jmin:jmax-jshift_obc_north-jshift_obc_south) :: hx_read
   REAL(KIND=rsh),DIMENSION(imin:imax,jmin:jmax) :: hx_tmp
   LOGICAL           :: l_filesubs

! Grid parameters
   INTEGER           :: imax_ww3
   REAL(KIND=rsh), PARAMETER :: dyval=24.0_rsh,dxval=24.0_rsh   ! mais pourquoi pas 20 m comme la ref ?
!   REAL(KIND=rsh), PARAMETER :: dyval=2.0_rsh,dxval=2.0_rsh
!   REAL(KIND=rsh), PARAMETER :: dyval=24.0_rsh,dxval=8.0_rsh
! Bathymetric parameters
   REAL(kind=rsh), PARAMETER :: a1 = 2.97_rsh, & ! [m]
                                a2 = 1.5_rsh,  & ! [m]
                                b1 = 0.075,    & ! beach slope close to the shore
                                b2 = 0.0064,   & ! slope offshore of the bar
                                g1 = TAN(b1)/TAN(b2), &
                                perturb = 0.1, & ! magnitude of the perturbation and measures the depth of rip channels
                                lambda = 256,  & ! wavelength of the perturbation (or spacing of rip channels)
                                xc = 80.0_rsh    ! alongshore bar location
   REAL(kind=rsh)            :: x,  &            ! distance [m] from the alongshore
                                y,z,minz
   REAL(KIND=rsh),ALLOCATABLE,DIMENSION(:,:) :: bat



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

   !PRINT*, 'enter CAS_RIP routine'

! general parametrization (each parameter can be reassigned)
! -----------------------

! dynamics
   CALL cas_init_param

! over-assignement
   l_maxu=.FALSE.
   h0fond=0.001_rsh
   hm=0.01_rsh

! tracers : among others, saliref=saliref_lin=35.5_rsh,
! tetaref=tetaref_lin=10.0_rsh
   CALL cas_init_tspar

! over-assignement
   icon=0
   l_stateeq_lin=.FALSE.


! Bathymetry  (Yu and Slinn, 2003)
! ----------

! assignement : input variables for cas_init_bathy routine
   depthcst = 0.0_rsh   ! analytical distribution
   ni_wall_w = 3  ! land over ni_wall_w cells (from imin column) at the western part of the domain


! Approximated beach profile measured at Duck, North Carolina
   j=jmin
   DO i=imin,imax
     x=i*dxval - ni_wall_w*dxval + dxval/2.0_rsh  ! Along shore at u location (en fait -1 cm et non 0 si dxval=2)
     h0_g(i,j) = (a1-a1/g1)*TANH(b1*x/a1) + b1*x/g1 - a2*EXP(-5.0_rsh*((x-xc)/xc)**2)
   END DO
   DO j=jmin+1,jmax
     h0_g(:,j)=h0_g(:,jmin)
   END DO
! Sinusoidal longshore perturbation
   DO j=jmin+jshift_obc_south+1,jmax-jshift_obc_north
     y=dyval*(j-jshift_obc_south)
     DO i=imin,imax
       x=i*dxval - ni_wall_w*dxval + dxval/2.0_rsh  ! Along shore at u location (en fait -1 cm et non 0)
       z=perturb * COS(2.0_rsh*pi*y/lambda) *  EXP(-5.0_rsh*((x-xc)/xc)**2)
       h0_g(i,j) = h0_g(i,j)*(1.0_rsh + z)
     END DO
   END DO
! Elimination de certains rips
   DO i=imin,imax
     x=i*dxval - ni_wall_w*dxval + dxval/2.0_rsh  ! Along shore at u location (en fait -1 cm et non 0)
     IF (x <= 300.0_rsh) THEN
       minz = 10000.0_rsh
       DO j=jmin+jshift_obc_south+1,jmax-jshift_obc_north
         minz = MIN(h0_g(i,j),minz)
       END DO
       DO j=jmin+jshift_obc_south+1,jmax-jshift_obc_north
         y=dyval*(j-jshift_obc_south)
         IF (y <= 3.0_rsh*lambda/2.0_rsh) THEN
           h0_g(i,j) = minz
         ELSE IF (y >= 5.0_rsh*lambda/2.0_rsh .AND. y <= 7.0_rsh*lambda/2.0_rsh) THEN
           h0_g(i,j) = minz
         ELSE IF (y >= 9.0_rsh*lambda/2.0_rsh) THEN
           h0_g(i,j) = minz
         END IF
       END DO
     END IF
   END DO
   DO j=jmin,jmin+jshift_obc_south
     h0_g(:,j) = h0_g(:,jshift_obc_south+1)
   END DO
   DO j=MIN0(jmax,jmax-jshift_obc_north+1),jmax
     h0_g(:,j) = h0_g(:,jmax-jshift_obc_north)
   END DO

! Changement de resolution
! ========================
! - 1 - Le fait d enlever minz fait que l on ne retrouve pas exactement les memes valeurs d
! une resolution a l autre
! - 2 - attention a avoir les memes distances pour comparer les veleurs entre 2 resolutions (8m OK)

!!!Faire autrement pour les bancs decouvrants
!!!#ifdef moyenne
   do i=imin,imax-1
   do j=jmin,jmax
   hx_g(i,j)=0.5_rsh*(h0_g(i,j)+h0_g(i+1,j))
   enddo
   enddo
   do i=imin,imax
   do j=jmin,jmax-1
   hy_g(i,j)=0.5_rsh*(h0_g(i,j)+h0_g(i,j+1)) 
   enddo
   enddo
   hy_g(:,jmax)=hy_g(:,jmax-1)
   hx_g(imax,:)=hx_g(imax-1,:)
!!!!#endif



! ajout d une ligne

! ATTENTION. Comment gere t on les bancs decouvrants si on augmente la
! resolution ?  Penser a fct aussi attention siggens

! analytical bathymetry, plus walls
   CALL cas_init_bathy(depthcst,ni_wall_w=ni_wall_w)

! over-assignement
! c est quoi ces dhx superieurs a 0.35 FRANCK ? ATTENTION
   dhx(:,:)=0.5_rsh
   dhy(:,:)=0.5_rsh
   dhx_g(:,:)=0.5_rsh
   dhy_g(:,:)=0.5_rsh


! grid parametrization
! --------------------

! assignement : input variables for cas_init_grid routine

! parametres de forme
!---------------------
   fisac=0.0_riolg
   gwac=0.0_riolg



   CALL cas_init_grid(dxval,dyval)

! over-assignement
   fu(:,:)=0.0_rsh
   fv(:,:)=0.0_rsh


! IO parametrization
! ------------------

! assignement : input variables for cas_init_io routine
   date_start=   '31/12/1900 12:00:00'
   date_end=   '03/01/1901 14:00:00'

   date_end=   '31/12/1900 14:00:00'

!date_end=   '01/02/1901 14:00:00'
!date_end=   '31/12/1900 12:10:00'
   date_start3d=  '31/12/1900 12:00:00'
   date_startagrif='31/12/1900 00:00:00'
   dtmin=2.0_rlg
   dtmax=2.0_rlg
   dt=2.0_rlg
   halfdt=dt*0.5_rlg

   suffix='qeuler'
   l_nc4par = .TRUE.
!   l_nc4par = .FALSE.

   nb_outfile=1

   CALL cas_init_io(date_start,date_end,date_start3d,date_startagrif,suffix,nb_file=nb_outfile,l_nc4par=l_nc4par)

! over-assignement

   pasor(:)=24.0_rlg*3600.0_rlg
   pasor(:)=4.0_rlg*60.0_rlg

   l_out_sal(:)=.FALSE.
   l_out_temp(:)=.FALSE.
   l_out_wave(:)=.TRUE.
   tobs=500.0
!   l_saverestart_bydate=.TRUE.
!   saverestart_step=0.00001
   iscreenlog=12


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

! no open boundary
   CALL cas_init_obc

! over-assignement
! geometry (where are the open boundaries ?)
   l_obc_east=.TRUE.
! used scheme to manage obc
   l_obc_mars = .TRUE.
!   l_obc_diri = .TRUE.
   l_obc_cycl = .TRUE.
   l_obc_cycl_y=.TRUE.
! physics introduced into the simulated area (tide or ogcm solution)
   l_obc_tide=.TRUE.


! viscosity
!-----------

! none
   CALL cas_init_visc

! over-assignement
   vish_xe(:,:,:) = 0.2_rsh
   vish_phi(:,:,:)= 0.2_rsh
   fvisc = 0.5173532   ! to get 0.2
   l_smagor=.FALSE.
   vismul=10.0_rsh
   l_spongw=.FALSE.
   sponge_nbcell=5 ! par defaut
   DO j=ljminp1,ljmax
     IF(l_spongw) THEN                        ! couche eponge a l ouest si ouverte
       IF(h0_g(ig(j),j).GT.-valmanq) THEN
         DO i=ig(j),MAX0(MIN0(ig(j)+sponge_nbcell,imax),imin)
           IF ((i >= limin).AND.(i <= limax)) THEN
!     vismuloc=(5.* exp(-real(i-1)/10.))
!     vismuloc=(i-(ig(j)+sponge_nbcell))*(1-vismul)/sponge_nbcell+1.0_rsh
           vismuloc=vismul
           vis2d(i,j)=vismuloc*vis2d(i,j)
           IF(.NOT.l_modele2d) THEN
             vish_phi(:,i,j)=vismuloc*vish_phi(:,i,j)
             vish_xe(:,i,j) =vismuloc*vish_xe(:,i,j)
           END IF
           ENDIF
         ENDDO
       ENDIF
     ENDIF
   END DO

   CALL ex_i_rsh(-1,0,kmax,liminm1,limax,ljminm1,ljmax,vish_phi(:,liminm1:limax,ljminm1:ljmax))
   CALL ex_j_rsh(-1,0,kmax,liminm1,limax,ljminm1,ljmax,vish_phi(:,liminm1:limax,ljminm1:ljmax))
   CALL ex_i_rsh(0,1,kmax,limin,limaxp1,ljmin,ljmaxp1,vish_xe(:,limin:limaxp1,ljmin:ljmaxp1))
   CALL ex_j_rsh(0,1,kmax,limin,limaxp1,ljmin,ljmaxp1,vish_xe(:,limin:limaxp1,ljmin:ljmaxp1))



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

   z0bottom   =0.005_rsh
   z0surface  =1.6_rsh   ! 1.6*hs
   nzinit0=0.0_rsh
   kzinit=0.0_rsh

! Gaspard parametrization
   CALL cas_init_turb(kzinit,nzinit0,z0surface,z0bottom)

! over-assignement
   turb_nbeq=2
   turb_2eq_option=2
   frot_lin=0.004_rsh
   ect(:,:,:)=1.0e-6
   bgdiff=1.0e-6
! attention ces coefs ont change de base (initturb)       c1turb=1.38_rsh
! attention ces coefs ont change de base (initturb)       c2turb=1.98_rsh
! attention ces coefs ont change de base (initturb)       cmu=0.5544_rsh


! meteorology
! -----------

! no forcing
   CALL cas_init_meteo



! wave
! ----

! no forcing
   CALL cas_init_wave

! over-assignement
   l_wave_rdir=.TRUE.;l_wave_rfoc=.TRUE.;l_wave_rubot=.FALSE.;l_wave_ruubot=.TRUE.;l_wave_rvubot=.TRUE.
   l_wave_rtp=.FALSE.;l_wave_rhs=.TRUE.;l_wave_rt0m1=.TRUE.;l_wave_rtauwx=.FALSE.;l_wave_rtauwy=.FALSE.
   l_wave_rtauox=.TRUE.;l_wave_rtauoy=.TRUE.;l_wave_rwlm=.TRUE.;l_wave_rbhd=.TRUE.
   l_wave_rtus=.TRUE.;l_wave_russ=.TRUE.;l_wave_rfbb=.TRUE.;l_wave_rtbb=.TRUE.
   file_wave="../inputs/ww3.190012.nc"
   file_wave="../inputs/results_mars.nc"
   file_wave="../inputs/ww3tomars.nc"
   name_out_wave_tauox="utwo";name_out_wave_tauoy="vtwo"



! River and outflows
! -------------------

! no river
   CALL cas_init_river


! OASIS coupling
! --------------

   CALL init_comoasis

! over-assignement


! initialization of the pronostic variables
! -----------------------------------------

! set to 0,tetaref,saliref
   CALL cas_init_progvar

! over-assignement


! Miscellaneous
! -------------

! Domain of calculation
! Percent of wet surface on each cell
   CALL cas_init_misc

! over-assignement
   jh(:)=jh(:)-jshift_obc_north
   jb(:)=jb(:)+jshift_obc_south


!fdumas ajout tracer cas test 060112

   nvpc=0; nvp=0; nv_cp=0; nv_ncp=0; nv_sorb=0; nv_fix=0; nv_dri=0; nv_int=0;
   nv_dis=1 
   l_filesubs=.FALSE.  ! TRUE si variable lu dans le fichier variable.dat
   IF (l_filesubs) THEN  ! we use "sub_read_var" instead of "cas_init_subs"
     name_filesubs = '../inputs/variable.dat'
   ELSE
! allocation et initialisation des tableaux cv_wat,f90,t90...
     CALL cas_init_subs(date_start)    ! nvp must be defined before
   ENDIF


! fin fdumas

! faire en fonction du choix de BThouvenin


   !PRINT*, 'exit CAS_RIP routine'

  END SUBROUTINE cas_rip
