Hi all, There are 365 days of soil moisture NC files and I am trying to find out how many days the values are below and above this certain threshold are repeated by R. However, I couldn't reach exactly what I wanted. For example, Daily soil moisture is below 0.3 without interrupting how many days in 365 days. NC file contains annual soil moisture values daily nctoarray <- function(ncfname, varid = NA) { nc <- nc_open(ncfname) a <- aperm(ncvar_get(nc), c(2,1,3)) nc_close(nc) a } function(x, threshold = 0.28, below = TRUE) { if (below) { y <- ifelse(x < threshold,1,0) } else y <- ifelse(x > threshold,1,0) y2 <- rle(y) sel <- which(y2$values == 1) max(y2$lengths[sel]) } m1 <- suppressWarnings(apply(a,c(1,2), consechours, 0.3, TRUE)) m2 <- suppressWarnings(apply(a,c(1,2), consechours, 0.4, FALSE)) [[alternative HTML version deleted]]
You need to make a small fake dataset that illustrates what you have and what you want out of it. Telling us you are not getting what you want is simply not useful. On August 6, 2020 8:58:09 AM PDT, "ahmet varl?" <varli61 at windowslive.com> wrote:>Hi all, > > >There are 365 days of soil moisture NC files and I am trying to find >out how many days the values are below and above this certain threshold >are repeated by R. However, I couldn't reach exactly what I wanted. For >example, Daily soil moisture is below 0.3 without interrupting how many >days in 365 days. NC file contains annual soil moisture values daily > >nctoarray <- function(ncfname, varid = NA) { nc <- nc_open(ncfname) > >a <- aperm(ncvar_get(nc), c(2,1,3)) nc_close(nc) a } > > > >function(x, threshold = 0.28, below = TRUE) { > > if (below) { > > y <- ifelse(x < threshold,1,0) > > } else y <- ifelse(x > threshold,1,0) > > > > y2 <- rle(y) > > sel <- which(y2$values == 1) > > max(y2$lengths[sel]) > >} > > > >m1 <- suppressWarnings(apply(a,c(1,2), consechours, 0.3, TRUE)) > > > >m2 <- suppressWarnings(apply(a,c(1,2), consechours, 0.4, FALSE)) > > > > > [[alternative HTML version deleted]] > >______________________________________________ >R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see >https://stat.ethz.ch/mailman/listinfo/r-help >PLEASE do read the posting guide >http://www.R-project.org/posting-guide.html >and provide commented, minimal, self-contained, reproducible code.-- Sent from my phone. Please excuse my brevity.
Hi Ahmet, I think what you are looking for can be done using run length encoding (rle). # make up some data soil_moisture<-sin(seq(0,4*pi,length.out=730))+1.1 dates<-as.Date(as.Date("2018-01-01"):as.Date("2019-12-31"), origin=as.Date("1970-01-01")) # get a logical vector for your condition under.28<-soil_moisture < 0.28 # show the soil moisture against time plot(dates,soil_moisture,pch=".",col=under.28+3,cex=2) abline(h=0.28) # use rle to get the runs of low soil moisture sm.rle<-rle(soil_moisture < 0.28) cat("Consecutive days below 0.28", paste(1:sum(sm.rle$values),sm.rle$lengths[sm.rle$values==TRUE],sep="-"), "\n") Jim On Fri, Aug 7, 2020 at 3:33 AM ahmet varl? <varli61 at windowslive.com> wrote:> > Hi all, > > > There are 365 days of soil moisture NC files and I am trying to find out how many days the values are below and above this certain threshold are repeated by R. However, I couldn't reach exactly what I wanted. For example, Daily soil moisture is below 0.3 without interrupting how many days in 365 days. NC file contains annual soil moisture values daily > > nctoarray <- function(ncfname, varid = NA) { nc <- nc_open(ncfname) > > a <- aperm(ncvar_get(nc), c(2,1,3)) nc_close(nc) a } > > > > function(x, threshold = 0.28, below = TRUE) { > > if (below) { > > y <- ifelse(x < threshold,1,0) > > } else y <- ifelse(x > threshold,1,0) > > > > y2 <- rle(y) > > sel <- which(y2$values == 1) > > max(y2$lengths[sel]) > > } > > > > m1 <- suppressWarnings(apply(a,c(1,2), consechours, 0.3, TRUE)) > > > > m2 <- suppressWarnings(apply(a,c(1,2), consechours, 0.4, FALSE)) > > > > > [[alternative HTML version deleted]] > > ______________________________________________ > R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code.
Hi Ahmet, Here is a way to get the result you ask for for one geographic grid cell. You may want more detail or something, but this is a "reproducible example". # retrieved from ftp://ftp2.psi.noaa.gov/Datasets/ncep.renalysis.dailyavgs/surface_gauss/soilw.1-10cm.gauss.1949.nc library(ncdf4) soilm<-nc_open("soilw.0-10cm.gauss.1949.nc") soil_moist<-ncvar_get(soilm) # find a suitable grid cell for(i in 1:192) { for(j in 1:94) { minmoist<-min(soil_moist[i,j,],na.rm=TRUE) if(minmoist < 0.3) cat(i,j,minmoist,"\n") } } # data is 3D numeric array use cell 159,66 soil_moisture<-soil_moist[159,66,] dates<-as.Date(as.Date("1949-01-01"):as.Date("1949-12-31"), origin=as.Date("1970-01-01")) # get a logical vector for your condition under.28<-soil_moisture < 0.28 plot(dates,soil_moisture, main="Soil moisture for grid cell 159,66 1949", col=under.28+3) abline(h=0.28) sm.rle<-rle(soil_moisture < 0.28) day_of_year<-1 cat("Consecutive days below 0.28\n") for(run in 1:length(sm.rle$lengths)) { if(sm.rle$values[run]) { cat("Day",day_of_year,"-",day_of_year+sm.rle$lengths[run], "-",sm.rle$lengths[run],"days\n") day_of_year<-day_of_year+sm.rle$lengths[run] } } Jim On Fri, Aug 7, 2020 at 11:54 AM ahmet varl? <varli61 at windowslive.com> wrote:> > Many thanks for your answer > > > nc > > File soilw_1949.nc (NC_FORMAT_NETCDF4_CLASSIC): > > > > 1 variables (excluding dimension variables): > > float soilw[lon,lat,time] > > long_name: mean Daily Volumetric Soil Moisture between 0-10 cm Below Ground Level > > units: fraction > > precision: 4 > > least_significant_digit: 3 > > GRIB_id: 144 > > GRIB_name: SOILW > > var_desc: Volumetric Soil Moisture > > dataset: NCEP Reanalysis Daily Averages > > level_desc: Between 0-10 cm BGL > > statistic: Mean > > parent_stat: Individual Obs > > missing_value: -9.96920996838687e+36 > > actual_range: 0.100000143051147 > > actual_range: 0.434000015258789 > > valid_range: 0 > > valid_range: 1 > > > > 3 dimensions: > > lon Size:192 > > units: degrees_east > > long_name: Longitude > > actual_range: 0 > > actual_range: 358.125 > > standard_name: longitude > > axis: X > > lat Size:94 > > units: degrees_north > > actual_range: 88.5419998168945 > > actual_range: -88.5419998168945 > > long_name: Latitude > > standard_name: latitude > > axis: Y > > time Size:365 *** is unlimited *** > > long_name: Time > > delta_t: 0000-00-01 00:00:00 > > avg_period: 0000-00-01 00:00:00 > > standard_name: time > > axis: T > > units: hours since 1800-01-01 00:00:0.0 > > actual_range: 1306104 > > actual_range: 1314840 > > > > 7 global attributes: > > Conventions: COARDS > > title: mean daily NMC reanalysis (1949) > > description: Data is from NMC initialized reanalysis > > (4x/day). It consists of T62 variables interpolated to > > pressure surfaces from model (sigma) surfaces. > > platform: Model > > history: created 99/05/29 by Hoop (netCDF2.3) > > Converted to chunked, deflated non-packed NetCDF4 2014/09 > > dataset_title: NCEP-NCAR Reanalysis 1 > > References: http://www.psl.noaa.gov/data/gridded/data.ncep.reanalysis.html > > > > I swiched nc to array to calculate threshold it is a 3d matrix and there is no date in files > > > > > > > dim(a) > > > > [1] 94 192 365 > > > > >a > > , , 1 > > > > [,169] [,170] [,171] [,172] [,173] [,174] [,175] [,176] [,177] [,178] [,179] [,180] [,181] [,182] > > [1,] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 > > [2,] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 > > [3,] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 > > [4,] 0.3580000 0.3580000 0.3580000 0.3580000 0.3580000 0.3580000 0.3580000 0.3577001 0.3577001 0.3577001 0.0000000 0.0000000 0.0000000 0.0000000 > > [5,] 0.3580000 0.3575001 0.3572001 0.3572001 0.3572001 0.3570001 0.3570001 0.3575001 0.3577001 0.3580000 0.3580000 0.3580000 0.3580000 0.3580000 > > > > > > > > [ reached getOption("max.print") -- omitted 89 row(s) and 364 matrix slice(s) ] > > > > > > Windows 10 i?in Posta ile g?nderildi > > > > Kimden: Jim Lemon > G?nderilme: 7 A?ustos 2020 Cuma 02:17 > Kime: ahmet varl?; r-help mailing list > Konu: Re: [R] find number of consecutive days in NC files by R > > > > Hi Ahmet, > I think what you are looking for can be done using run length encoding (rle). > > # make up some data > soil_moisture<-sin(seq(0,4*pi,length.out=730))+1.1 > dates<-as.Date(as.Date("2018-01-01"):as.Date("2019-12-31"), > origin=as.Date("1970-01-01")) > # get a logical vector for your condition > under.28<-soil_moisture < 0.28 > # show the soil moisture against time > plot(dates,soil_moisture,pch=".",col=under.28+3,cex=2) > abline(h=0.28) > # use rle to get the runs of low soil moisture > sm.rle<-rle(soil_moisture < 0.28) > cat("Consecutive days below 0.28", > paste(1:sum(sm.rle$values),sm.rle$lengths[sm.rle$values==TRUE],sep="-"), > "\n") > > Jim > > On Fri, Aug 7, 2020 at 3:33 AM ahmet varl? <varli61 at windowslive.com> wrote: > > > > Hi all, > > > > > > There are 365 days of soil moisture NC files and I am trying to find out how many days the values are below and above this certain threshold are repeated by R. However, I couldn't reach exactly what I wanted. For example, Daily soil moisture is below 0.3 without interrupting how many days in 365 days. NC file contains annual soil moisture values daily > > > > nctoarray <- function(ncfname, varid = NA) { nc <- nc_open(ncfname) > > > > a <- aperm(ncvar_get(nc), c(2,1,3)) nc_close(nc) a } > > > > > > > > function(x, threshold = 0.28, below = TRUE) { > > > > if (below) { > > > > y <- ifelse(x < threshold,1,0) > > > > } else y <- ifelse(x > threshold,1,0) > > > > > > > > y2 <- rle(y) > > > > sel <- which(y2$values == 1) > > > > max(y2$lengths[sel]) > > > > } > > > > > > > > m1 <- suppressWarnings(apply(a,c(1,2), consechours, 0.3, TRUE)) > > > > > > > > m2 <- suppressWarnings(apply(a,c(1,2), consechours, 0.4, FALSE)) > > > > > > > > > > [[alternative HTML version deleted]] > > > > ______________________________________________ > > R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see > > https://stat.ethz.ch/mailman/listinfo/r-help > > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > > and provide commented, minimal, self-contained, reproducible code. > >
Hi Ahmet, My apologies, the final for loop should read: for(run in 1:length(sm.rle$lengths)) { if(sm.rle$values[run]) { cat("Day",day_of_year,"-",day_of_year+sm.rle$lengths[run], "-",sm.rle$lengths[run],"days\n") } day_of_year<-day_of_year+sm.rle$lengths[run] } Jim On Fri, Aug 7, 2020 at 3:17 PM Jim Lemon <drjimlemon at gmail.com> wrote:> > Hi Ahmet, > Here is a way to get the result you ask for for one geographic grid > cell. You may want more detail or something, but this is a > "reproducible example". > > # retrieved from > ftp://ftp2.psi.noaa.gov/Datasets/ncep.renalysis.dailyavgs/surface_gauss/soilw.1-10cm.gauss.1949.nc > library(ncdf4) > soilm<-nc_open("soilw.0-10cm.gauss.1949.nc") > soil_moist<-ncvar_get(soilm) > # find a suitable grid cell > for(i in 1:192) { > for(j in 1:94) { > minmoist<-min(soil_moist[i,j,],na.rm=TRUE) > if(minmoist < 0.3) cat(i,j,minmoist,"\n") > } > } > # data is 3D numeric array use cell 159,66 > soil_moisture<-soil_moist[159,66,] > dates<-as.Date(as.Date("1949-01-01"):as.Date("1949-12-31"), > origin=as.Date("1970-01-01")) > # get a logical vector for your condition > under.28<-soil_moisture < 0.28 > plot(dates,soil_moisture, > main="Soil moisture for grid cell 159,66 1949", > col=under.28+3) > abline(h=0.28) > sm.rle<-rle(soil_moisture < 0.28) > day_of_year<-1 > cat("Consecutive days below 0.28\n") > for(run in 1:length(sm.rle$lengths)) { > if(sm.rle$values[run]) { > cat("Day",day_of_year,"-",day_of_year+sm.rle$lengths[run], > "-",sm.rle$lengths[run],"days\n") > day_of_year<-day_of_year+sm.rle$lengths[run] > } > } > > Jim > > On Fri, Aug 7, 2020 at 11:54 AM ahmet varl? <varli61 at windowslive.com> wrote: > > > > Many thanks for your answer > > > > > nc > > > > File soilw_1949.nc (NC_FORMAT_NETCDF4_CLASSIC): > > > > > > > > 1 variables (excluding dimension variables): > > > > float soilw[lon,lat,time] > > > > long_name: mean Daily Volumetric Soil Moisture between 0-10 cm Below Ground Level > > > > units: fraction > > > > precision: 4 > > > > least_significant_digit: 3 > > > > GRIB_id: 144 > > > > GRIB_name: SOILW > > > > var_desc: Volumetric Soil Moisture > > > > dataset: NCEP Reanalysis Daily Averages > > > > level_desc: Between 0-10 cm BGL > > > > statistic: Mean > > > > parent_stat: Individual Obs > > > > missing_value: -9.96920996838687e+36 > > > > actual_range: 0.100000143051147 > > > > actual_range: 0.434000015258789 > > > > valid_range: 0 > > > > valid_range: 1 > > > > > > > > 3 dimensions: > > > > lon Size:192 > > > > units: degrees_east > > > > long_name: Longitude > > > > actual_range: 0 > > > > actual_range: 358.125 > > > > standard_name: longitude > > > > axis: X > > > > lat Size:94 > > > > units: degrees_north > > > > actual_range: 88.5419998168945 > > > > actual_range: -88.5419998168945 > > > > long_name: Latitude > > > > standard_name: latitude > > > > axis: Y > > > > time Size:365 *** is unlimited *** > > > > long_name: Time > > > > delta_t: 0000-00-01 00:00:00 > > > > avg_period: 0000-00-01 00:00:00 > > > > standard_name: time > > > > axis: T > > > > units: hours since 1800-01-01 00:00:0.0 > > > > actual_range: 1306104 > > > > actual_range: 1314840 > > > > > > > > 7 global attributes: > > > > Conventions: COARDS > > > > title: mean daily NMC reanalysis (1949) > > > > description: Data is from NMC initialized reanalysis > > > > (4x/day). It consists of T62 variables interpolated to > > > > pressure surfaces from model (sigma) surfaces. > > > > platform: Model > > > > history: created 99/05/29 by Hoop (netCDF2.3) > > > > Converted to chunked, deflated non-packed NetCDF4 2014/09 > > > > dataset_title: NCEP-NCAR Reanalysis 1 > > > > References: http://www.psl.noaa.gov/data/gridded/data.ncep.reanalysis.html > > > > > > > > I swiched nc to array to calculate threshold it is a 3d matrix and there is no date in files > > > > > > > > > > > > > dim(a) > > > > > > > > [1] 94 192 365 > > > > > > > > >a > > > > , , 1 > > > > > > > > [,169] [,170] [,171] [,172] [,173] [,174] [,175] [,176] [,177] [,178] [,179] [,180] [,181] [,182] > > > > [1,] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 > > > > [2,] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 > > > > [3,] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 > > > > [4,] 0.3580000 0.3580000 0.3580000 0.3580000 0.3580000 0.3580000 0.3580000 0.3577001 0.3577001 0.3577001 0.0000000 0.0000000 0.0000000 0.0000000 > > > > [5,] 0.3580000 0.3575001 0.3572001 0.3572001 0.3572001 0.3570001 0.3570001 0.3575001 0.3577001 0.3580000 0.3580000 0.3580000 0.3580000 0.3580000 > > > > > > > > > > > > > > > > [ reached getOption("max.print") -- omitted 89 row(s) and 364 matrix slice(s) ] > > > > > > > > > > > > Windows 10 i?in Posta ile g?nderildi > > > > > > > > Kimden: Jim Lemon > > G?nderilme: 7 A?ustos 2020 Cuma 02:17 > > Kime: ahmet varl?; r-help mailing list > > Konu: Re: [R] find number of consecutive days in NC files by R > > > > > > > > Hi Ahmet, > > I think what you are looking for can be done using run length encoding (rle). > > > > # make up some data > > soil_moisture<-sin(seq(0,4*pi,length.out=730))+1.1 > > dates<-as.Date(as.Date("2018-01-01"):as.Date("2019-12-31"), > > origin=as.Date("1970-01-01")) > > # get a logical vector for your condition > > under.28<-soil_moisture < 0.28 > > # show the soil moisture against time > > plot(dates,soil_moisture,pch=".",col=under.28+3,cex=2) > > abline(h=0.28) > > # use rle to get the runs of low soil moisture > > sm.rle<-rle(soil_moisture < 0.28) > > cat("Consecutive days below 0.28", > > paste(1:sum(sm.rle$values),sm.rle$lengths[sm.rle$values==TRUE],sep="-"), > > "\n") > > > > Jim > > > > On Fri, Aug 7, 2020 at 3:33 AM ahmet varl? <varli61 at windowslive.com> wrote: > > > > > > Hi all, > > > > > > > > > There are 365 days of soil moisture NC files and I am trying to find out how many days the values are below and above this certain threshold are repeated by R. However, I couldn't reach exactly what I wanted. For example, Daily soil moisture is below 0.3 without interrupting how many days in 365 days. NC file contains annual soil moisture values daily > > > > > > nctoarray <- function(ncfname, varid = NA) { nc <- nc_open(ncfname) > > > > > > a <- aperm(ncvar_get(nc), c(2,1,3)) nc_close(nc) a } > > > > > > > > > > > > function(x, threshold = 0.28, below = TRUE) { > > > > > > if (below) { > > > > > > y <- ifelse(x < threshold,1,0) > > > > > > } else y <- ifelse(x > threshold,1,0) > > > > > > > > > > > > y2 <- rle(y) > > > > > > sel <- which(y2$values == 1) > > > > > > max(y2$lengths[sel]) > > > > > > } > > > > > > > > > > > > m1 <- suppressWarnings(apply(a,c(1,2), consechours, 0.3, TRUE)) > > > > > > > > > > > > m2 <- suppressWarnings(apply(a,c(1,2), consechours, 0.4, FALSE)) > > > > > > > > > > > > > > > [[alternative HTML version deleted]] > > > > > > ______________________________________________ > > > R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see > > > https://stat.ethz.ch/mailman/listinfo/r-help > > > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > > > and provide commented, minimal, self-contained, reproducible code. > > > >
Hi Ahmet, An easy way is this: library(ncdf4) soilm<-nc_open("soilw.0-10cm.gauss.1949.nc") soil_moist<-ncvar_get(soilm) smdim<-dim(soil_moist) # identify NA grid cells sm_NA_count<-matrix(NA,nrow=smdim[1],ncol=smdim[2]) for(i in 1:smdim[1]) { for(j in 1:smdim[2]) { sm_NA_count[i,j]<-sum(!is.na(soil_moist[i,j,])) } } The resulting matrix contains the counts of valid (not NA) values in each 365 day series in the array. It looks to me as though there are 5914 complete series and the rest are all NA. This does not tell you why some files (the third dimension) are all NA. Probably the best guess is that the soil moisture content is not measurable for some reason. Here is the explanation from NOAA: Missing Data: There is no missing data though the ocean has 0's. There is a file with the percent of the grid that is land. Another file has simply 1 and 0's for land/ocean. Grids where the percent of land is zero are "missing". You can get a feel for the geographic coverage like this (white cells are not NA): library(maps) library(plotrix) color2D.matplot(t(sm_NA_count)) Jim On Mon, Aug 10, 2020 at 4:09 AM ahmet varl? <varli61 at windowslive.com> wrote:> > Hi Jim, > > > > Could you help me to remove NA values which are water values ? > > > > > > Kimden: Jim Lemon > G?nderilme: 7 A?ustos 2020 Cuma 22:53 > Kime: ahmet varl? > Konu: Re: [R] find number of consecutive days in NC files by R > > > > There are 17848 grid cells in the file I downloaded for 1949. Many of > them only contain NA values, probably because they are from a > geographic grid that is covered by water. In the code there is a > section that prints out a list of the grid cells that contain minimum > values less than 0.3. Since I don't know which grid cell you are > using, I had to find one that would produce interpretable results for > the problem you are trying to solve. > > Jim > > On Fri, Aug 7, 2020 at 11:03 PM ahmet varl? <varli61 at windowslive.com> wrote: > > > > I am greatfull for your helps and ? just want to ask why did you use cell 159,66 > > > >
Hi Ahmet, I really can't work out what your problem is. I don't have access to the data you are using and so cannot inspect "a" to see what might be in it. Jim On Wed, Sep 2, 2020 at 8:54 AM ahmet varl? <varli61 at windowslive.com> wrote:> > Hi jim, > > > > I have a new question. I have 71 years raster data from 1949 to 2019 and ? am trying to convert these 71 yeears raster an array to calculate a liner regration. > > > > > > library(raster) > > r<-raster("C:/Teaching/MSCprojects/2020/Ahmet/soilm/max_consecutive_days/max_consecutive_days_1949.tif") > > a<-array(NA,dim=c(dim(r)[1:2],70)) > > i <- 1 > > for (year in 1949:2019) { > > fi<-paste0("C:/Teaching/MSCprojects/2020/Ahmet/soilm/max_consecutive_days_",year,".tif") > > r<-raster(fi) > > a[,,i]<-getValues(r,format="matrix") > > i<-i+1 > > } > > Best wishes, > > Windows 10 i?in Posta ile g?nderildi > > > > Kimden: Jim Lemon > G?nderilme: 10 A?ustos 2020 Pazartesi 06:28 > Kime: ahmet varl? > Konu: Re: [R] find number of consecutive days in NC files by R > > > > That's right. If you want to get the result for all cells with valid > readings, step through the cells, texting for valid data. The result > will be quite large, so I suggest sending the output to a file: > > sink("soil_moisture_result.txt") > for(i in 1:nrows(soil_moist)) { > for(j in 1:ncols(soil_moist)) { > if(sum(!is.na(soil_moist[i,j,]) > 0) { > # process the cell here and print out the result > } > } > sink() > > caution: untested > > Jim > > On Mon, Aug 10, 2020 at 2:02 PM ahmet varl? <varli61 at windowslive.com> wrote: > > > > Hi Jim, > > > > > > > > I would like to find out how many consecutive days each cell is under the specific value for a certain date range. ?f I am right your code is for just one cell > > > > > > > > > > > > Kimden: Jim Lemon > > G?nderilme: 10 A?ustos 2020 Pazartesi 04:23 > > Kime: ahmet varl?; r-help mailing list > > Konu: Re: [R] find number of consecutive days in NC files by R > > > > > > > > Hi Ahmet, > > An easy way is this: > > > > library(ncdf4) > > soilm<-nc_open("soilw.0-10cm.gauss.1949.nc") > > soil_moist<-ncvar_get(soilm) > > smdim<-dim(soil_moist) > > # identify NA grid cells > > sm_NA_count<-matrix(NA,nrow=smdim[1],ncol=smdim[2]) > > for(i in 1:smdim[1]) { > > for(j in 1:smdim[2]) { > > sm_NA_count[i,j]<-sum(!is.na(soil_moist[i,j,])) > > } > > } > > > > The resulting matrix contains the counts of valid (not NA) values in > > each 365 day series in the array. It looks to me as though there are > > 5914 complete series and the rest are all NA. This does not tell you > > why some files (the third dimension) are all NA. Probably the best > > guess is that the soil moisture content is not measurable for some > > reason. Here is the explanation from NOAA: > > > > Missing Data: > > > > There is no missing data though the ocean has 0's. There is a file > > with the percent of the grid that is land. Another file has simply 1 > > and 0's for land/ocean. Grids where the percent of land is zero are > > "missing". > > > > You can get a feel for the geographic coverage like this (white cells > > are not NA): > > > > library(maps) > > library(plotrix) > > color2D.matplot(t(sm_NA_count)) > > > > Jim > > > > On Mon, Aug 10, 2020 at 4:09 AM ahmet varl? <varli61 at windowslive.com> wrote: > > > > > > Hi Jim, > > > > > > > > > > > > Could you help me to remove NA values which are water values ? > > > > > > > > > > > > > > > > > > Kimden: Jim Lemon > > > G?nderilme: 7 A?ustos 2020 Cuma 22:53 > > > Kime: ahmet varl? > > > Konu: Re: [R] find number of consecutive days in NC files by R > > > > > > > > > > > > There are 17848 grid cells in the file I downloaded for 1949. Many of > > > them only contain NA values, probably because they are from a > > > geographic grid that is covered by water. In the code there is a > > > section that prints out a list of the grid cells that contain minimum > > > values less than 0.3. Since I don't know which grid cell you are > > > using, I had to find one that would produce interpretable results for > > > the problem you are trying to solve. > > > > > > Jim > > > > > > On Fri, Aug 7, 2020 at 11:03 PM ahmet varl? <varli61 at windowslive.com> wrote: > > > > > > > > I am greatfull for your helps and ? just want to ask why did you use cell 159,66 > > > > > > > > > > > > > > > >
Martin Maechler
2020-Sep-04 07:19 UTC
[R] find number of consecutive days in NC files by R
>>>>> Jeff Newmiller >>>>> on Thu, 06 Aug 2020 10:49:50 -0700 writes:> You need to make a small fake dataset that illustrates > what you have and what you want out of it. Telling us you > are not getting what you want is simply not useful. Indeed. In addition: Do *not* use suppressWarnings( . ) lightly ! Warnings are there for a good reason, and you should think hard and may be ask for help before "blindly" using suppressWarnings(). Whoever told you to do that routinely has not been a good teacher of R .. Best regards, Martin Maechler ETH Zurich and R Core team > On August 6, 2020 8:58:09 AM PDT, "ahmet varl?" > <varli61 at windowslive.com> wrote: >> Hi all, >> >> >> There are 365 days of soil moisture NC files and I am >> trying to find out how many days the values are below and >> above this certain threshold are repeated by R. However, >> I couldn't reach exactly what I wanted. For example, >> Daily soil moisture is below 0.3 without interrupting how >> many days in 365 days. NC file contains annual soil >> moisture values daily >> >> nctoarray <- function(ncfname, varid = NA) { nc <- >> nc_open(ncfname) >> >> a <- aperm(ncvar_get(nc), c(2,1,3)) nc_close(nc) a } >> >> >> >> function(x, threshold = 0.28, below = TRUE) { >> >> if (below) { >> >> y <- ifelse(x < threshold,1,0) >> >> } else y <- ifelse(x > threshold,1,0) >> >> >> >> y2 <- rle(y) >> >> sel <- which(y2$values == 1) >> >> max(y2$lengths[sel]) >> >> } >> >> >> >> m1 <- suppressWarnings(apply(a,c(1,2), consechours, 0.3, >> TRUE)) >> >> >> >> m2 <- suppressWarnings(apply(a,c(1,2), consechours, 0.4, >> FALSE)) >> >> >> >> >> [[alternative HTML version deleted]] >> >> ______________________________________________ >> R-help at r-project.org mailing list -- To UNSUBSCRIBE and >> more, see https://stat.ethz.ch/mailman/listinfo/r-help >> PLEASE do read the posting guide >> http://www.R-project.org/posting-guide.html and provide >> commented, minimal, self-contained, reproducible code. > -- > Sent from my phone. Please excuse my brevity. > ______________________________________________ > R-help at r-project.org mailing list -- To UNSUBSCRIBE and > more, see https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html and provide > commented, minimal, self-contained, reproducible code.