Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Crop (an also when crop & mask) to the output file does not change the datatype #1128

Closed
Rapsodia86 opened this issue Apr 21, 2023 · 4 comments

Comments

@Rapsodia86
Copy link

Hello,
Was there recently any change in crop to the file?
I do have an MCD12Q2 EVI2 file that when loaded using rast() has a scale factor applied automatically. That is fine.
The problem is when I do crop to file in one line. It saves the files with the original datatype ('INT2S') and I end up with all values being 0.
When I do it in two steps: first crop() to variable and then the variable is saved using writeRaster() then it changes the datatype according to the values and everything is ok.
I am using 1.7-24 dev version on Windows.

Not long ago, everything was ok, and the output from one linear was ok!
MCD12Q2.061_EVI_Minimum_0_doy2015001_proj.zip

@rhijmans
Copy link
Member

I do not see that. This is what I see:

library(terra)
#terra 1.7.30
f <- "MCD12Q2.061_EVI_Minimum_0_doy2015001_proj.tif"
r <- rast(f)
e <- ext(160000, 400000, 620000, 950000 )

minmax(crop(r, e))
#    MCD12Q2.061_EVI_Minimum_0_doy2015001_proj
#min                                    0.1200
#max                                    0.4518
 
x = crop(r, e, filename="test.tif", overwrite=TRUE)
minmax(x)
#    MCD12Q2.061_EVI_Minimum_0_doy2015001_proj
#min                                    0.1200
#max                                    0.4518

And while the input file has a scale and offset, the output does not.

describe(f)[62]
#[1] "  Offset: 0,   Scale:0.0001"

I think that first a change was made such that crop maintains the data-type, and later an additional change, to not do that when there is a scale and offset (the alternative would have been to also write the scale and offset, but I find that too much complexity to be worth it). So perhaps your version is in-between these two changes? Can you please try with the current CRAN or devel version?

@Rapsodia86
Copy link
Author

Rapsodia86 commented Apr 25, 2023

Ok, I see the problem is due to the mask (a shapefile loaded via vect() ) within the crop function:

x = crop(r, xmask,mask=TRUE, filename="test.tif", overwrite=TRUE)
minmax(x)
#    MCD12Q2.061_EVI_Minimum_0_doy2015001_proj
#min                                         0
#max                                         0
 
x2 = crop(r, xmask,mask=TRUE)
minmax(x2)
#    MCD12Q2.061_EVI_Minimum_0_doy2015001_proj
#min                                    0.1200
#max                                    0.5029

@rhijmans
Copy link
Member

Yes, that was it. This example now also works:

library(terra)
#terra 1.7.30
r = rast("MCD12Q2.061_EVI_Minimum_0_doy2015001_proj.tif")
m = matrix(c(11000, 488000, -105000, 11000, 1065000, 669000, 669000, 1065000), ncol=2)
v = vect(m, "polygons", crs=crs(r))
 
x = crop(r, v, mask=T, filename="test.tif", overwrite=TRUE)
minmax(x)
#    MCD12Q2.061_EVI_Minimum_0_doy2015001_proj
#min                                    0.1200
#max                                    0.5029
``` 

@Rapsodia86
Copy link
Author

Thank you!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants