# Waste management #

# VS: volatile solid per animal category (kg VS / animal / year) (IPCC 2006)
volatileSolid <- function(VSrate, LW){
  VS <- (VSrate * LW/1000) * 365 # kg VS / animal / year
  
  return(VS)
}

# Wasting management total emission (Mg CH4)
# EF: vector of "fator de emissão para a população de animais definida" (g CH4)
# N: vector of "número de cabeças por espécie animal (por categoria T no país)"
wasteTotalEmission <- function(EF, N){
  CH4emission <- (EF * N) / 1e+3    # Mg CH4 / year

  return(CH4emission)
} 

# B0 (kg CH4 / VS)
# MCF (%)
# MS (%)
# [EF] = VS * kg CH4/VS = kg CH4
emissionFactor <- function(VS, B0, MCF, MS){
  EF <- VS * B0 * 0.67 * MCF/100 * MS  # kg CH4 / year
  
  return(EF)
}

# Excreção anual de N por animais do rebanho #
# Nrate: taxa diária de excreção de N (kg N / kg animal mass / day)
# LW: média dos pesos vivos de cada categoria (kg)
NExcretion <- function(Nrate, LW){
  Nex <- (Nrate * LW / 1000) * 365  # kg N / animal / year
    
  return(Nex)
}

# Emissões diretas de N2O pelo manejo de dejetos #
# Nex : kg N / animal / year
# EF: Fator de emissão direta de N (kg N2O-N / kg N)
# Ncdg: annual nitrogen input via co-digestate (kg N / year)
directN2OEmissions <- function(N, Nex, MS, Ncdg, EF){
  N2Od <- (N * Nex * MS + Ncdg) * EF * (44 / 28)  # kg N2O / year
  N2Od <- N2Od / 1e+3                             # Mg N2O / year
  
  return(N2Od)
}

# Estimativa do N volatilizado com o manejo de dejetos #
# Nex : kg N / animal / year
# N: números de animais por categoria
volatilizationN <- function(N, Nex, MS, Ncdg, Frac){
  Nvol <- (N * Nex * MS + Ncdg) * Frac  # kg N / year
  
  return(Nvol)
}

# Emissões indiretas de N2O pela reposição do N volatilizado do manejo de dejetos # 
# EF: Fator de emissão indireta de N (kg N2O-N / (kg NH3-N + NOx-N volatilised)) (IPCC 2006)
# Nvol: N volatilizado (kg N / year)
indirectN2OEmissions <- function(Nvol, EF){
  N2Oi <- (Nvol * EF) * 44 / 28   # kg N2O / year
  N2Oi <- N2Oi / 1e+3             # Mg N2O / year
  
  return(N2Oi)
}

herdWasteManagement <- function(allInp, allHerdStructure){
  
  resList <- list()
  
  for(hd in c("Equilibrium", "Current")){
    
    allMS  <- allInp$MS[,allInp$climateRegion]
    allMCF <- allInp$MCF[,allInp$climateRegion]
    
    catEmissions <- lapply(names(allInp$Category), function(nm){
        N  <- allHerdStructure[allHerdStructure$Category == nm, hd]
        
        if(N != 0){
          LW <- allInp$averageLW[[nm]]
          manageType <- allInp$Category[[nm]][["wastingManageType"]]
          MS     <- allMS[[manageType]]
          MCF    <- allMCF[[manageType]]
          Ncdg   <- allInp$Ncdg[[manageType]]
          EFN2Od <- allInp$EFN2Od[[manageType]]
          Frac   <- allInp$Frac[[manageType]]

          Nex   <- NExcretion(Nrate = allInp$Nrate, LW = LW)
          VS    <- volatileSolid(VSrate = allInp$VSrate, LW = LW)
          Nvol  <- volatilizationN(N = N, Nex = Nex, MS = MS, Ncdg = Ncdg, Frac = Frac)
          EFCH4 <- emissionFactor(VS = VS, B0 = allInp$B0, MCF = MCF, MS = MS)
          CH4   <- wasteTotalEmission(EF = EFCH4, N = N)
          dN2O  <- directN2OEmissions(N = N, Nex = Nex, MS = MS, Ncdg = Ncdg, EF = EFN2Od)
          iN2O  <- indirectN2OEmissions(Nvol = Nvol, EF = allInp$EFN2Oi)
          
          resDF <- data.frame(Category             = nm,
                              NExcretion           = Nex,
                              volatileSolid        = VS,
                              volatilizedN         = Nvol,
                              emissionFactor       = EFCH4,
                              CH4emission          = CH4,
                              directN2OEmissions   = dN2O,
                              indirectN2OEmissions = iN2O)
        }else{
          resDF <- data.frame(Category             = nm,
                              NExcretion           = 0,
                              volatileSolid        = 0,
                              volatilizedN         = 0,
                              emissionFactor       = 0,
                              CH4emission          = 0,
                              directN2OEmissions   = 0,
                              indirectN2OEmissions = 0)
        }
        
        return(resDF)
      }
    )
    
    resList[[hd]] <- do.call(rbind, catEmissions)
  }

  return(resList)
}