mdbartos / pysheds

:earth_americas: Simple and fast watershed delineation in python.
GNU General Public License v3.0
718 stars 195 forks source link

trying to plot dem #83

Closed tommylees112 closed 4 years ago

tommylees112 commented 5 years ago

I read in the dem file from raster as shown on the README.md file

grid = Grid.from_raster(BASE_DATA_DIR / 'hydrosheds_dem'/ 'n00e005_con/', data_name='dem')
grid.read_raster(BASE_DATA_DIR / 'hydrosheds_dem'/ 'n00e005_con/'/ 'n00e005_con/', data_name='dir')

The file was downloaded from here. It's BIL format.

The file looks like this:

n00e005_con
├── info
│   ├── arc0000.dat
│   ├── arc0000.nit
│   ├── arc0001.dat
│   ├── arc0001.nit
│   ├── arc0002.dat
│   ├── arc0002.nit
│   ├── arc0002r.001
│   └── arc.dir
└── n00e005_con
    ├── dblbnd.adf
    ├── hdr.adf
    ├── metadata.htm
    ├── metadata.xml
    ├── prj.adf
    ├── sta.adf
    ├── vat.adf
    ├── w001001.adf
    └── w001001x.adf

2 directories, 17 files

Then I try and plot as in the quickstart notebook

plt.imshow(grid.dem, extent=grid.extent, cmap='cubehelix', zorder=1)
plt.colorbar(label='Elevation (m)')
plt.grid(zorder=0)
plt.title('Digital elevation map')
plt.xlabel('Longitude')
plt.ylabel('Latitude')
plt.tight_layout()
plt.savefig('img/conditioned_dem.png', bbox_inches='tight')

But I get an error message:

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-74-dccabe6b56be> in <module>
      1 plt.imshow(grid.dem, extent=grid.extent, cmap='cubehelix', zorder=1)
----> 2 plt.colorbar(label='Elevation (m)')
      3 plt.grid(zorder=0)
      4 plt.title('Digital elevation map')
      5 plt.xlabel('Longitude')

~/.conda/envs/crp/lib/python3.7/site-packages/matplotlib/pyplot.py in colorbar(mappable, cax, ax, **kw)
   2098         ax = gca()
   2099
-> 2100     ret = gcf().colorbar(mappable, cax = cax, ax=ax, **kw)
   2101     return ret
   2102 colorbar.__doc__ = matplotlib.colorbar.colorbar_doc

~/.conda/envs/crp/lib/python3.7/site-packages/matplotlib/figure.py in colorbar(self, mappable, cax, ax, use_gridspec, **kw)
   2127                              'panchor']
   2128         cb_kw = {k: v for k, v in kw.items() if k not in NON_COLORBAR_KEYS}
-> 2129         cb = cbar.colorbar_factory(cax, mappable, **cb_kw)
   2130
   2131         self.sca(current_ax)

~/.conda/envs/crp/lib/python3.7/site-packages/matplotlib/colorbar.py in colorbar_factory(cax, mappable, **kwargs)
   1565         cb = ColorbarPatch(cax, mappable, **kwargs)
   1566     else:
-> 1567         cb = Colorbar(cax, mappable, **kwargs)
   1568
   1569     cid = mappable.callbacksSM.connect('changed', cb.on_mappable_changed)

~/.conda/envs/crp/lib/python3.7/site-packages/matplotlib/colorbar.py in __init__(self, ax, mappable, **kw)
   1096                 kw['alpha'] = mappable.get_alpha()
   1097
-> 1098             ColorbarBase.__init__(self, ax, **kw)
   1099
   1100     def on_mappable_changed(self, mappable):

~/.conda/envs/crp/lib/python3.7/site-packages/matplotlib/colorbar.py in __init__(self, ax, cmap, norm, alpha, values, boundaries, orientation, ticklocation, extend, spacing, ticks, format, drawedges, filled, extendfrac, extendrect, label)
    412             self.formatter = format  # Assume it is a Formatter
    413         # The rest is in a method so we can recalculate when clim changes.
--> 414         self.draw_all()
    415
    416     def _extend_lower(self):

~/.conda/envs/crp/lib/python3.7/site-packages/matplotlib/colorbar.py in draw_all(self)
    436         # sets self._boundaries and self._values in real data units.
    437         # takes into account extend values:
--> 438         self._process_values()
    439         # sets self.vmin and vmax in data units, but just for
    440         # the part of the colorbar that is not part of the extend

~/.conda/envs/crp/lib/python3.7/site-packages/matplotlib/colorbar.py in _process_values(self, b)
    842                 self.norm.vmin,
    843                 self.norm.vmax,
--> 844                 expander=0.1)
    845
    846             b = self.norm.inverse(self._uniform_y(self.cmap.N + 1))

~/.conda/envs/crp/lib/python3.7/site-packages/matplotlib/transforms.py in nonsingular(vmin, vmax, expander, tiny, increasing)
   2905             vmax = expander
   2906         else:
-> 2907             vmin -= expander*abs(vmin)
   2908             vmax += expander*abs(vmax)
   2909

~/.conda/envs/crp/lib/python3.7/site-packages/numpy/ma/core.py in __isub__(self, other)
   4195             self._mask += m
   4196         self._data.__isub__(np.where(self._mask, self.dtype.type(0),
-> 4197                                      getdata(other)))
   4198         return self
   4199

TypeError: Cannot cast ufunc subtract output from dtype('float64') to dtype('int16') with casting rule 'same_kind'

And this plot: Screenshot 2019-03-26 at 18 41 49

mdbartos commented 5 years ago

Greetings,

I don't think I've ever tried using a BIL file before. All the tests use the GRID-type files instead.

kassouna commented 4 years ago

I am facing the Same issue using GRID!

any solution! thanks in advance

mdbartos commented 4 years ago

If the underlying data is int16, you may need to cast to float to make the colorbar work correctly. So, when using imshow, try:

plt.imshow(grid.dem.astype(float), extent=grid.extent, cmap='cubehelix', zorder=1)

kassouna commented 4 years ago

it works ! thanks