0

Below is an example of a larger process (downloading MODIS data-- multiple products for 23 cities between 2004-2016) that is giving me a bit of a hassle:

#load libraries
import ee
import datetime

#Init earth engine and load a collection
ee.Initialize()
imgs = ee.ImageCollection('MODIS/006/MOD11A2')

#filer by some dates and extract a single "band" and extract one image. I'm sure there is a cleaner way to do this
dat = imgs.filterDate('2004-01-01','2004-01-02').select('QC_Day').toList(1).get(0)

#print some info
print dat.getInfo()

#define a bounding box
bbox = (38.94993166800003, -7.127943276999929, 39.44200277300007, -6.552704129999938)

#Create the bounding box as an earth engine object thing
target_projection = ee.Image(imgs.select('QC_Day').first()).projection()
thecity = ee.Geometry.Rectangle(bbox).transform(target_projection,250) #hoping the 250 is in meters

#For some reason ee sees dat as a computed object and so clip won't work
#trick it
dat = ee.Image.cat([dat])

#print the info again
print dat.getInfo()

#clip the image to the city
dat = dat.clip(thecity)

#make the ask to EE
ask = ee.batch.Export.image.toDrive(image = dat, description = 'hello')

#start the ask
ask.start()

When the download completes and I move to R (my preference) to visualize the image:

#Import the raster package, analgous to import raster
library('raster')

#load the raster
ras = ras("C:/Users/Me/Downloads/hello.tif")
#Set all values that equal 0 to 1 and everything else to 0
zeros = ras == 0
#make a crappy map
plot(zeros)

I get the following map. My problem is the green area on the sides that looks like some combination of a projection + no data issue. The raw MODIS data (downloaded from LPDAAC) does not have this pattern. Thoughts on how to make the green go away?

R Visualization of the EE export

  • Why are you downloading a single MODIS image from GEE in the first place instead of using the canoncical source? It looks like you are displaying the data in UTM, so the no-data edges are tobe expected when the input is in an equal area, Sinusoidal projection. – Kersten Jan 17 '18 at 11:35
  • @Kersten, This is a small snippet of a larger process. In reality, I'm trying to download MODIS imagery for 23 cities over 12 years. I've done some work using the "raw" MODIS data, but the storage and memory overhead is annoying (I only need tiny parts of tiles). The underlying Proj.4 string of hello.tif (using gdalinfo) is '+proj=sinu +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs '. ArcGIS, which is hopefully better than R at accounting for only the fly projections produces a similar map. The projection matches: spatialreference.org/ref/sr-org/6965 – dothatrumba Jan 17 '18 at 16:32
  • Do the edges persist if you use a GIS to display the image in its native CRS? Reprojecting MODIS imagery these "skewed" images are the norm. – Kersten Jan 18 '18 at 8:50
  • @Kersten: If I override the projection info a la: gdal_translate -a_srs "+proj=sinu +lon_0=0 +x_0=0 +y_0=0 +a=6371007.181 +b=6371007.181 +units=m +no_defs" hello.tif helloreproj.tif the underlying map still looks the same (at least in R). – dothatrumba Jan 18 '18 at 16:38
0

The problem seems to be during the dat.clip(the_city) step. Removing that line and setting the region option of the export seems to return an image without the green edge areas. I'm not sure what the underlying problem is, but my guess is has something to do with the translation of the polygon coordinates with respect to the image's projection.

ask = ee.batch.Export.image.toDrive(image = dat, description = 'hello', region = region= thecity.getInfo().get('coordinates'))

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.