Check for CRS on geopandas object

[bradyrx]

I think it would be nice to have a simple helper function to add a CRS to an input geopandas object that doesn't have one.

import geopandas as gp
import pooch
import xagg as xa
import xarray as xr

# Load in example data and shapefiles
ds = xr.tutorial.open_dataset("air_temperature")
file = pooch.retrieve(
    "https://pubs.usgs.gov/of/2006/1187/basemaps/continents/continents.zip", None
)
continents = gp.read_file("zip://" + file)

# Assume this geopandas object does not have a CRS set
continents.crs = None

# Calculate weights
wm = xa.pixel_overlaps(ds, continents)

This returns the following traceback:

ValueError: Must pass either crs or epsg.
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-30-52fa646d6f28> in <module>
----> 1 wm = xa.pixel_overlaps(ds, continents)

~/miniconda3/envs/analysis/lib/python3.8/site-packages/xagg/wrappers.py in pixel_overlaps(ds, gdf_in, weights, weights_target, subset_bbox)
     52     # Get overlaps between these pixel polygons and the gdf_in polygons
     53     print('calculating overlaps between pixels and output polygons...')
---> 54     wm_out = get_pixel_overlaps(gdf_in,pix_agg)
     55     print('success!')
     56 

~/miniconda3/envs/analysis/lib/python3.8/site-packages/xagg/core.py in get_pixel_overlaps(gdf_in, pix_agg)
    249 
    250     # Match up CRSes
--> 251     pix_agg['gdf_pixels'] = pix_agg['gdf_pixels'].to_crs(gdf_in.crs)
    252 
    253     # Get GeoDataFrame of the overlaps between every pixel and the polygons

~/miniconda3/envs/analysis/lib/python3.8/site-packages/geopandas/geodataframe.py in to_crs(self, crs, epsg, inplace)
    814         else:
    815             df = self.copy()
--> 816         geom = df.geometry.to_crs(crs=crs, epsg=epsg)
    817         df.geometry = geom
    818         df.crs = geom.crs

~/miniconda3/envs/analysis/lib/python3.8/site-packages/geopandas/geoseries.py in to_crs(self, crs, epsg)
    533             crs = CRS.from_epsg(epsg)
    534         else:
--> 535             raise ValueError("Must pass either crs or epsg.")
    536 
    537         # skip if the input CRS and output CRS are the exact same

ValueError: Must pass either crs or epsg.

Maybe this is super obvious, but I was a bit confused by the error message. The pixel_overlap function doesn't allow one to pass in a CRS, so I'd say either:

1. Update the error message to say ValueError: Please add a CRS to your Geopandas object with gdf_in.set_crs("EPSG:4326")
2. Add a quick line in the function that is something like:

if gdf_in.crs is None:
    gdf_in.set_crs("ESPG:4326")

[ks905383]

Thanks for that example!

I think personally I'm in favor of just putting in a better error message; xagg probably shouldn't assume a CRS of any input to avoid any odder silent errors. Thoughts?

I'll see about putting this in the next release!

[bradyrx]

I think personally I'm in favor of just putting in a better error message; xagg probably shouldn't assume a CRS of any input to avoid any odder silent errors. Thoughts?

Yes this is probably a good philosophy and a really easy implementation. I would just offer a suggestion in the error message of how to fix this issue, since it's relatively straight-forward, but might not be for someone newer to working with shapefiles.

[huazai-beyond121]

great packages!! i have "install" this xagg packages with amount of efforts! however, when i run the example as along (basic function as weightmap = xa.pixel_overlaps(ds,gdf,weights=ds_pop.pop)), i still got this errors "ValueError: Must pass either crs or epsg"! I'm stuck here. I don't know what to do. Any help to a newer ? Thanks very much!

full errors as:

creating polygons for each pixel...
regridding weights to data grid...

/home/miniconda3/envs/cmip_env/lib/python3.8/site-packages/xarray/core/dataarray.py:853: FutureWarning: elementwise 
comparison failed; returning scalar instead, but in the future will perform elementwise comparison
  return key in self.data

calculating overlaps between pixels and output polygons...

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Input In [7], in <cell line: 2>()
      1 # Calculate overlaps
----> 2 weightmap = xa.pixel_overlaps(ds,gdf,weights=ds_pop.pop)

File ~/miniconda3/envs/cmip_env/lib/python3.8/site-packages/xagg/wrappers.py:84, in pixel_overlaps(ds, gdf_in, weights, weights_target, subset_bbox, impl)
     82 # Get overlaps between these pixel polygons and the gdf_in polygons
     83 print('calculating overlaps between pixels and output polygons...')
---> 84 wm_out = get_pixel_overlaps(gdf_in,pix_agg,impl=impl)
     85 print('success!')
     87 return wm_out

File ~/miniconda3/envs/cmip_env/lib/python3.8/site-packages/xagg/core.py:372, in get_pixel_overlaps(gdf_in, pix_agg, impl)
    369 gdf_in['poly_idx'] = np.arange(0,len(gdf_in))
    371 # Match up CRSes
--> 372 pix_agg['gdf_pixels'] = pix_agg['gdf_pixels'].to_crs(gdf_in.crs)
    374 # Get GeoDataFrame of the overlaps between every pixel and the polygons
    375 # (using the EASE grid https://nsidc.org/data/ease)
    376 if np.all(gdf_in.total_bounds[[1,3]]>0):
    377     # If min/max lat are both in NH, use North grid
    378     #epsg_set = {'init':'EPSG:6931'} (change to below bc of depreciation of {'init':...} format in geopandas)

File ~/miniconda3/envs/cmip_env/lib/python3.8/site-packages/geopandas/geodataframe.py:1364, in GeoDataFrame.to_crs(self, crs, epsg, inplace)
   1362 else:
   1363     df = self.copy()
-> 1364 geom = df.geometry.to_crs(crs=crs, epsg=epsg)
   1365 df.geometry = geom
   1366 if not inplace:

File ~/miniconda3/envs/cmip_env/lib/python3.8/site-packages/geopandas/geoseries.py:1124, in GeoSeries.to_crs(self, crs, epsg)
   1047 def to_crs(self, crs=None, epsg=None):
   1048     """Returns a ``GeoSeries`` with all geometries transformed to a new
   1049     coordinate reference system.
   1050 
   (...)
   1121 
   1122     """
   1123     return GeoSeries(
-> 1124         self.values.to_crs(crs=crs, epsg=epsg), index=self.index, name=self.name
   1125     )

File ~/miniconda3/envs/cmip_env/lib/python3.8/site-packages/geopandas/array.py:771, in GeometryArray.to_crs(self, crs, epsg)
    769     crs = CRS.from_epsg(epsg)
    770 else:
--> 771     raise ValueError("Must pass either crs or epsg.")
    773 # skip if the input CRS and output CRS are the exact same
    774 if self.crs.is_exact_same(crs):

ValueError: Must pass either crs or epsg.