diff --git a/NAMESPACE b/NAMESPACE
index af92a8f..ab1c433 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -1,6 +1,7 @@
 # Generated by roxygen2: do not edit by hand
 
 export("%>%")
+export(Compute_BRF)
 export(Jfunc1)
 export(Jfunc2)
 export(Jfunc3)
@@ -18,9 +19,9 @@ export(apply_sensor_characteristics)
 export(campbell)
 export(check_brown_lop)
 export(check_spectral_sampling)
-export(compute_BRF)
 export(compute_SRF)
 export(compute_albedo)
+export(compute_brf)
 export(compute_fAPAR)
 export(conservative_scattering)
 export(cost_function_RMSE_PROSAIL)
diff --git a/R/compute_BRF.R b/R/compute_BRF.R
index 03cdd69..4fa6b26 100644
--- a/R/compute_BRF.R
+++ b/R/compute_BRF.R
@@ -10,23 +10,30 @@
 #' @param rdot numeric. Hemispherical-directional refl factor viewing direction
 #' @param rsot numeric. Bi-directional reflectance factor
 #' @param tts numeric. Solar zenith angle
-#' @param SpecATM_Sensor list. direct and diffuse radiation for clear conditions
+#' @param spec_atm_sensor list. direct and diffuse radiation for clear conditions
 #' @param skyl numeric. values for skyl
 #' @return BRF numeric. Bidirectional reflectance factor
 #' @export
-compute_BRF  <- function(rdot, rsot, tts, SpecATM_Sensor, skyl = NULL){
-
-  ############################## #
-  ##	direct / diffuse light	##
-  ############################## #
-  Es <- SpecATM_Sensor$Direct_Light
-  Ed <- SpecATM_Sensor$Diffuse_Light
-  rd <- pi/180
+compute_brf <- function(rdot, rsot, tts, spec_atm_sensor, skyl = NULL) {
+  ##############################
+  ## 	direct / diffuse light	##
+  ##############################
+  e_s <- spec_atm_sensor$Direct_Light
+  e_d <- spec_atm_sensor$Diffuse_Light
+  rd <- pi / 180
   # diffuse radiation (Francois et al., 2002)
-  if (is.null(skyl))
-    skyl <- 0.847-1.61*sin((90-tts)*rd) + 1.04*sin((90-tts)*rd)*sin((90-tts)*rd)
-  PARdiro <- (1-skyl)*Es
-  PARdifo <- skyl*Ed
-  BRF <- (rdot*PARdifo+rsot*PARdiro)/(PARdiro+PARdifo)
-  return(data.frame('BRF' = BRF))
+  if (is.null(skyl)) {
+    skyl <- 0.847 - 1.61 * sin((90 - tts) * rd) + 1.04 * sin((90 - tts) * rd) * sin((90 - tts) * rd)
+  }
+  par_dir_o <- (1 - skyl) * e_s
+  par_dif_o <- skyl * e_d
+  brf <- (rdot * par_dif_o + rsot * par_dir_o) / (par_dir_o + par_dif_o)
+  return(data.frame("BRF" = brf))
+}
+
+#' @rdname prosail-deprecated
+#' @export
+Compute_BRF <- function(rdot, rsot, tts, SpecATM_Sensor, skyl = NULL) {
+  .Deprecated("compute_brf")
+  compute_brf(rdot, rsot, tts, SpecATM_Sensor, skyl)
 }
diff --git a/R/generate_LUT_4SAIL.R b/R/generate_LUT_4SAIL.R
index 8b1e140..7156dff 100644
--- a/R/generate_LUT_4SAIL.R
+++ b/R/generate_LUT_4SAIL.R
@@ -18,65 +18,73 @@
 #' @export
 
 generate_LUT_4SAIL <- function(input_prosail, SpecPROSPECT, SpecSOIL, SpecATM,
-                               band_names = NULL, SAILversion ='4SAIL',
-                               brown_lop = NULL){
-
+                               band_names = NULL, SAILversion = "4SAIL",
+                               brown_lop = NULL) {
   nb_samples <- length(input_prosail[[1]])
-  rdot <- rsot <- rsdt <- rddt <-  BRF <- list()
-  Split <- round(nb_samples/10)
+  rdot <- rsot <- rsdt <- rddt <- BRF <- list()
+  Split <- round(nb_samples / 10)
   pb <- progress_bar$new(
     format = "Generate LUT [:bar] :percent in :elapsed",
-    total = 10, clear = FALSE, width= 100)
-  for (i in seq_len(nb_samples)){
-    if (i%%Split==0 & nb_samples>100) pb$tick()
-    rsoil <- input_prosail[i,]$psoil*SpecSOIL$Dry_Soil +
-      (1-input_prosail[i,]$psoil)*SpecSOIL$Wet_Soil
+    total = 10, clear = FALSE, width = 100
+  )
+  for (i in seq_len(nb_samples)) {
+    if (i %% Split == 0 & nb_samples > 100) pb$tick()
+    rsoil <- input_prosail[i, ]$psoil * SpecSOIL$Dry_Soil +
+      (1 - input_prosail[i, ]$psoil) * SpecSOIL$Wet_Soil
     # if 4SAIL
-    if (SAILversion=='4SAIL'){
-      RefSAIL <- PRO4SAIL(Spec_Sensor = SpecPROSPECT,
-                          input_prospect = input_prosail[i,],
-                          TypeLidf = input_prosail[i,]$TypeLidf,
-                          LIDFa = input_prosail[i,]$LIDFa,
-                          LIDFb = input_prosail[i,]$LIDFb,
-                          lai = input_prosail[i,]$lai,
-                          q = input_prosail[i,]$q,
-                          tts = input_prosail[i,]$tts,
-                          tto = input_prosail[i,]$tto,
-                          psi = input_prosail[i,]$psi,
-                          rsoil = rsoil)
-    } else if (SAILversion=='4SAIL2'){
-      RefSAIL <- PRO4SAIL(Spec_Sensor = SpecPROSPECT,
-                          input_prospect = input_prosail[i,],
-                          TypeLidf = input_prosail[i,]$TypeLidf,
-                          LIDFa = input_prosail[i,]$LIDFa,
-                          LIDFb = input_prosail[i,]$LIDFb,
-                          lai = input_prosail[i,]$lai,
-                          q = input_prosail[i,]$q,
-                          tts = input_prosail[i,]$tts,
-                          tto = input_prosail[i,]$tto,
-                          psi = input_prosail[i,]$psi, rsoil = rsoil,
-                          SAILversion = '4SAIL2',
-                          fraction_brown = input_prosail[i,]$fraction_brown,
-                          diss = input_prosail[i,]$diss,
-                          Cv = input_prosail[i,]$Cv,
-                          Zeta = input_prosail[i,]$Zeta, brown_lop = brown_lop)
+    if (SAILversion == "4SAIL") {
+      RefSAIL <- PRO4SAIL(
+        Spec_Sensor = SpecPROSPECT,
+        input_prospect = input_prosail[i, ],
+        TypeLidf = input_prosail[i, ]$TypeLidf,
+        LIDFa = input_prosail[i, ]$LIDFa,
+        LIDFb = input_prosail[i, ]$LIDFb,
+        lai = input_prosail[i, ]$lai,
+        q = input_prosail[i, ]$q,
+        tts = input_prosail[i, ]$tts,
+        tto = input_prosail[i, ]$tto,
+        psi = input_prosail[i, ]$psi,
+        rsoil = rsoil
+      )
+    } else if (SAILversion == "4SAIL2") {
+      RefSAIL <- PRO4SAIL(
+        Spec_Sensor = SpecPROSPECT,
+        input_prospect = input_prosail[i, ],
+        TypeLidf = input_prosail[i, ]$TypeLidf,
+        LIDFa = input_prosail[i, ]$LIDFa,
+        LIDFb = input_prosail[i, ]$LIDFb,
+        lai = input_prosail[i, ]$lai,
+        q = input_prosail[i, ]$q,
+        tts = input_prosail[i, ]$tts,
+        tto = input_prosail[i, ]$tto,
+        psi = input_prosail[i, ]$psi, rsoil = rsoil,
+        SAILversion = "4SAIL2",
+        fraction_brown = input_prosail[i, ]$fraction_brown,
+        diss = input_prosail[i, ]$diss,
+        Cv = input_prosail[i, ]$Cv,
+        Zeta = input_prosail[i, ]$Zeta, brown_lop = brown_lop
+      )
     }
     rdot[[i]] <- RefSAIL$rdot
     rsot[[i]] <- RefSAIL$rsot
     rsdt[[i]] <- RefSAIL$rsdt
     rddt[[i]] <- RefSAIL$rddt
     # Computes BRF based on outputs from PROSAIL and sun position
-    BRF[[i]] <- compute_BRF(rdot = RefSAIL$rdot, rsot = RefSAIL$rsot,
-                            tts = input_prosail$tts[[i]],
-                            SpecATM_Sensor = SpecATM)
+    BRF[[i]] <- compute_brf(
+      rdot = RefSAIL$rdot, rsot = RefSAIL$rsot,
+      tts = input_prosail$tts[[i]],
+      spec_atm_sensor = SpecATM
+    )
   }
-  BRF <- do.call(cbind,BRF)
-  rdot <- do.call(cbind,rdot)
-  rsot <- do.call(cbind,rsot)
-  rsdt <- do.call(cbind,rsdt)
-  rddt <- do.call(cbind,rddt)
+  BRF <- do.call(cbind, BRF)
+  rdot <- do.call(cbind, rdot)
+  rsot <- do.call(cbind, rsot)
+  rsdt <- do.call(cbind, rsdt)
+  rddt <- do.call(cbind, rddt)
   row.names(BRF) <- row.names(rdot) <- row.names(rsot) <- band_names
   row.names(rsdt) <- row.names(rddt) <- band_names
-  return(list('BRF' = BRF, 'rdot' = rdot, 'rsot' = rsot,
-              'rsdt' = rsdt, 'rddt' = rddt))
+  return(list(
+    "BRF" = BRF, "rdot" = rdot, "rsot" = rsot,
+    "rsdt" = rsdt, "rddt" = rddt
+  ))
 }
diff --git a/R/generate_LUT_BRF.R b/R/generate_LUT_BRF.R
index a584d7f..6fc06a7 100644
--- a/R/generate_LUT_BRF.R
+++ b/R/generate_LUT_BRF.R
@@ -61,10 +61,10 @@ generate_LUT_BRF <- function(input_prosail, SpecPROSPECT, SpecSOIL, SpecATM,
                           brown_lop = brown_lop)
     }
     # Computes BRF based on outputs from PROSAIL and sun position
-    BRF[[i]] <- compute_BRF(rdot = RefSAIL$rdot,
+    BRF[[i]] <- compute_brf(rdot = RefSAIL$rdot,
                             rsot = RefSAIL$rsot,
                             tts = input_prosail$tts[[i]],
-                            SpecATM_Sensor = SpecATM)
+                            spec_atm_sensor = SpecATM)
   }
   BRF <- do.call(cbind,BRF)
   row.names(BRF) <- band_names
diff --git a/R/generate_LUT_PROSAIL.R b/R/generate_LUT_PROSAIL.R
index 899f875..527c6ff 100644
--- a/R/generate_LUT_PROSAIL.R
+++ b/R/generate_LUT_PROSAIL.R
@@ -64,10 +64,10 @@ generate_LUT_PROSAIL <- function(input_prosail, SpecPROSPECT, SpecSOIL, SpecATM,
                           brown_lop = brown_lop)
     }
     # Computes BRF based on outputs from PROSAIL and sun position
-    BRF[[i]] <- compute_BRF(rdot = RefSAIL$rdot,
+    BRF[[i]] <- compute_brf(rdot = RefSAIL$rdot,
                             rsot = RefSAIL$rsot,
                             tts = input_prosail$tts[[i]],
-                            SpecATM_Sensor = SpecATM)
+                            spec_atm_sensor = SpecATM)
     fCover[i] <- RefSAIL$fCover
     fAPAR[i] <- compute_fAPAR(abs_dir = RefSAIL$abs_dir,
                               abs_hem = RefSAIL$abs_hem,
diff --git a/R/merit_RMSE_PROSAIL.R b/R/merit_RMSE_PROSAIL.R
index 30b1fe8..22b8634 100644
--- a/R/merit_RMSE_PROSAIL.R
+++ b/R/merit_RMSE_PROSAIL.R
@@ -37,9 +37,9 @@ merit_RMSE_PROSAIL <- function(xinit, parms_xinit, brf_mes, SpecPROSPECT_Sensor,
                   tts = input_prosail$tts, tto = input_prosail$tto,
                   psi = input_prosail$psi, rsoil = rsoil)
   # Computes BRF based on outputs from PROSAIL and sun position
-  brf_mod <- compute_BRF(rdot = Ref$rdot, rsot = Ref$rsot,
+  brf_mod <- compute_brf(rdot = Ref$rdot, rsot = Ref$rsot,
                          tts = input_prosail$tts,
-                         SpecATM_Sensor = SpecATM_Sensor)
+                         spec_atm_sensor = SpecATM_Sensor)
   # compute cost
   fc <- cost_function_RMSE_PROSAIL(brf_mes = brf_mes, brf_mod$BRF, xprior,
                                    prior_info = prior_info)
diff --git a/R/prosail-depcated.R b/R/prosail-depcated.R
new file mode 100644
index 0000000..f2d6482
--- /dev/null
+++ b/R/prosail-depcated.R
@@ -0,0 +1,12 @@
+#' Deprecated functions in prosail
+#'
+#' These functions still work but will be removed (defunct) in the next version.
+#'
+#' \itemize{
+#'  \item \code{\link{Compute_BRF}}: This function is deprecated, and will
+#'  be removed in the next version of this package.
+#'  Use \code{\link{compute_brf}} instead.
+#' }
+#'
+#' @name prosail-deprecated
+NULL
diff --git a/man/prosail-deprecated.Rd b/man/prosail-deprecated.Rd
new file mode 100644
index 0000000..47d763a
--- /dev/null
+++ b/man/prosail-deprecated.Rd
@@ -0,0 +1,19 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/compute_BRF.R, R/prosail-depcated.R
+\name{Compute_BRF}
+\alias{Compute_BRF}
+\alias{prosail-deprecated}
+\title{Deprecated functions in prosail}
+\usage{
+Compute_BRF(rdot, rsot, tts, SpecATM_Sensor, skyl = NULL)
+}
+\description{
+These functions still work but will be removed (defunct) in the next version.
+}
+\details{
+\itemize{
+\item \code{\link{Compute_BRF}}: This function is deprecated, and will
+be removed in the next version of this package.
+Use \code{\link{compute_brf}} instead.
+}
+}
diff --git a/paper/paper.md b/paper/paper.md
index 4630d41..0170c64 100644
--- a/paper/paper.md
+++ b/paper/paper.md
@@ -329,12 +329,12 @@ sun zenith angle and assuming clear sky conditions.
 
 ```r
 # Compute BRF with known skyl
-BRF_4SAIL <- compute_BRF(rdot = Ref_4SAIL$rdot, rsot = Ref_4SAIL$rsot,
-                         skyl = 0.23, SpecATM_Sensor = SpecATM)
+BRF_4SAIL <- compute_brf(rdot = Ref_4SAIL$rdot, rsot = Ref_4SAIL$rsot,
+                         skyl = 0.23, spec_atm_sensor = SpecATM)
 
 # Compute BRF assuming clear sky conditions and using sun zenith angle
-BRF_4SAIL <- compute_BRF(rdot = Ref_4SAIL$rdot, rsot = Ref_4SAIL$rsot,
-                         tts = 40, SpecATM_Sensor = SpecATM)
+BRF_4SAIL <- compute_brf(rdot = Ref_4SAIL$rdot, rsot = Ref_4SAIL$rsot,
+                         tts = 40, spec_atm_sensor = SpecATM)
 
 ```
 
@@ -449,8 +449,8 @@ Refl_1nm <- PRO4SAIL(N = truth$N, CHL = truth$CHL, CAR = truth$CAR,
                      lai = truth$lai, q = parm_set$q, LIDFa = Init$LIDFa, 
                      rsoil = rsoil, tts = parm_set$tts, tto = parm_set$tto, 
                      psi = parm_set$psi)
-brf_1nm <- compute_BRF(rdot = Refl_1nm$rdot, rsot = Refl_1nm$rsot,
-                       tts = parm_set$tts, SpecATM_Sensor = SpecATM)
+brf_1nm <- compute_brf(rdot = Refl_1nm$rdot, rsot = Refl_1nm$rsot,
+                       tts = parm_set$tts, spec_atm_sensor = SpecATM)
 
 # invert PROSAIL on BRF with 1 nm spectral sampling
 est_1nm <- invert_PROSAIL(brf_mes = brf_1nm$BRF, initialization = Init,
diff --git a/tests/testthat/test-PRO4SAIL_iterative_opt.R b/tests/testthat/test-PRO4SAIL_iterative_opt.R
index a1a0c64..5549729 100644
--- a/tests/testthat/test-PRO4SAIL_iterative_opt.R
+++ b/tests/testthat/test-PRO4SAIL_iterative_opt.R
@@ -20,10 +20,10 @@ test_that("PROSAIL iterative optimization", {
                        lai = truth$lai, q = parm_set$q, LIDFa = parm_set$LIDFa,
                        rsoil = rsoil, tts = parm_set$tts, tto = parm_set$tto,
                        psi = parm_set$psi)
-  brf_1nm <- prosail::compute_BRF(rdot = Refl_1nm$rdot,
+  brf_1nm <- prosail::compute_brf(rdot = Refl_1nm$rdot,
                                   rsot = Refl_1nm$rsot,
                                   tts = parm_set$tts,
-                                  SpecATM_Sensor = SpecATM)
+                                  spec_atm_sensor = SpecATM)
   # invert 1 nm data
   est <- invert_PROSAIL(brf_mes = brf_1nm$BRF,
                         initialization = Init,
diff --git a/tests/testthat/test-compute_BRF.R b/tests/testthat/test-compute_BRF.R
index 2f934c2..f55b55a 100644
--- a/tests/testthat/test-compute_BRF.R
+++ b/tests/testthat/test-compute_BRF.R
@@ -1,10 +1,10 @@
 test_that("PRO4AIL and PRO4SAIL2 produce physically possible BRF values", {
   # run PROSAIL with 4SAIL2
   Refl <- PRO4SAIL()
-  BRF_4SAIL <- compute_BRF(rdot = Refl$rdot,
+  BRF_4SAIL <- compute_brf(rdot = Refl$rdot,
                            rsot = Refl$rsot,
                            tts = 30,
-                           SpecATM_Sensor = prosail::SpecATM)
+                           spec_atm_sensor = prosail::SpecATM)
 
   # run PROSAIL with 4SAIL2
   input_prospect <- data.frame('CHL' = c(40, 5), 'CAR' = c(8, 4),
@@ -15,10 +15,10 @@ test_that("PRO4AIL and PRO4SAIL2 produce physically possible BRF values", {
                          TypeLidf = 2, LIDFa = 30, lai = 5, q = 0.1, tts = 30,
                          tto = 10, psi = 90, rsoil = prosail::SpecSOIL$Dry_Soil,
                          fraction_brown = 0.5, diss = 0.5, Cv = 1, Zeta = 1)
-  BRF_4SAIL2 <- compute_BRF(rdot = Ref_4SAIL2$rdot,
+  BRF_4SAIL2 <- compute_brf(rdot = Ref_4SAIL2$rdot,
                             rsot = Ref_4SAIL2$rsot,
                             tts = 30,
-                            SpecATM_Sensor = prosail::SpecATM)
+                            spec_atm_sensor = prosail::SpecATM)
 
   expect_true(all(BRF_4SAIL$BRF >= 0))
   expect_true(all(BRF_4SAIL2$BRF >= 0))
diff --git a/vignettes/prosail2.Rmd b/vignettes/prosail2.Rmd
index ca0a98f..fc6d083 100644
--- a/vignettes/prosail2.Rmd
+++ b/vignettes/prosail2.Rmd
@@ -165,7 +165,7 @@ Ref_4SAIL2 <- PRO4SAIL(SAILversion = '4SAIL2',
 
 ## Compute simplified bidirectional reflectance factor (BRF) under the assumption of clear conditions conditions
 
-The function `compute_BRF` computes the bi-directional reflectance factor in the 
+The function `compute_brf` computes the bi-directional reflectance factor in the 
 direction of the observer, by combining both hemispherical-directional and 
 bi-directional reflectance factors.
 The direct and diffuse light are taken into account as described by 
@@ -177,14 +177,14 @@ by [Spitters et al., 1986](https://www.sciencedirect.com/science/article/pii/016
 # Ref_4SAIL is the variable obtained when running PRO4SAIL as in the previous 
 # illustration.
 # SpecATM corresponds to the direct and diffuse radiation solar spectra
-BRF_4SAIL <- compute_BRF(rdot = Ref_4SAIL$rdot,
+BRF_4SAIL <- compute_brf(rdot = Ref_4SAIL$rdot,
                          rsot = Ref_4SAIL$rsot,
                          tts = tts,
-                         SpecATM_Sensor = SpecATM)
-BRF_4SAIL2 <- compute_BRF(rdot = Ref_4SAIL2$rdot,
+                         spec_atm_sensor = SpecATM)
+BRF_4SAIL2 <- compute_brf(rdot = Ref_4SAIL2$rdot,
                           rsot = Ref_4SAIL2$rsot,
                           tts = tts,
-                          SpecATM_Sensor = SpecATM)
+                          spec_atm_sensor = SpecATM)
 
 ```
 
@@ -198,7 +198,7 @@ absorptance for hemispherical diffuse incident flux.
 fAPAR_4SAIL <- compute_fAPAR(abs_dir = Ref_4SAIL$abs_dir,
                              abs_hem = Ref_4SAIL$abs_hem,
                              tts = tts,
-                             SpecATM_Sensor = SpecATM)
+                             spec_atm_sensor = SpecATM)
 ```
 
 ## Compute albedo
diff --git a/vignettes/prosail5.Rmd b/vignettes/prosail5.Rmd
index 63db4fd..4363e1b 100644
--- a/vignettes/prosail5.Rmd
+++ b/vignettes/prosail5.Rmd
@@ -73,10 +73,10 @@ Ref_Sensor <- PRO4SAIL(Spec_Sensor = Spec_Sensor,
                        rsoil = SpecSOIL_Sensor$Dry_Soil)
 
 # Computes BRF based on outputs from PROSAIL and sun position
-BRF_Sensor <- compute_BRF(rdot = Ref_Sensor$rdot, 
+BRF_Sensor <- compute_brf(rdot = Ref_Sensor$rdot, 
                           rsot = Ref_Sensor$rsot, 
                           tts = tts,
-                          SpecATM_Sensor = SpecATM_Sensor)
+                          spec_atm_sensor = SpecATM_Sensor)
 
 # Computes fAPAR based on outputs from PROSAIL and sun position
 fAPAR_Sensor <- compute_fAPAR(abs_dir = Ref_Sensor$abs_dir,