-1

I ordered netCDF chlorophyll data from 2002 to 2010 of the Caspian Sea region from http://hermes.acri.fr. Then, I got 726 .nc files listed in a folder and downloaded them from the FTP site. So now I have 726 time series chlorophyll data in my working directory (time series data and the size is 61 MB). I am trying to open all of the data (because I want to study changes through time) by ncdf4package and extract chlorophyll-a values (pixels):

mydata<- nc_open('L3m_20020423-20020430__631633423_4_AV-SWF_CHL1_8D_00')
mydata_extract<- ncvar_get(mydata)

But the problem is that I have 726 netCDF file and I have no idea how I can read and extract all extracted chlorophyll-a values into one dataframe or matrix. In fact I need a table rows with longitude, latitude, time and chlorophyll data from all of these 726 files in one table to apply K-means clustering and classify chlorophyll data based on time and region. Anyone has a clue?

  • 1
    Why can't you loop over the file names? How do you want to combine the 726 files? Are they at different times, or do they mosaic over a wider space? – Spacedman Oct 30 '17 at 21:53
  • Actually all of 726 files are time series and also mosaic data from a lake. – mohammad gerami Oct 31 '17 at 19:39
  • So what do you want the result to look like? – Spacedman Oct 31 '17 at 19:51
  • I want all of the 726 files into a datafram for further analysis. But I do not know how to write a script to read all of it into a matrix – mohammad gerami Oct 31 '17 at 20:04
  • What columns do you want in your result data frame? What's in mydata_extract after your code in your question? Can you point us to a sample data file? – Spacedman Oct 31 '17 at 20:40
1

I just downloaded a random bunch of netCDF files from http://hermes.acri.fr/index.php?class=demonstration_products - I picked a couple of years and drew a region off the coast of the UK. I got 25 files listed and downloaded them from the FTP site.

Let's open one up and see what's in it:

> library(ncdf4)
> CHL1 <- nc_open('./L3m_20020901-20020930__802386687_25_AVW-MERMODSWF_CHL-CIA_MO_00.nc')
> names(CHL1$var)
[1] "CHL-CIA_mean"  "CHL-CIA_flags"

There are two variables. The first one is:

> CHL1$var[[1]]$varsize
[1] 134  48

a grid of 134x48 values. I can get that with:

> CHLmean = ncvar_get(CHL1, CHL1$var[[1]])
> dim(CHLmean)
[1] 134  48

I can get the X and Y coordinates of the grid from the variable metadata:

> (longitude = CHL1$var[[1]]$dim[[1]]$vals)
  [1] -30.875 -30.625 -30.375 -30.125 -29.875 -29.625 -29.375 -29.125 -28.875
 [10] -28.625 -28.375 -28.125 -27.875 -27.625 -27.375 -27.125 -26.875 -26.625
 [19] -26.375 -26.125 -25.875 -25.625 -25.375 -25.125 -24.875 -24.625 -24.375
 [28] -24.125 -23.875 -23.625 -23.375 -23.125 -22.875 -22.625 -22.375 -22.125
> (latitude = CHL1$var[[1]]$dim[[2]]$vals)
 [1] 59.625 59.375 59.125 58.875 58.625 58.375 58.125 57.875 57.625 57.375
[11] 57.125 56.875 56.625 56.375 56.125 55.875 55.625 55.375 55.125 54.875
[21] 54.625 54.375 54.125 53.875 53.625 53.375 53.125 52.875 52.625 52.375

I can plot this as an image but I have to flip the coordinates:

> image(list(x=longitude, y=rev(latitude), z=CHLmean[,length(latitude):1]))

enter image description here

Another easier option is to read the data in directly with the raster package:

> CHLr = raster("./L3m_20020901-20020930__802386687_25_AVW-MERMODSWF_CHL-CIA_MO_00.nc")
Warning message:
In .varName(nc, varname, warn = warn) : varname used is: CHL-CIA_mean
If that is not correct, you can set it to one of: CHL-CIA_mean, CHL-CIA_flags

> CHLr
class       : RasterLayer 
dimensions  : 48, 134, 6432  (nrow, ncol, ncell)
resolution  : 0.25, 0.25  (x, y)
extent      : -31, 2.5, 47.75, 59.75  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
data source : /nobackup/rowlings/Downloads/Hermes/L3m_20020901-20020930__802386687_25_AVW-MERMODSWF_CHL-CIA_MO_00.nc 
names       : CHL.CIA...Mean.of.the.binned.pixels 
zvar        : CHL-CIA_mean 


> plot(CHLr)

enter image description here

That warns that it has only read the first variable, but if you only have one variable in each NetCDF then that'll do.

So, that's how you read data in from one of them, or at least from one file that might have a good change of being a bit like the files you've got (but you've supplied no data or metadata so its a bit of a guess).

What do you want to do with 726 of them? How can you turn the things into a data frame? Do you want rows with X, Y, Time, and CHL1_mean or whatever variables you've got? This remains unclear. You should at least now edit your question to show us more of the metadata from the NetCDF you have, including the variable names, sizes, and dimensions.

  • Question has been edited @Spacedman – mohammad gerami Nov 4 '17 at 9:45
  • 1
    Given a raster read in as above, rasterToPoints(r) will get you a data frame of x,y, CHL_1mean, and then you can cbind the time on it, repeat over all your files. – Spacedman Nov 4 '17 at 15:15
  • How do you open all of the 25 files list of the coast of UK simultaneously? – mohammad gerami Nov 4 '17 at 20:35
  • 1
    You don't open them all simultaneously, you write a loop, open one at a time, then merge the results. That's basic R programming. – Spacedman Nov 5 '17 at 11:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.