MODULE module_sf_sstskin

CONTAINS

   SUBROUTINE sst_skin_update(xland,glw,gsw,hfx,qfx,tsk,ust,emiss,  &
                dtw1,sstsk,dt,stbolt,                                &
                ids, ide, jds, jde, kds, kde,                       &
                ims, ime, jms, jme, kms, kme,                       &
                its, ite, jts, jte, kts, kte                       )


   USE module_wrf_error
   IMPLICIT NONE



   INTEGER , INTENT(IN)           :: ids, ide, jds, jde, kds, kde,   &
                                     ims, ime, jms, jme, kms, kme,   &
                                     its, ite, jts, jte, kts, kte


   REAL,     DIMENSION( ims:ime , jms:jme ) , INTENT(IN   ) :: xland, glw, gsw
   REAL,     DIMENSION( ims:ime , jms:jme ) , INTENT(IN   ) :: hfx, qfx, tsk
   REAL,     DIMENSION( ims:ime , jms:jme ) , INTENT(IN   ) :: ust, emiss
   REAL,     DIMENSION( ims:ime , jms:jme ) , INTENT(INOUT  ) :: dtw1         
   REAL,     DIMENSION( ims:ime , jms:jme ) , INTENT(INOUT  ) :: sstsk        
   REAL,                                      INTENT(IN )   ::   DT           
   REAL,                                      INTENT(IN )   ::   STBOLT       


   REAL :: lw, sw, q, qn, zeta, dep, dtw3, skinmax, skinmin
   REAL :: fs, con1, con2, con3, con4, con5, zlan, q2, ts, phi, qn1
   REAL :: usw, qo, swo, us, tb, dtc, dtw, alw, dtwo, delt, f1
   INTEGER :: i, j, k

   INTEGER , PARAMETER :: n=1152
   REAL , PARAMETER :: z1=3.,an=.3,zk=.4,rho=1.2,rhow=1025.,cw=4190.
   REAL , PARAMETER :: g=9.8,znuw=1.e-6,zkw=1.4e-7,sdate=1201.6667
























      skinmax=-9999.
      skinmin=9999.
      do i=its,ite
      do j=jts,jte

      if(xland(i,j).ge.1.5) then
      qo=glw(i,j)-emiss(i,j)*stbolt*(sstsk(i,j)**4)-2.5e6*qfx(i,j)-hfx(i,j)
      swo=gsw(i,j)
      us=MAX(ust(i,j), 0.01)
      tb=tsk(i,j)-273.15
      dtwo=dtw1(i,j)
      delt=dt

      q=qo/(rhow*cw)
      sw=swo/(rhow*cw)


      f1=1.                   -0.27*exp(-2.8*z1)-0.45*exp(-0.07*z1)

      dtc=0.0

      alw=1.e-5*max(tb,1.)
      con4=16.*g*alw*znuw**3/zkw**2
      usw=sqrt(rho/rhow)*us
      con5=con4/usw**4


      q2=max(1./(rhow*cw),-q)
      zlan=6./(1.+(con5*q2)**0.75)**0.333
      dep=zlan*znuw/usw                    
      fs=0.065+11.*dep-(6.6e-5/dep)*(1.-exp(-dep/8.e-4))
      fs=max(fs,0.01)          
      dtc=dep*(q+sw*fs)/zkw            
      dtc=min(dtc,0.)

      dtw=0.0

      alw=1.e-5*max(tb,1.)
      con1=sqrt(5.*z1*g*alw/an)
      con2=zk*g*alw
      qn=q+sw*f1
      usw=sqrt(rho/rhow)*us

      if(dtwo.gt.0..and.qn.lt.0.) then
         qn1=sqrt(dtwo)*usw**2/con1
         qn=max(qn,qn1)
      endif
      zeta=z1*con2*qn/usw**3
      if(zeta.gt.0.) then
         phi=1.+5.*zeta
      else
         phi=1./sqrt(1.-16.*zeta)
      endif
      con3=zk*usw/(z1*phi)

      dtw=(dtwo+(an+1.)/an*(q+sw*f1)*                             &
                          delt/z1)/(1.+(an+1.)*con3*delt)
      dtw=max(0.,dtw)
      dtwo=dtw
      ts = tb + dtw + dtc

      skinmax=amax1(skinmax,ts-tb)
      skinmin=amin1(skinmin,ts-tb)
      sstsk(i,j)=ts+273.15      
      dtw1(i,j)=dtw              
      endif

      end do
      end do


      return

   END SUBROUTINE sst_skin_update


END MODULE module_sf_sstskin