!*******************************************************************************
!$RCSfile$
!$Revision$
!$Date$
!*******************************************************************************
SUBROUTINE CountChemKindMC( ID,nkind )

USE chem_fi
USE error_fi

IMPLICIT NONE

INTEGER, INTENT( IN  ) :: ID
INTEGER, INTENT( OUT ) :: nkind

INTEGER nmc, ndos, ndep, i

nmc = chemMC(ID)%nFast + chemMC(ID)%nSlow + chemMC(ID)%nParticle + chemMC(ID)%nEquilibrium
ndos = 0
ndep = 0

DO i = 1,nmc
  IF( chemMC(ID)%species(i)%ldos )ndos = ndos + 1
  IF( chemMC(ID)%species(i)%ldep )ndep = ndep + 2 ! Dry and Wet dep
END DO

nkind = nmc   !CONCENTRATION plots
IF( chemMC(ID)%nOutGroup > 0 )THEN
  nkind = nkind + ndep + 2*chemMC(ID)%nOutGroup  !Group dep output + individual species
ELSE
  IF( ndep < nmc )nkind = nkind + ndep           !Otherwise, point to concentration kinds
END IF
IF( chemMC(ID)%nOutGroup > 0 )THEN
  nkind = nkind + ndos + chemMC(ID)%nOutGroup  !Group dose output + individual species
ELSE
  IF( ndos < nmc )nkind = nkind + ndos         !Otherwise, point to concentration kinds
END IF

RETURN
END

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

SUBROUTINE AddChemKindMC( ID,ikind,KindString,is,nkind )

USE chem_fi
USE error_fi

IMPLICIT NONE

INTEGER,                        INTENT( IN    ) :: ID
INTEGER,                        INTENT( INOUT ) :: ikind
CHARACTER(64), DIMENSION(*),    INTENT( INOUT ) :: KindString
INTEGER,       DIMENSION(4),    INTENT( OUT   ) :: is
INTEGER,       DIMENSION(4),    INTENT( OUT   ) :: nkind

INTEGER nmc, ndos, ndep, i

nmc = chemMC(ID)%nFast + chemMC(ID)%nSlow + chemMC(ID)%nParticle + chemMC(ID)%nEquilibrium

!---- Concentration plots

is(1)    = ikind + 1
nkind(1) = nmc

!---- Integrated Concentration plots

is(2)    = 0
nkind(2) = 0

!---- Add component names to KindString and check for surface plots

ndos = 0
ndep = 0

DO i = 1,nmc

  ikind = ikind + 1
  KindString(ikind) = TRIM(chemMC(ID)%species(i)%name)

  IF( chemMC(ID)%species(i)%ldos )ndos = ndos + 1
  IF( chemMC(ID)%species(i)%ldep )ndep = ndep + 1

END DO

!---- Deposition plots

nkind(3) = 0
IF( ndep == 0 )THEN
  is(3) = 0
ELSE IF( ndep == nmc .AND. chemMC(ID)%nOutGroup == 0 )THEN
  is(3)    = is(1)
  nkind(3) = nmc
ELSE
  is(3) = ikind + 1
  DO i = 1,nmc
    IF( chemMC(ID)%species(i)%ldep )THEN
      nkind(3) = nkind(3) + 1
      ikind    = ikind + 1
      KindString(ikind) = TRIM(chemMC(ID)%species(i)%Name)//" (DRY)"
      nkind(3) = nkind(3) + 1
      ikind    = ikind + 1
      KindString(ikind) = TRIM(chemMC(ID)%species(i)%Name)//" (WET)"
    END IF
  END DO
END IF

IF( chemMC(ID)%nOutGroup > 0 )THEN
  IF( ndep == 0 )is(3) = ikind + 1
  DO i = 1,chemMC(ID)%nOutGroup
    nkind(3) = nkind(3) + 1
    ikind    = ikind + 1
    KindString(ikind) = TRIM(chemMC(ID)%OutGroup(i)%Name)//" (DRY)"
    nkind(3) = nkind(3) + 1
    ikind    = ikind + 1
    KindString(ikind) = TRIM(chemMC(ID)%OutGroup(i)%Name)//" (WET)"
  END DO
END IF

!---- Dosage plots

nkind(4) = 0
IF( ndos == 0 )THEN
  is(4) = 0
ELSE IF( ndos == nmc .AND. chemMC(ID)%nOutGroup == 0 )THEN
  is(4)    = is(1)
  nkind(4) = nmc
ELSE
  is(4) = ikind + 1
  DO i = 1,nmc
    IF( chemMC(ID)%species(i)%ldos )THEN
      nkind(4) = nkind(4) + 1
      ikind    = ikind + 1
      KindString(ikind) = TRIM(chemMC(ID)%species(i)%name)
    END IF
  END DO
END IF

IF( chemMC(ID)%nOutGroup > 0 )THEN
  IF( ndos == 0 )is(4) = ikind + 1
  DO i = 1,chemMC(ID)%nOutGroup
    nkind(4) = nkind(4) + 1
    ikind    = ikind + 1
    KindString(ikind) = TRIM(chemMC(ID)%OutGroup(i)%Name)
  END DO
END IF

RETURN
END

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

SUBROUTINE InitChemSliceSurfMC( ID,icomp,grdI )

USE chem_fi
USE error_fi
USE sagstr_fd
USE scipuff_fi
USE PtrGrdStrItf

IMPLICIT NONE

INTEGER, INTENT( IN ) :: ID    ! Multicomponent ID
INTEGER, INTENT( IN ) :: grdI  ! index of slice SAG grid structure
INTEGER, INTENT( IN ) :: icomp ! Component number

INTEGER :: i, ix, iy
REAL    :: x, y, z

TYPE( SAGgrid_str ), POINTER :: srf

INTERFACE
  SUBROUTINE SetChemAmbient( xIn,yIn,zIn,t,leq,lmet,lplot )
    REAL,    INTENT( IN ) :: xIn, yIn, zIn, t
    LOGICAL, INTENT( IN ) :: leq
    LOGICAL, INTENT( IN ) :: lmet
    LOGICAL, OPTIONAL, INTENT( IN ) :: lplot
  END SUBROUTINE
END INTERFACE

!----- set SAG grid structure

srf => SAG_PtrGrdStr( grdI )         ! Associate "local" grid structure pointer
IF( .NOT.ASSOCIATED( srf ))THEN
  nError = UK_ERROR
  eRoutine = 'InitChemSliceSurfMC'
  GOTO 9999
END IF

!--- Initialize top level grid to ambient value

chem  => chemMC(ID)
z     = 0.

!--- Initialize top level grid to ambient value
IF( chem%lAmbFile )THEN
  DO ix = 1,srf%nx
    x = srf%xmin + (FLOAT(ix-1)+0.5)*srf%dx
    DO iy = 1,srf%ny
      i = (iy-1)*srf%nx + ix
      y = srf%ymin + (FLOAT(iy-1)+0.5)*srf%dy
      IF( chem%lAddAmb )THEN
        CALL SetChemAmbient( x,y,z,t,.FALSE.,.FALSE.,.TRUE. )
      ELSE
        chem%species(icomp)%amb = 0.0
      ENDIF
      IF( nError /= NO_ERROR )GOTO 9999
      srf%ipdat(i) = chem%species(icomp)%amb
    END DO
  END DO
ELSE IF( chem%lStepAmb )THEN  !use stepped ambient for single point
  srf%ipdat(:) = chem%Ambient%stepAmb(1,1,icomp)
ELSE
  IF( chem%lAddAmb )THEN
    DO i = 1,srf%nx*srf%ny
      srf%ipdat(i) = chem%species(icomp)%ambient
    END DO
  ELSE
    DO i = 1,srf%nx*srf%ny
      srf%ipdat(i) = 0.0
    END DO
  ENDIF
END IF

9999 CONTINUE

RETURN
END

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

SUBROUTINE InitChemSliceHorizMC( ID,icomp,grdI,slice )

USE met_fi
USE chem_fi
USE error_fi
USE slice_fd
USE sagstr_fd
USE scipuff_fi
USE PtrGrdStrItf

IMPLICIT NONE

INTEGER,           INTENT( IN ) :: ID    ! Multicomponent ID
INTEGER,           INTENT( IN ) :: grdI  ! index of slice SAG grid structure
INTEGER,           INTENT( IN ) :: icomp ! Component number
TYPE( slice_str ), INTENT( IN ) :: slice ! slice definition structure

INTEGER :: i, ix, iy
REAL    :: x, y, z
REAL    :: hx, hy, hz

TYPE( SAGgrid_str ), POINTER :: srf

INTERFACE
  SUBROUTINE SetChemAmbient( xIn,yIn,zIn,t,leq,lmet,lplot )
    REAL,    INTENT( IN ) :: xIn, yIn, zIn, t
    LOGICAL, INTENT( IN ) :: leq
    LOGICAL, INTENT( IN ) :: lmet
    LOGICAL, OPTIONAL, INTENT( IN ) :: lplot
  END SUBROUTINE
END INTERFACE

!----- set SAG grid structure

srf => SAG_PtrGrdStr( grdI )         ! Associate "local" grid structure pointer
IF( .NOT.ASSOCIATED( srf ))THEN
  nError = UK_ERROR
  eRoutine = 'InitChemSliceHorizMC'
  GOTO 9999
END IF

chem  => chemMC(ID)
z     = slice%data(SD_HEIGHT) - hmin

!--- Initialize top level grid to ambient value
IF( chem%lAmbFile )THEN
  DO iy = 1,srf%ny
    y = srf%ymin + (FLOAT(iy-1)+0.5)*srf%dy
    DO ix = 1,srf%nx
      i = (iy-1)*srf%nx + ix
      x = srf%xmin + (FLOAT(ix-1)+0.5)*srf%dx
      IF( lter )THEN
        CALL get_topogPrj( x,y,hz,hx,hy )
      ELSE
        hz = 0.0
      END IF
      IF( chem%lAddAmb )THEN
        CALL SetChemAmbient( x,y,z,t,.FALSE.,.FALSE.,.TRUE. )
      ELSE
        chem%species(icomp)%amb = 0.0
      ENDIF
      IF( nError /= NO_ERROR )GOTO 9999
      srf%ipdat(i) = chem%species(icomp)%amb
    END DO
  END DO
ELSE IF( chem%lStepAmb )THEN !use stepped ambient for single point
  srf%ipdat(:) = chem%Ambient%stepAmb(1,1,icomp)
ELSE
  IF( chem%lAddAmb )THEN
    DO i = 1,srf%nx*srf%ny
      srf%ipdat(i) = chem%species(icomp)%ambient
    END DO
  ELSE
    DO i = 1,srf%nx*srf%ny
      srf%ipdat(i) = 0.0
    END DO
  ENDIF
END IF

9999 CONTINUE

RETURN
END

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

SUBROUTINE InitChemSliceVertMC( ID,icomp,grdI,slice )

USE met_fi
USE chem_fi
USE error_fi
USE slice_fd
USE sagstr_fd
USE scipuff_fi
USE PtrGrdStrItf

IMPLICIT NONE

INTEGER,           INTENT( IN ) :: ID    ! Multicomponent ID
INTEGER,           INTENT( IN ) :: grdI  ! index of slice SAG grid structure
INTEGER,           INTENT( IN ) :: icomp ! Component number
TYPE( slice_str ), INTENT( IN ) :: slice ! slice definition structure

INTEGER :: i, ix, iy
REAL    :: x, y, z
REAL    :: xt, dxt, dyt
REAL    :: hx, hy, hz

TYPE( SAGgrid_str ), POINTER :: srf

INTERFACE
  SUBROUTINE SetChemAmbient( xIn,yIn,zIn,t,leq,lmet,lplot )
    REAL,    INTENT( IN ) :: xIn, yIn, zIn, t
    LOGICAL, INTENT( IN ) :: leq
    LOGICAL, INTENT( IN ) :: lmet
    LOGICAL, OPTIONAL, INTENT( IN ) :: lplot
  END SUBROUTINE
END INTERFACE

!----- set SAG grid structure

srf => SAG_PtrGrdStr( grdI )         ! Associate "local" grid structure pointer
IF( .NOT.ASSOCIATED( srf ))THEN
  nError = UK_ERROR
  eRoutine = 'InitChemSliceVertMC'
  GOTO 9999
END IF

chem  => chemMC(ID)
dxt   = slice%xmax - slice%xmin
dyt   = slice%ymax - slice%ymin

!--- Initialize top level grid to ambient value

IF( chem%lAmbFile )THEN
  DO ix = 1,srf%nx
    xt = srf%xmin + (FLOAT(ix-1)+0.5)*srf%dx
    x  = slice%xmin + xt*dxt
    y  = slice%ymin + xt*dyt
    DO iy = 1,srf%ny
      i = (iy-1)*srf%nx + ix
      IF( lter )THEN
        CALL get_topogPrj( x,y,hz,hx,hy )
      ELSE
        hz = 0.0
      END IF
      z = srf%ymin + (FLOAT(iy-1)+0.5)*srf%dy - hmin - hz
      IF( chem%lAddAmb )THEN
        CALL SetChemAmbient( x,y,z,t,.FALSE.,.FALSE.,.TRUE. )
      ELSE
        chem%species(icomp)%amb = 0.0
      ENDIF
      IF( nError /= NO_ERROR )GOTO 9999
      srf%ipdat(i) = chem%species(icomp)%amb
    END DO
  END DO
ELSE IF( chem%lStepAmb )THEN !use stepped ambient for single point
  srf%ipdat(:) = chem%Ambient%stepAmb(1,1,icomp)
ELSE
  IF( chem%lAddAmb )THEN
  DO i = 1,srf%nx*srf%ny
    srf%ipdat(i) = chem%species(icomp)%ambient
  END DO
  ELSE
    DO i = 1,srf%nx*srf%ny
      srf%ipdat(i) = 0.0
    END DO
  END IF
END IF
9999 CONTINUE

RETURN
END

!=======================================================================

SUBROUTINE SetChemAmbDos( imat,ifld0,tfac )

USE scipuff_fi
USE error_fi
USE chem_fi

IMPLICIT NONE

INTEGER, INTENT( IN ) :: imat
INTEGER, INTENT( IN ) :: ifld0
REAL,    INTENT( IN ) :: tfac

INTEGER mcID, ID, nmc, nsd, i, alloc_stat

INTEGER, DIMENSION(:), ALLOCATABLE :: ispecie

mcID = material(imat)%mcID

SELECT CASE( mat_mc%type(mcID) )
  CASE( MC_CHEM )
    ID  = mat_mc%ID(mcID)
    nmc = chemMC(ID)%nFast + chemMC(ID)%nSlow + chemMC(ID)%nParticle + chemMC(ID)%nEquilibrium
    nsd = 0

    DO i = 1,nmc
      IF( chemMC(ID)%species(i)%ldos )THEN
        nsd = nsd + 1
      END IF
    END DO

    IF( nsd > 0 )THEN

      ALLOCATE( ispecie(nsd),STAT=alloc_stat )
      IF( alloc_stat /= 0 )THEN
        nError = UK_ERROR
        eRoutine = 'SetChemAmbDos'
        eMessage = 'Error allocating list of multi-component species for ambient dose contribution'
        GOTO 9999
      END IF
      nsd = 0
      DO i = 1,nmc
        IF( chemMC(ID)%species(i)%ldos )THEN
          nsd = nsd + 1
          ispecie(nsd) = i
        END IF
      END DO

    ELSE

      ALLOCATE( ispecie(1),STAT=alloc_stat )
      IF( alloc_stat /= 0 )THEN
        nError = UK_ERROR
        eRoutine = 'SetChemAmbDos'
        eMessage = 'Error allocating list of multi-component species for ambient dose contribution'
        GOTO 9999
      END IF
      ispecie(1) = 0

    END IF

    IF( chemMC(ID)%nOutGroup > 0 .OR. nsd > 0 )THEN
      CALL AddChemAmbientDoseMC(ID,ifld0,tfac,nsd,ispecie )
    END IF

  CASE DEFAULT
    nError   = UK_ERROR
    eRoutine = 'SetChemAmbDosRes'
    eMessage = 'Multicomponent error'
    WRITE(eInform,'(A,I6)') 'Unknown multicomponent type :',mat_mc%type(mcID)
    GOTO 9999
END SELECT

9999 CONTINUE

IF( ALLOCATED(ispecie) )DEALLOCATE(ispecie,STAT=alloc_stat)

RETURN
END

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

SUBROUTINE AddChemAmbientDoseMC( ID,ifld0,tfac,nsd,ispecie )

!------ Add ambient contribution to top level of dose field

USE chem_fi
USE error_fi
USE sagstr_fd
USE scipuff_fi
USE surface_fi
USE PtrGrdStrItf

IMPLICIT NONE

INTEGER, INTENT( IN ) :: ID      !Multicomponent ID
INTEGER, INTENT( IN ) :: ifld0   !Output field offset
REAL,    INTENT( IN ) :: tfac    !Time factor (s)
INTEGER, INTENT( IN ) :: nsd     !No. of individual specie doses
INTEGER, DIMENSION(*), &
         INTENT( IN ) :: ispecie !Indices of species

INTEGER :: i, ix, iy, i0, ig, j, k, l, n, isp
REAL    :: x, y, z, f, ppm2ugm3

TYPE( SAGgrid_str ), POINTER :: srf

INTERFACE
  SUBROUTINE SetChemAmbient( xIn,yIn,zIn,t,leq,lmet,lplot )
    REAL,    INTENT( IN ) :: xIn, yIn, zIn, t
    LOGICAL, INTENT( IN ) :: leq
    LOGICAL, INTENT( IN ) :: lmet
    LOGICAL, OPTIONAL, INTENT( IN ) :: lplot
  END SUBROUTINE
END INTERFACE

!----- Set SAG grid structure

srf => SAG_PtrGrdStr( srfdos )
IF( .NOT.ASSOCIATED( srf ))THEN
  nError = UK_ERROR
  eRoutine = 'AddChemAmbientDoseMC'
  GOTO 9999
END IF

!----- Loop over top level grid

chem  => chemMC(ID)

ppm2ugm3 = 1013.25*750.06/(62.4*298.15)   ! *MW later (STP: pb mbar to mmHg, T in K, R = 62.4 and mg to ug)

z = z_dosage

IF( chem%lAmbFile )THEN

  DO ix = 1,srf%nx
    x = srf%xmin + (FLOAT(ix-1)+0.5)*srf%dx
    DO iy = 1,srf%ny
      i = (iy-1)*srf%nx + ix
      y = srf%ymin + (FLOAT(iy-1)+0.5)*srf%dy
      CALL SetChemAmbient( x,y,z,t,.FALSE.,.FALSE.,.TRUE. )
      IF( nError /= NO_ERROR )GOTO 9999

      ig = ifld0

      IF( nsd > 0 )THEN
        DO j = 1,nsd
          ig = ig + 1
          i0 = (ig-1)*srf%mxgrd
          isp = ispecie(j)
          ! Assume ambient non particle species concentrations are in ppm. Also, convert to dosage to concentration.
          IF( chem%oUnits == UNIT_UGM3 .AND. chem%species(isp)%class /= ID_SPECIES_PARTICLE )THEN
            srf%ipdat(i0+i) = srf%ipdat(i0+i) + chem%species(isp)%amb*tfac/dt_save*ppm2ugm3*chem%species(isp)%MW
          ELSE
            srf%ipdat(i0+i) = srf%ipdat(i0+i) + chem%species(isp)%amb*tfac/dt_save
          END IF
        END DO
      END IF

      DO k = 1,chem%nOutGroup
        ig = ig + 1
        i0 = (ig-1)*srf%mxgrd
        DO j = 1,chem%OutGroup(k)%nComp
          l = chem%OutGroup(k)%iComp(j)
          f = chem%OutGroup(k)%fComp(j)
          ! Assume ambient non particle species concentrations are in ppm. Also, convert to dosage to concentration.
          IF( chem%oUnits == UNIT_UGM3 .AND. chem%species(l)%class /= ID_SPECIES_PARTICLE )&
            f = f*ppm2ugm3*chem%species(l)%MW
          srf%ipdat(i0+i) = srf%ipdat(i0+i) + f*chem%species(l)%amb*tfac/dt_save
        END DO
      END DO

    END DO
  END DO

ELSE IF( chem%lStepAmb )THEN  !use stepped ambient for single point. Currently disabled

  n  = srf%nx*srf%ny
  ig = ifld0

  IF( nsd > 0 )THEN
    DO j = 1,nsd
      ig = ig + 1
      i0 = (ig-1)*srf%mxgrd
      srf%ipdat(i0+1:i0+n) = srf%ipdat(i0+1:i0+n) + chem%Ambient%stepAmb(1,1,ispecie(j)) * tfac
    END DO
  END IF

  DO k = 1,chem%nOutGroup
    ig = ig + 1
    i0 = (ig-1)*srf%mxgrd
    DO j = 1,chem%OutGroup(k)%nComp
      l = chem%OutGroup(k)%iComp(j)
      f = chem%OutGroup(k)%fComp(j)
      srf%ipdat(i0+1:i0+n) = srf%ipdat(i0+1:i0+n) + f*chem%Ambient%stepAmb(1,1,l) * tfac
    END DO
  END DO

ELSE

  n  = srf%nx*srf%ny
  ig = ifld0

  IF( nsd > 0 )THEN
    DO j = 1,nsd
      ig  = ig + 1
      i0  = (ig-1)*srf%mxgrd
      isp = ispecie(j)
      ! Assume ambient non particle species concentrations are in ppm. Also, convert to dosage to concentration.
      IF( chem%oUnits == UNIT_UGM3 .AND. chem%species(isp)%class /= ID_SPECIES_PARTICLE )THEN
        srf%ipdat(i0+1:i0+n) = srf%ipdat(i0+1:i0+n) + chem%species(isp)%ambient*tfac/dt_save*ppm2ugm3*chem%species(isp)%MW
      ELSE
        srf%ipdat(i0+1:i0+n) = srf%ipdat(i0+1:i0+n) + chem%species(isp)%ambient*tfac/dt_save
      END IF
    END DO
  END IF

  DO k = 1,chem%nOutGroup
    ig = ig + 1
    i0 = (ig-1)*srf%mxgrd
    DO j = 1,chem%OutGroup(k)%nComp
      l = chem%OutGroup(k)%iComp(j)
      f = chem%OutGroup(k)%fComp(j)
      ! Assume ambient non particle species concentrations are in ppm. Also, convert to dosage to concentration.
      IF( chem%oUnits == UNIT_UGM3 .AND. chem%species(l)%class /= ID_SPECIES_PARTICLE )&
        f = f*ppm2ugm3*chem%species(l)%MW
      srf%ipdat(i0+1:i0+n) = srf%ipdat(i0+1:i0+n) + f*chem%species(l)%ambient*tfac/dt_save
    END DO
  END DO

END IF

9999 CONTINUE

RETURN
END


