1

I originally wanted to ask the following question:

"In the FeaturesCollection r_features, why do the two features have such different areas? They should have the same areas because they are both buffers with a 10km radius around two different points. (Note that the points are "close enough" so that the earth's curvature shouldn't be creating "huge" differences.)"

But, no one would be able to help me with the question, because I generated these features in R, exported them to .shp, and imported them into EE.

So I decided to generate two buffers in EE around the same points, and compute their area. The answer is what I expected -- the areas are the same.

So now my question is:

Assuming that the two features that I imported in EE were constructed in the way I just described (they are both 10-km buffers around two different points), can anyone tell me why my calculation in EE is returning different areas? Perhaps something with the export-import process? But what exactly?

To be clear, I would expect some (although small) differences in the areas calculated in R and in GEE, just because of projections and machines' precisions. But I find bizarre that the two geometries have the same area in R, but different areas in GEE.

Here is my code -- and here's the link to run in in GEE, where you can also load the imported features:

// load features imported from R (public access)
var r_features = ee.FeatureCollection("users/dimitrijoe/r_features");

// create two points
var pt1 = ee.Geometry.Point([-64.64, -9.27]);
var pt2 = ee.Geometry.Point([-51.77, -3.13]);

// create buffers
var bf1 = pt1.buffer(10 * 1000);
var bf2 = pt2.buffer(10 * 1000);

// check areas of the buffers I just created (equal to the 9th decimal)
print('Area of bf1 in km2:', bf1.area().divide(1000 * 1000));
> Area of bf1 in km2: 310.4060072301467

print('Area of bf2 in km2:', bf2.area().divide(1000 * 1000));
> Area of bf2 in km2: 310.40600723013716

// calculate areas of the buffers I imported
var addArea = function(feature) {
  return feature.set({areaSqkm2: feature.geometry().area().divide(1000 * 1000)});
};

// map the above function over the FeatureCollection.
var fc = ee.FeatureCollection(r_features);
var fc = fc.map(addArea);

// Print the first feature from the collection with the added property.
var area1 = fc.filter(ee.Filter.eq('id', 206)).first().get('areaSqkm2');
var area2 = fc.filter(ee.Filter.eq('id', 1)).first().get('areaSqkm2');

print('Area of Buffer #1 in km2:', area1);
> Area of Buffer #1 in km2: 305.13

print('Area of Buffer #2 in km2:', area2);
> Area of Buffer #2 in km2: 310.15

Update: here's the R code that generates the same two features:

library(sp)

# projections    
WGS84 = "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
SAD69_POLYCONIC = "+proj=poly +lat_0=0 +lon_0=-54 +x_0=5000000 +y_0=10000000 +ellps=aust_SA +towgs84=-57,1,-41,0,0,0,0 +units=m +no_defs"

# data frame with poitns
points = data.frame(
            long = c(-64.64, -51.77), 
            lat = c(-9.27, -3.13), 
            id = c(206, 1)
        )

# create spatial object
sp = SpatialPointsDataFrame(
        cbind(points$long, points$lat), points, proj4string = CRS(WGS84))

# trabsform to planar projections to create buffer
sp = spTransform(sp, CRS(SAD69_POLYCONIC))

# create buffers
buffer = gBuffer(sp, byid = T, width = (10 * 1000))

# check area in km2
gArea(buffer, byid = T)/(1000 * 1000)
>  1        2 
>  309.017  309.017 

# export to shapefile
writeOGR(obj = buffer, dsn = TMPDIR, layer = "r_features", driver = "ESRI Shapefile")
  • 2
    Coding questions here are expected to contain code, not just a link to be followed. If a graphic context is necessary, you should include image(s) in the question. Please Edit the question. – Vince Jul 6 at 14:12
  • I included the code. There's not "graphic context" needed -- but the question is clear that one may need to look at an external object, the link for which is provided. – djas Jul 6 at 15:12
  • 2
    I think it might help this question if you included the R code you used to generate the features. That way the entire process is reproducible, and if there is something up with the export and import it'll possibly be easier to identify. – Kevin Reid Jul 6 at 15:46
  • 1
    Thanks for the suggestion, @KevinReid. I've just posted the R code necessary to make this problem reproducible. – djas Jul 6 at 18:29

Your Answer

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

Browse other questions tagged or ask your own question.