Export time series rasters to QGis

[emmalijohansson]

Hi,

I'm wondering if there is a smooth way to export the output obtained through plotSpatial() as a raster file containing many layers? I would like to look at annual changes in carbon pools in QGis (annual data from 1987 to 2015). thanks!

Emma

[MagicForrest]

Hi Emma,

That should be easy, you can write out the Field that you plotted with plotSpatial() as a netCDF file using the writeNetCDF() function.

Alternatively you can convert the Field to a raster using as.Raster(), and then write out using the raster::writeRaster() function (you can choose many file types there).

[emmalijohansson]

Yes, I've tried that, but that just gives me one raster file to import to QGis, and I'd like to have all years and layers in one file.

Emma

Den ons 18 sep. 2019 15:36MagicForrest notifications@github.com skrev:

Hi Emma,

That should be easy, you can write out the Field that you plotted with plotSpatial() as a netCDF file using the writeNetCDF() function.

Alternatively you can convert the Field to a raster using as.Raster(), and then write out using the raster::writeRaster() function (you can choose many file types there).

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/MagicForrest/DGVMTools/issues/36?email_source=notifications&email_token=AK6KZIO6ANRCHPCTSRC5R33QKIVEXA5CNFSM4IX5QPL2YY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOD7ACTCQ#issuecomment-532687242, or mute the thread https://github.com/notifications/unsubscribe-auth/AK6KZIPESAVQZPN5FLSE2OLQKIVEXANCNFSM4IX5QPLQ .

[MagicForrest]

I am not exactly sure what format QGis expects the files to be in. But the first method (writeNetCDF() function) should put all layers and years in one file.

In the second method, the as.Raster() command might not work properly if you have multiple years and multiple variables because in such a case it is not immediately clear what one layer of the raster should represent - is it a year of a different variable.

But post the code/script that you are trying and I can take a look.

[emmalijohansson]

THE DYNAMIC MODEL

run.dynamic <- defineSource(name = "dynamic", dir = "v3/1987-2015/", format = "LPJ-GUESS")

model.field.dynamic <- getField(source = run.dynamic, var = "cpool", first.year = 1987, last.year = 2015)

plotSpatial(model.field.dynamic, layers= "Total") # this gives me maps of the total carbon pool for each year from 1987-2015. Now I want to export those raster maps.

dynamic <- as.Raster(model.field.dynamic) crs(dynamic) <- CRS('+init=EPSG:4326')

writeRaster(dynamic, "dynamic_87_15.tif", bylayer=T) # I am only interested in the total cpool, for each year.

[MagicForrest]

Okay, try slipping in a selectLayer() command to keep only the "Total" layer before in the Field before you use the Field to make the raster, like so:

model.field.dynamic.Totalonly <- selectLayer(model.field.dynamic, "Total") dynamic <- as.Raster(model.field.dynamic.Totalonly)

And then as before. That might get you further.

But also using writeNetCDF() should also work, did you try that? It should keep the layers and years intact. I guess QGis should be able to read the resulting netCDF file.

[emmalijohansson]

How would I write it using writeNetCDF()?

model.field.dynamic.total <- selectLayers(model.field.dynamic, "Total") dynamic <- as.Raster(model.field.dynamic.total) crs(dynamic) <- CRS('+init=EPSG:4326')

writeNetCDF(x = dynamic, filename = "dynamic_87_15.nc")?

[MagicForrest]

You don't need to convert it to a raster. The writeNetCDF() function writes a Field directly as a netCDF file, hopefully with all the right meta-data (such as years etc..). So do something like:

model.field.dynamic.total <- selectLayers(model.field.dynamic, "Total") writeNetCDF(x = model.field.dynamic.total, filename = "dynamic_87_15.nc")

The arguments on writeNetCDF() should be pretty well documentented.

[emmalijohansson]

I need to identify the start.date to make it a raster stack. But there are no examples for how to write that in the code, so I think the problem is that I'm not sure how to tell the code how to find the data. Sorry, but I'm not so familiar with how to work with this type of data in R.

I manage to export the data with both methods, but I only get one raster layer.

Emma

On Thu, Sep 19, 2019 at 11:02 AM MagicForrest notifications@github.com wrote:

You don't need to convert it to a raster. The writeNetCDF() function writes a Field directly as a netCDF file, hopefully with all the right meta-data (such as years etc..). So do something like:

model.field.dynamic.total <- selectLayers(model.field.dynamic, "Total") writeNetCDF(x = model.field.dynamic.total, filename = "dynamic_87_15.nc")

The arguments on writeNetCDF() should be pretty well documentented.

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/MagicForrest/DGVMTools/issues/36?email_source=notifications&email_token=AK6KZIMEPMPOFW24PWMIC7LQKM53NA5CNFSM4IX5QPL2YY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOD7CYLKQ#issuecomment-533038506, or mute the thread https://github.com/notifications/unsubscribe-auth/AK6KZION7DNWPKTZRFRPMP3QKM53NANCNFSM4IX5QPLQ .

[MagicForrest]

Hmm... it should be making multiple layers. How do you know it is making only one? Because you only see one in QGis? Or using a different method, like in R? And what file format are you using .tiff?

In any case, I will do some testing and get back to you, hopefully with a little script that works!

[emmalijohansson]

Yes, when I import the file to QGis, it only gives me one layer (both as .nc, and as .tif). Also when I use as.Raster it doesn't stack them but makes it to one layer.

Emma

On Fri, Sep 20, 2019 at 11:15 AM MagicForrest notifications@github.com wrote:

Hmm... it should be making multiple layers. How do you know it is making only one? Because you only see one in QGis? Or using a different method, like in R? And what file format are you using .tiff?

In any case, I will do some testing and get back to you, hopefully with a little script that works!

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/MagicForrest/DGVMTools/issues/36?email_source=notifications&email_token=AK6KZIJ7KX4JPKIFRDOYTODQKSICJA5CNFSM4IX5QPL2YY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOD7GDJPQ#issuecomment-533476542, or mute the thread https://github.com/notifications/unsubscribe-auth/AK6KZIKU65GYSSSED44QU6TQKSICJANCNFSM4IX5QPLQ .

[MagicForrest]

Hmm... it should be making multiple layers. How do you know it is making only one? Because you only see one in QGis? Or using a different method, like in R? And what file format are you using .tiff?

In any case, I will do some testing and get back to you, hopefully with a little script that works!

[emmalijohansson]

I get this message when I try to export as netCDF

writeNetCDF(model.field.dynamic.total, filename = "v3_2/Figures/dynamic_87_17.nc") [1] "ncvar_put: warning: you asked to write 5780 values, but the passed data array has 179180 entries!"

[MagicForrest]

Hmm, that definitely seems like a bug. Can you send me the the output of print(model.field.dynamic.total) please?

In general I am up against a couple of deadlines, but I'll try to sort this today.

[MagicForrest]

Hi Emma. There was a bug. I fixed it. You can use the fixed code either by installing the 'factorise-facetting' branch from github, or install (thought Rstudio) the package file that I just sent you via email...

And the attached script should do what you want! ;-)

script_emma_example.txt

[emmalijohansson]

Thanks a LOT!

[MagicForrest]

So, did it work out? Shall I close the issue?

[emmalijohansson]

It works, but I cannot read the NC file in QGis, but that's another issue :) I have created separate raster files for each year instead. Thanks for all the fixing!

Emma

[MagicForrest]

Okay, sounds like a QGIS, format thing, so I'll go ahead and close this issue.