final changes for 0.2.0

IREA-CNR-MI · Jun 8, 2020 · c8acb65 · c8acb65
1 parent 3e70a6d
commit c8acb65
Show file tree

Hide file tree

Showing 57 changed files with 934 additions and 352 deletions.
diff --git a/DESCRIPTION b/DESCRIPTION
@@ -12,7 +12,7 @@ Description: The goal of prismaread is allowing to easily import PRISMA (http://
 License: GPL-3
 Encoding: UTF-8
 LazyData: true
-URL: https://github.com/lbusett/prismaread
+URL: https://lbusett.github.io/prismaread/
 BugReports: https://github.com/lbusett/prismaread/issues
 Imports: 
     hdf5r,
@@ -30,7 +30,8 @@ Imports:
 Suggests: testthat, 
     spelling,
     knitr,
-    rmarkdown
+    rmarkdown, 
+    ggplot2
 Remotes: 
     isciences/exactextractr
 RoxygenNote: 7.1.0
diff --git a/R/convert_prisma.R b/R/convert_prisma.R
@@ -208,7 +208,7 @@ convert_prisma <- function(in_file,
   acqtime <- hdf5r::h5attr(f, "Product_StartTime")
 
   out_file_angles <- paste0(tools::file_path_sans_ext(out_file), "_", source,
-                            "_ANGLES.txt")
+                            ".ang")
   utils::write.table(data.frame(date = acqtime,
                                 sunzen   = sunzen,
                                 sunaz = sunaz, stringsAsFactors = FALSE),
@@ -247,11 +247,18 @@ convert_prisma <- function(in_file,
     }
   }
 
+  if (!is.null(selbands_vnir)){
+    seqbands_vnir <- unlist(lapply(selbands_vnir, FUN = function(x) which.min(abs(wl_vnir - x))))
 
+  } else {
+    seqbands_vnir <- (1:66)[wl_vnir != 0]
+  }
+  wl_vnir   <- wl_vnir[seqbands_vnir]
+  fwhm_vnir <- fwhm_vnir[seqbands_vnir]
   # be sure to remove zeroes also if VNIR already present to avoid  ----
   # errors on creation of FULL
-  fwhm_vnir <- fwhm_vnir[wl_vnir != 0]
-  wl_vnir <- wl_vnir[wl_vnir != 0]
+  # fwhm_vnir <- fwhm_vnir[wl_vnir != 0]
+  # wl_vnir <- wl_vnir[wl_vnir != 0]
 
   # get SWIR data cube and convert to raster ----
 
@@ -285,23 +292,17 @@ convert_prisma <- function(in_file,
                          in_L2_file = in_L2_file)
     }
   }
-
-  # be sure to remove zeroes also if swir already present to avoid  ----
-  # errors on creation of FULL if recycling an existing cube from file
-  fwhm_swir <- fwhm_swir[wl_swir != 0]
-  wl_swir   <- wl_swir[wl_swir != 0]
-
-  if (!is.null(selbands_vnir)){
-    seqbands_vnir <- unlist(lapply(selbands_vnir, FUN = function(x) which.min(abs(wl_vnir - x))))
-    wl_vnir   <- wl_vnir[seqbands_vnir]
-    fwhm_vnir <- fwhm_vnir[seqbands_vnir]
-  }
-
   if (!is.null(selbands_swir)){
     seqbands_swir <- unlist(lapply(selbands_swir, FUN = function(x) which.min(abs(wl_swir - x))))
-    wl_swir   <- wl_swir[seqbands_swir]
-    fwhm_swir <- fwhm_swir[seqbands_swir]
+  } else {
+    seqbands_swir <- (1:173)[wl_swir != 0]
   }
+  wl_swir   <- wl_swir[seqbands_swir]
+  fwhm_swir <- fwhm_swir[seqbands_swir]
+  # be sure to remove zeroes also if swir already present to avoid  ----
+  # errors on creation of FULL if recycling an existing cube from file
+  # fwhm_swir <- fwhm_swir[wl_swir != 0]
+  # wl_swir   <- wl_swir[wl_swir != 0]
 
   # create FULL data cube and convert to raster ----
   out_file_full <- paste0(tools::file_path_sans_ext(out_file), "_", source,
@@ -320,16 +321,19 @@ convert_prisma <- function(in_file,
       if (file.exists(out_file_vnir) && file.exists(out_file_swir)) {
         rast_vnir <- raster::stack(out_file_vnir)
         rast_swir <- raster::stack(out_file_swir)
+
         if (join_priority == "SWIR") {
           rast_tot <- raster::stack(rast_vnir[[which(wl_vnir < min(wl_swir))]], rast_swir)
           wl_tot   <- c(wl_vnir[which(wl_vnir < min(wl_swir))], wl_swir)
           fwhm_tot <- c(fwhm_vnir[which(wl_vnir < min(wl_swir))], fwhm_swir)
-          names(rast_tot) <- paste0("wl_", round(wl_tot, digits = 4))
+          names(rast_tot) <- c(paste0("b", seqbands_vnir[which(wl_vnir < min(wl_swir))], "_v"),
+                               paste0("b", seqbands_swir, "_s"))
         } else {
           rast_tot <- raster::stack(rast_vnir, rast_swir[[which(wl_swir > max(wl_vnir))]])
           wl_tot   <- c(wl_vnir, wl_swir[which(wl_swir > max(wl_vnir))])
           fwhm_tot <- c(fwhm_vnir, fwhm_swir[which(wl_swir > max(wl_vnir))])
-          names(rast_tot) <- paste0("wl_", round(wl_tot, digits = 4))
+          names(rast_tot) <- c(paste0("b", seqbands_vnir, "_v"),
+                               paste0("b", seqbands_swir[which(wl_swir > max(wl_vnir))], "_s"))
         }
 
         rm(rast_vnir)
@@ -341,20 +345,28 @@ convert_prisma <- function(in_file,
         rastwrite_lines(rast_tot, out_file_full, out_format, proc_lev, join = TRUE)
 
         if (out_format == "ENVI") {
+          cat("band names = {", paste(names(rast_tot),collapse=","), "}", "\n",
+              file=raster::extension(out_file_full, "hdr"), append=TRUE)
           out_hdr <- paste0(tools::file_path_sans_ext(out_file_full), ".hdr")
           write(c("wavelength = {",
                   paste(round(wl_tot, digits = 4), collapse = ","), "}"),
                 out_hdr, append = TRUE)
           write(c("fwhm = {",
                   paste(round(fwhm_tot, digits = 4), collapse = ","), "}"),
                 out_hdr, append = TRUE)
+          write("wavelength units = Nanometers")
+          write("sensor type = PRISMA")
         }
-        out_file_txt <- paste0(tools::file_path_sans_ext(out_file_full), "_wavelengths.txt")
-        utils::write.table(data.frame(band = 1:length(wl_tot),
+
+        out_file_txt <- paste0(tools::file_path_sans_ext(out_file_full), ".wvl")
+        utils::write.table(data.frame(band   = 1:length(wl_tot),
+                                      orband = substring(names(rast_tot) , 2),
                                       wl   = wl_tot,
-                                      fwhm = fwhm_tot, stringsAsFactors = FALSE),
+                                      fwhm = fwhm_tot,
+                                      stringsAsFactors = FALSE),
                            file = out_file_txt, row.names = FALSE)
 
+
         rm(rast_tot)
         gc()
       } else {
@@ -395,7 +407,7 @@ convert_prisma <- function(in_file,
                             paste0(out_file_latlon, ".envi"))
 
   if (file.exists(out_file_latlon) & !overwrite) {
-    message("CLD file already exists - use overwrite = TRUE or change output file name to reprocess")
+    message("LATLON file already exists - use overwrite = TRUE or change output file name to reprocess")
   } else {
     if (LATLON) {
       prisma_create_latlon(f,

diff --git a/R/prisma_create_additional.R b/R/prisma_create_additional.R
@@ -13,11 +13,12 @@ prisma_create_additional <- function(f,
                                      out_format,
                                      base_georef,
                                      fill_gaps,
+                                     fix_geo,
                                      in_L2_file = NULL){
 
     message(" - Accessing ", type, " dataset - ")
 
-    # Get geo info ----
+    in_L2_file# Get geo info ----
     geo <- prisma_get_geoloc(f, "1", "HCO", "VNIR", in_L2_file)
 
     if (type == "CLD") {

diff --git a/R/prisma_create_angles.R b/R/prisma_create_angles.R
@@ -80,4 +80,8 @@ prisma_create_angles <- function(f,
     gc()
     message(" - Writing ANGLES raster - ")
     rastwrite_lines(rastang, out_file, out_format)
+    if (out_format == "ENVI") {
+        cat("band names = {", paste(names(rastang),collapse=","), "}", "\n",
+            file=raster::extension(out_file, "hdr"), append=TRUE)
+    }
 }
diff --git a/R/prisma_create_latlon.R b/R/prisma_create_latlon.R
@@ -55,4 +55,8 @@ prisma_create_latlon <- function(f,
     gc()
     message(" - Writing LATLON raster - ")
     rastwrite_lines(rastang, out_file, out_format)
+    if (out_format == "ENVI") {
+        cat("band names = {", paste(names(rastang),collapse=","), "}", "\n",
+            file=raster::extension(out_file, "hdr"), append=TRUE)
+    }
 }
diff --git a/R/prisma_create_pan.R b/R/prisma_create_pan.R
@@ -102,7 +102,6 @@ prisma_create_pan <- function(f,
     rm(pan_cube)
     gc()
 
-
     message("- Writing PAN raster -")
 
     rastwrite_lines(rast_pan, out_file_pan, out_format, proc_lev,

diff --git a/R/prisma_create_swir.R b/R/prisma_create_swir.R
@@ -70,7 +70,8 @@ prisma_create_swir <- function(f,
                     band <- prisma_basegeo(band, lon, lat, fill_gaps)
                 } else {
                     message("Importing Band: ", band_swir, " of: 173")
-                    band <- raster::raster((swir_cube[,order_swir[band_swir], ]))
+                    band <- raster::raster(swir_cube[,order_swir[band_swir], ],
+                                           crs = "+proj=longlat +datum=WGS84")
                     if (proc_lev == "1") {
                         band <- (band / swir_scale) - swir_offset
                     }
@@ -107,47 +108,51 @@ prisma_create_swir <- function(f,
     }
 
     # Write the cube ----
-    wl_swir   <- wl_swir[wl_swir != 0]
-    fwhm_swir <- fwhm_swir[fwhm_swir != 0]
     if (is.null(selbands_swir)) {
-        names(rast_swir) <- paste0("wl_", round(wl_swir, digits = 4))
+        # names(rast_vnir) <- paste0("b", seqbands, "_", round(wl_vnir, digits = 3))
+        orbands <- seqbands[wl_swir != 0]
+        names(rast_swir) <- paste0("b", orbands)
+        wl_sub   <- wl_swir[wl_swir != 0]
+        fwhm_sub <- fwhm_swir[wl_swir != 0]
     } else {
-        names(rast_swir) <- paste0("wl_", round(wl_swir[seqbands], digits = 4))
+        # names(rast_vnir) <- paste0("b", seqbands, "_", round(wl_vnir[seqbands], digits = 3))
+        orbands <- seqbands
+        names(rast_swir) <- paste0("b", orbands)
+        wl_sub   <- wl_swir[seqbands]
+        fwhm_sub <- fwhm_swir[seqbands]
     }
     rm(swir_cube)
     rm(band)
     gc()
     message("- Writing SWIR raster -")
+
     rastwrite_lines(rast_swir,
                     out_file_swir,
                     out_format,
                     proc_lev,
                     scale_min = swir_min,
                     scale_max = swir_max)
 
-    if (!is.null(selbands_swir)) {
-        wl_sub   <- wl_swir[seqbands]
-        fwhm_sub <- fwhm_swir[seqbands]
-    } else {
-        wl_sub   <- wl_swir
-        fwhm_sub <- fwhm_swir
-    }
-
     if (out_format == "ENVI") {
+        cat("band names = {", paste(names(rast_swir),collapse=","), "}", "\n",
+            file=raster::extension(out_file_swir, "hdr"), append=TRUE)
         out_hdr <- paste0(tools::file_path_sans_ext(out_file_swir), ".hdr")
-
         write(c("wavelength = {",
                 paste(round(wl_sub, digits = 4), collapse = ","), "}"),
               out_hdr, append = TRUE)
         write(c("fwhm = {",
                 paste(round(fwhm_sub, digits = 4), collapse = ","), "}"),
               out_hdr, append = TRUE)
+        write("wavelength units = Nanometers")
+        write("sensor type = PRISMA")
     }
 
-    out_file_txt <- paste0(tools::file_path_sans_ext(out_file_swir), "_meta.txt")
+    out_file_txt <- paste0(tools::file_path_sans_ext(out_file_swir), ".wvl")
     utils::write.table(data.frame(band = 1:length(wl_sub),
+                                  orband = orbands,
                                   wl   = wl_sub,
-                                  fwhm = fwhm_sub, stringsAsFactors = FALSE),
+                                  fwhm = fwhm_sub,
+                                  stringsAsFactors = FALSE),
                        file = out_file_txt, row.names = FALSE)
 
 }
diff --git a/R/prisma_create_vnir.R b/R/prisma_create_vnir.R
@@ -61,6 +61,7 @@ prisma_create_vnir <- function(f,
     } else {
         seqbands <- unlist(lapply(selbands_vnir, FUN = function(x) which.min(abs(wl_vnir - x))))
     }
+
     for (band_vnir in seqbands) {
 
         if (wl_vnir[band_vnir] != 0) {
@@ -115,46 +116,39 @@ prisma_create_vnir <- function(f,
 
     # Write the cube ----
 
-    wl_vnir   <- wl_vnir[wl_vnir != 0]
-    fwhm_vnir <- fwhm_vnir[fwhm_vnir != 0]
-    # browser()
     if (is.null(selbands_vnir)) {
-        names(rast_vnir) <- paste0("wl_", round(wl_vnir, digits = 4))
+        # names(rast_vnir) <- paste0("b", seqbands, "_", round(wl_vnir, digits = 3))
+        orbands <- seqbands[wl_vnir != 0]
+        names(rast_vnir) <- paste0("b", orbands)
+        wl_sub   <- wl_vnir[wl_vnir != 0]
+        fwhm_sub <- fwhm_vnir[wl_vnir != 0]
     } else {
-        names(rast_vnir) <- paste0("wl_", round(wl_vnir[seqbands], digits = 4))
+        # names(rast_vnir) <- paste0("b", seqbands, "_", round(wl_vnir[seqbands], digits = 3))
+        orbands <- seqbands
+        names(rast_vnir) <- paste0("b", orbands)
+        wl_sub   <- wl_vnir[seqbands]
+        fwhm_sub <- fwhm_vnir[seqbands]
     }
+
+    # wl_vnir   <- wl_vnir[wl_vnir != 0]
+    # fwhm_vnir <- fwhm_vnir[fwhm_vnir != 0]
     rm(vnir_cube)
     rm(band)
     gc()
 
     message("- Writing VNIR raster -")
 
-    # if(proc_lev %in% c("1", "2B", "2C") && !base_georef) {
-    #     browser()
-    # }
-
     rastwrite_lines(rast_vnir,
                     out_file_vnir,
                     out_format,
                     proc_lev,
                     scale_min = vnir_min,
                     scale_max = vnir_max)
 
-    if (!is.null(selbands_vnir)) {
-        wl_sub   <- wl_vnir[seqbands]
-        fwhm_sub <- fwhm_vnir[seqbands]
-    } else {
-        wl_sub   <- wl_vnir
-        fwhm_sub <- fwhm_vnir
-    }
-
-
     if (out_format == "ENVI") {
-        # browser()
+        cat("band names = {", paste(names(rast_vnir),collapse=","), "}", "\n",
+            file=raster::extension(out_file_vnir, "hdr"), append=TRUE)
         out_hdr <- paste0(tools::file_path_sans_ext(out_file_vnir), ".hdr")
-        write(c("wavelength = {",
-                paste(round(wl_sub, digits = 4), collapse = ","), "}"),
-              out_hdr, append = TRUE)
         write(c("wavelength = {",
                 paste(round(wl_sub, digits = 4), collapse = ","), "}"),
               out_hdr, append = TRUE)
@@ -163,10 +157,23 @@ prisma_create_vnir <- function(f,
               out_hdr, append = TRUE)
         write("wavelength units = Nanometers")
         write("sensor type = PRISMA")
+        # hdr_in  <- readLines(out_hdr)
+        # start_line_bnames <- grep("band names", hdr_in) + 1
+        # end_line_bnames <- start_line_bnames + length(wl_sub) - 1
+        # for (bb in start_line_bnames:(end_line_bnames - 1)) {
+        #     hdr_in[bb] <- paste0(names(rast_vnir)[bb-start_line_bnames+1], ",")
+        # }
+        # hdr_in[end_line_bnames] <- paste0(tail(names(rast_vnir), 1), "}")
+        # fileConn <-file(out_hdr)
+        # writeLines(hdr_in, fileConn)
+        # close(fileConn)
+        # writeLines(hdr_in, )
+
     }
 
-    out_file_txt <- paste0(tools::file_path_sans_ext(out_file_vnir), "_meta.txt")
+    out_file_txt <- paste0(tools::file_path_sans_ext(out_file_vnir), ".wvl")
     utils::write.table(data.frame(band = 1:length(wl_sub),
+                                  orband = orbands,
                                   wl   = wl_sub,
                                   fwhm = fwhm_sub,
                                   stringsAsFactors = FALSE),
-Original file line number
+Diff line change
@@ Expand Up / @@ -102,7 +102,6 @@ prisma_create_pan <- function(f, @@
         rm(pan_cube)
         gc()
         message("- Writing PAN raster -")
         rastwrite_lines(rast_pan, out_file_pan, out_format, proc_lev,
@@ Expand Down @@