module m_ECOSM_riverloads
      use mod_xc
      contains

!**********************************************************************

      subroutine ECOSM_readnutrloads()
      implicit none

      call ECOSM_nitrloads()
      call ECOSM_phosloads()
      call ECOSM_sililoads()
!AS      call ECOSM_donoloads()

      end subroutine ECOSM_readnutrloads

!**********************************************************************

      subroutine ECOSM_rivloads(n)
      use mod_necessary_ecovars
c      use mod_ECOSM_bioini
      implicit none
      include 'common_blocks.h'
      include 'biovar.h'
      integer,intent(in) :: n

!AS      print*, 'maxnit', maxval(rivnit)
!AS      print*, 'maxpho', maxval(rivpho)
!AS      print*, 'maxsil', maxval(rivsil)
!AS      print*, 'maxdon', maxval(rivdon)

c        print*,'test REDF ',REDF

!AS divide by dp-layer depth when using TRIP-type nutrients
! distribute over the upper 5 layers
      do i=1-margin,ii+margin
      do j=1-margin,jj+margin
      do k=1,5
      if (sum(dp(i,j,1:5,n))>0) then
         
         bio(i,j,k,n,init)=bio(i,j,k,n,init)+
c     &         (rivnit(i,j)*REDF(18)*REDF(1)*REDF(6))
     &         (rivnit(i,j)*REDF(18)*REDF(1)*REDF(6))
     &          *biodt_sou*onem/sum(dp(i,j,1:5,n))

         bio(i,j,k,n,ipho)=bio(i,j,k,n,ipho)+
     &               (rivpho(i,j)*REDF(19)*REDF(2)*REDF(6))
     &               *biodt_sou*onem/sum(dp(i,j,1:5,n))

         bio(i,j,k,n,isil)=bio(i,j,k,n,isil)+
     &               (rivsil(i,j)*REDF(20)*REDF(3)*REDF(6))
     &               *biodt_sou*onem/sum(dp(i,j,1:5,n))

!AS  take out because it causes detp to be negative
!AS         bio(i,j,k,n,idet)=bio(i,j,k,n,idet)+
!AS     &               rivdon(i,j)*biodt_sou*onem/sum(dp(i,j,1:5,n))
      end if
      end do
      end do
      end do

      end subroutine ECOSM_rivloads

!**********************************************************************

      subroutine ECOSM_nitrloads()
       
      use mod_necessary_ecovars
      use mod_xc
      use mod_za
      use mod_year_info, only : year_info, rt

      implicit none
      include 'common_blocks.h'
      
      logical lexist
      character(len=20)::filename 
      character(len=4) :: char8
      character(len=24) :: txt
      character preambl(5)*79
      real :: xmin,xmax
      real :: xmin2,xmax2
      integer :: month
      integer :: nop1,nop2,nop3
      integer :: nfields,ifld
      integer :: k
      
!--- Reading nitrate forcing field
      nop1=710
      nop2=2710

      filename='forcing.rivnitr'
      inquire(file=trim(filename)//'.a',exist=lexist) 
      if(.not.lexist) then
         write(lp,*)
         write(lp,*) '***** no nitrate riverload *****'
         write(lp,*)
      else
        if(mnproc==1) then
           write(lp,'(a,i4)')'read nitrate riverloads'
        endif
      endif
      call zaiopf(trim(filename)//'.a','old',nop1)
      open(nop2,file=trim(filename)//'.b',status='old',action='read')
      read(nop2,'(a79)') preambl
      call preambl_print(preambl)
c
c --- read desired month
c
      month=rt%imm
      if (mnproc.eq.1) print*, 'month', rt
      if (mnproc.eq.1) print*, 'month', month
      do ifld=1,month
        read(nop2,'(a24,i2.2,2e16.8)')  
     &                  txt,k,xmin2, xmax2
        call zaiord(rivnit(1-nbdy,1-nbdy),ip,.false.,
     &                  xmin,xmax,nop1)
          if     (mnproc.eq.1) then
            write(lp,'(/ a / a,1p3e14.6 / a,1p3e14.6 /)')
     &           '.a and .b files max  min:',
     &           '.a,.b min = ',xmin,xmin2 ,xmin2-xmin ,
     &           '.a,.b max = ',xmax,xmax2 ,xmax2-xmax
            print *,'Month :',ifld
          endif
        ! Check xmin/xmax
        if (abs(xmin-xmin2)>abs(xmin)*1e-4 .or.
     &      abs(xmax-xmax2)>abs(xmax)*1e-4) then
          if     (mnproc.eq.1) then
            write(lp,'(/ a / a,1p3e14.6 / a,1p3e14.6 /)')
     &           'error - .a and .b files not consistent:',
     &           '.a,.b min = ',xmin,xmin2 ,xmin2-xmin ,
     &           '.a,.b max = ',xmax,xmax2 ,xmax2-xmax
            print *,'Month :',ifld
          endif
          call xcstop('(ECOSM_nitrloads)')
          stop '(ECOSM_nitrloads)'
        end if
      enddo
      call zaiocl(nop1)
      close(nop2)

      end subroutine ECOSM_nitrloads

!**********************************************************************

      subroutine ECOSM_phosloads()

      use mod_necessary_ecovars
      use mod_xc
      use mod_za
      use mod_year_info, only : year_info, rt

      implicit none
      include 'common_blocks.h'

      logical lexist
      character(len=20)::filename 
      character(len=4) :: char8
      character(len=24) :: txt
      character preambl(5)*79
      real :: xmin,xmax
      real :: xmin2,xmax2
      integer :: month
      integer :: nop1,nop2,nop3
      integer :: nfields,ifld
      integer :: k
      integer :: i,j
      
!--- Reading phosphate forcing field
      nop1=710
      nop2=2710

      filename='forcing.rivphos'
      inquire(file=trim(filename)//'.a',exist=lexist) 
      if(.not.lexist) then
         write(lp,*)
         write(lp,*) '***** no phosphate riverload *****'
         write(lp,*)
      else
         if(mnproc==1) then
            write(lp,'(a,i4)')'read phosphate riverloads'
         endif
      endif
      call zaiopf(trim(filename)//'.a','old',nop1)
      open(nop2,file=trim(filename)//'.b',status='old',action='read')
      read(nop2,'(a79)') preambl
      call preambl_print(preambl)
c
c --- read desired month
c
      month=rt%imm
      if (mnproc.eq.1) print*, 'month', month
      do ifld=1,month
        read(nop2,'(a24,i2.2,2e16.8)')  
     &                  txt,k,xmin2, xmax2
        call zaiord(rivpho(1-nbdy,1-nbdy),ip,.false.,
     &                  xmin,xmax,nop1)
        ! Check xmin/xmax
        if (abs(xmin-xmin2)>abs(xmin)*1e-4 .or.
     &      abs(xmax-xmax2)>abs(xmax)*1e-4) then
          if     (mnproc.eq.1) then
            write(lp,'(/ a / a,1p3e14.6 / a,1p3e14.6 /)')
     &           'error - .a and .b files not consistent:',
     &           '.a,.b min = ',xmin,xmin2 ,xmin2-xmin ,
     &           '.a,.b max = ',xmax,xmax2 ,xmax2-xmax
            print *,'Month :',ifld
          endif
          call xcstop('(ECOSM_phosloads)')
          stop '(ECOSM_phosloads)'
        end if
      enddo
      call zaiocl(nop1)
      close(nop2)

      end subroutine ECOSM_phosloads

!**********************************************************************
c      
      subroutine ECOSM_sililoads()
     
      use mod_necessary_ecovars
      use mod_xc
      use mod_za
      use mod_year_info, only : year_info, rt

      implicit none
      include 'common_blocks.h'
      
      logical lexist
      character(len=20)::filename 
      character(len=4) :: char8
      character(len=24) :: txt
      character preambl(5)*79
      real :: xmin,xmax
      real :: xmin2,xmax2
      integer :: month
      integer :: nop1,nop2,nop3
      integer :: nfields,ifld
      integer :: k
      integer :: i,j
      
!--- Reading silicate forcing field
      nop1=710
      nop2=2710

      filename='forcing.rivsili'
      inquire(file=trim(filename)//'.a',exist=lexist) 
      if(.not.lexist) then
         write(lp,*)
         write(lp,*) '***** no silicate riverload *****'
         write(lp,*)
      else
         if(mnproc==1) then
            write(lp,'(a,i4)')'read silicate riverloads'
         endif
      endif
      call zaiopf(trim(filename)//'.a','old',nop1)
      open(nop2,file=trim(filename)//'.b',status='old',action='read')
      read(nop2,'(a79)') preambl
      call preambl_print(preambl)
c
c --- read desired month
c
      month=rt%imm
      if (mnproc.eq.1) print*, 'month', month
      do ifld=1,month
        read(nop2,'(a24,i2.2,2e16.8)')  
     &                  txt,k,xmin2, xmax2
        call zaiord(rivsil(1-nbdy,1-nbdy),ip,.false.,
     &                  xmin,xmax,nop1)
        ! Check xmin/xmax
        if (abs(xmin-xmin2)>abs(xmin)*1e-4 .or.
     &      abs(xmax-xmax2)>abs(xmax)*1e-4) then
          if     (mnproc.eq.1) then
            write(lp,'(/ a / a,1p3e14.6 / a,1p3e14.6 /)')
     &           'error - .a and .b files not consistent:',
     &           '.a,.b min = ',xmin,xmin2 ,xmin2-xmin ,
     &           '.a,.b max = ',xmax,xmax2 ,xmax2-xmax
            print *,'Month :',ifld
          endif
          call xcstop('(ECOSM_sililoads)')
          stop '(ECOSM_sililoads)'
        end if
      enddo
      call zaiocl(nop1)
      close(nop2)
 
      end subroutine ECOSM_sililoads

!**********************************************************************
cc      
c      subroutine ECOSM_donoloads()
c     
c      use mod_necessary_ecovars
c      use mod_xc
c      use mod_za
c      use mod_year_info, only : year_info, rt
c
c      implicit none
c      include 'common_blocks.h'
c      
c      logical lexist
c      character(len=20)::filename 
c      character(len=4) :: char8
c      character(len=24) :: txt
c      character preambl(5)*79
c      real :: xmin,xmax
c      real :: xmin2,xmax2
c      integer :: month
c      integer :: nop1,nop2,nop3
c      integer :: nfields,ifld
c      integer :: k
c      integer :: i,j
c      
c!--- Reading DON forcing field
c      nop1=710
c      nop2=2710
c
c      filename='forcing.rivdono'
c      inquire(file=trim(filename)//'.a',exist=lexist) 
c      if(.not.lexist) then
c         write(lp,*)
c         write(lp,*) '***** no DON riverload *****'
c         write(lp,*)
c      else
c         if(mnproc==1) then
c            write(lp,'(a,i4)')'read DON riverloads'
c         endif
c      endif
c      call zaiopf(trim(filename)//'.a','old',nop1)
c      open(nop2,file=trim(filename)//'.b',status='old',action='read')
c      read(nop2,'(a79)') preambl
c      call preambl_print(preambl)
cc
cc --- read desired month
cc
c      month=rt%imm
c      if (mnproc.eq.1) print*, 'month', month
c      do ifld=1,month
c        read(nop2,'(a24,i2.2,2e16.8)')  
c     &                  txt,k,xmin2, xmax2
c        call zaiord(rivdon(1-nbdy,1-nbdy),ip,.false.,
c     &                  xmin,xmax,nop1)
c        ! Check xmin/xmax
c        if (abs(xmin-xmin2)>abs(xmin)*1e-4 .or.
c     &      abs(xmax-xmax2)>abs(xmax)*1e-4) then
c          if     (mnproc.eq.1) then
c            write(lp,'(/ a / a,1p3e14.6 / a,1p3e14.6 /)')
c     &           'error - .a and .b files not consistent:',
c     &           '.a,.b min = ',xmin,xmin2 ,xmin2-xmin ,
c     &           '.a,.b max = ',xmax,xmax2 ,xmax2-xmax
c            print *,'Month :',ifld
c          endif
c          call xcstop('(ECOSM_donoloads)')
c          stop '(ECOSM_donoloads)'
c        end if
c      enddo
c      call zaiocl(nop1)
c      close(nop2)
c 
c      end subroutine ECOSM_donoloads

!**********************************************************************

      end module m_ECOSM_riverloads