Module pytools4dart.tools.voxreader
This module contains tools to read .vox file (output from AMAPvox), intersect with polygons and export to data.frame ready for DART simulation plots file.
Classes
voxel
class voxel(
)
AMAPvox data reader and transformer
Class methods
from_data
def from_data(
i=0,
j=0,
k=0,
pad=1.0,
min_corner=[0.0, 0.0, 0.0],
res=[1.0],
lad='Spherical',
pad_max='NA'
)
Create a voxel object from data
Parameters
----------
i: int
x index of cell
j: int
y index of cell
k: int
z index of cell
pad: float
Plant area density (m2/m3)
min_corner:
coordinates of minimum corner
res:
resolution
lad: str
leaf angle distribution
pad_max:
Maximum PAD value used in AMAPVox. If not known, it can be left to 'NA'.
Examples
--------
# create a 3D diagonal of 3 voxels of crescent density.
>>> import pytools4dart as ptd
>>> ijk = [1, 2, 3]
>>> vox = ptd.voxreader.voxel().from_data(i=ijk, j=ijk, k=ijk, pad=ijk)
>>> vox.data
i j k pad
0 1 1 1 1
1 2 2 2 2
2 3 3 3 3
>>> vox.header
{'min_corner': [0.0, 0.0, 0.0], 'max_corner': array([4., 4., 4.]), 'split': [4, 4, 4], 'type': 'pytools4dart', 'res': [1.0], 'MAX_PAD': 'NA', 'LAD_TYPE': 'Spherical'}
from_vox
def from_vox(
filename
)
Load an AMAPVox file.
Parameters
----------
filename: str
Path to an AMAPVox .vox file
Returns
-------
Examples
--------
# Load the example voxel file
>>> import pytools4dart as ptd
>>> from path import Path
>>> voxfile = Path(ptd.__file__) / 'data' /'forest.vox'
>>> vox = ptd.voxreader.voxel.from_vox(voxfile)
# See use-case 3 for a simulation case
Instance variables
extent
Returns
tuple (xmin, ymin, xmax, ymax)
Methods
affine_transform
def affine_transform(
self,
matrix,
inplace=False
)
Apply an affine transformation to the voxel grid, like with shapely.affinity.affine_transform.
The coefficient matrix is provided as a list or tuple with 6 items
for 2D transformations.
For 2D affine transformations, the 6 parameter matrix is::
[a, b, d, e, xoff, yoff]
which represents the augmented matrix::
[x'] | a b xoff | [x]
[y'] = | d e yoff | [y]
[1 ] | 0 0 1 | [1]
or the equations for the transformed coordinates::
x' = a * x + b * y + xoff
y' = d * x + e * y + yoff
Parameters
----------
matrix: list or tuple or numpy.ndarray
[a, b, d, e, xoff, yoff]
inplace: bool
If True, the grid is updated by reference and the transformation is added to header.
Otherwise the transformed grid is returned.
Returns
-------
geopandas.GeoDataFrame
Examples
--------
>>> import pytools4dart as ptd
>>> import numpy as np
>>> from path import Path
>>> voxfile = Path(ptd.__file__).parent / 'data' / 'forest.vox'
>>> vox = ptd.voxreader.voxel().from_vox(voxfile)
Get the xy coordinates of the first cell of voxel grid
>>> vox.extent
(10.0, 20.0, 30.0, 40.0)
Make a translation of the grid coordinates x-10 and y-20
>>> vox.affine_transform((1, 0, 0, 1, -10, -20), inplace=True)
>>> vox.extent
(0.0, 0.0, 20.0, 20.0)
>>> vox.header['transforms']
[(1, 0, 0, 1, -10, -20)]
intersect
def intersect(
self,
x,
columns=None,
inplace=False
)
Intersection of voxel grid with shapefile or a raster
Parameters
----------
x: geopandas.GeoDataFrame or rasterio.DatasetReader or str
Path to object that can be read by geopandas or rasterio.
GeoDataFrame is expected with polygon geometries to intersect voxel grid with.
columns: list of str
Only used with raster.
Names of the band columns when inserted in the voxel GeoDataFrame.
If None, bands are named band_{i}.
inplace: bool
If True adds intersecting ID and attributes to data, or raster bands otherwise returns dataframe.
Examples
--------
The following example makes the intersection between voxel data and an individual tree crown shapefile
containing specific leaf chemical properties, see data/README.md for details on the content.
>>> import pytools4dart as ptd
>>> import geopandas as gpd
>>> data_dir = Path(ptd.__file__).parent / 'data'
>>> vox_file = data_dir / 'forest.vox'
>>> crowns_file = data_dir / 'crowns.shp'
>>> vox = ptd.voxreader.voxel().from_vox(vox_file)
>>> vox.intersect(crowns_file) # doctest: +ELLIPSIS +NORMALIZE_WHITESPACE
i j k pad angleMean ... Can Cab CBrown Car id_crown
0 0 0 19 0.0 9.154615 ... 8.0 5.0 0.0 5.0 3.0
1 0 0 20 0.0 16.293077 ... 8.0 5.0 0.0 5.0 3.0
2 0 0 21 0.0 16.293077 ... 8.0 5.0 0.0 5.0 3.0
3 0 0 22 0.0 16.293077 ... 8.0 5.0 0.0 5.0 3.0
4 0 0 23 0.0 9.154615 ... 8.0 5.0 0.0 5.0 3.0
... .. .. .. ... ... ... ... ... ... ... ...
22395 19 19 10 0.0 28.191290 ... NaN NaN NaN NaN NaN
22396 19 19 11 0.0 28.191290 ... NaN NaN NaN NaN NaN
22397 19 19 12 0.0 22.092730 ... NaN NaN NaN NaN NaN
22398 19 19 13 0.0 23.520610 ... NaN NaN NaN NaN NaN
22399 19 19 14 0.0 23.520610 ... NaN NaN NaN NaN NaN
<BLANKLINE>
[22400 rows x 24 columns]
>>> raster_file = Path(ptd.__file__).parent / 'data' / 'Can_Cab_Car_CBrown.tif'
>>> band_names = raster_file.stem.split('_')
>>> vox.intersect(raster_file, columns=band_names) # doctest: +ELLIPSIS +NORMALIZE_WHITESPACE
i j k pad ... Can Cab Car CBrown
0 0 0 19 0.0 ... 1.752778 33.722435 7.891806 0.271965
1 0 0 20 0.0 ... 1.752778 33.722435 7.891806 0.271965
2 0 0 21 0.0 ... 1.752778 33.722435 7.891806 0.271965
3 0 0 22 0.0 ... 1.752778 33.722435 7.891806 0.271965
4 0 0 23 0.0 ... 1.752778 33.722435 7.891806 0.271965
... .. .. .. ... ... ... ... ... ...
22395 19 19 10 0.0 ... 1.015507 35.692669 9.584843 0.205697
22396 19 19 11 0.0 ... 1.015507 35.692669 9.584843 0.205697
22397 19 19 12 0.0 ... 1.015507 35.692669 9.584843 0.205697
22398 19 19 13 0.0 ... 1.015507 35.692669 9.584843 0.205697
22399 19 19 14 0.0 ... 1.015507 35.692669 9.584843 0.205697
<BLANKLINE>
[22400 rows x 20 columns]
Returns
-------
geopandas.GeoDataFrame
reduce_xy
def reduce_xy(
self,
inplace=False
)
Shift the grid minimum corner to x,y=(0,0).
Parameters
----------
inplace: bool
If True, the grid is updated by reference and the transformation is added to header.
Otherwise the transformed grid is returned.
Returns
-------
None or (geopandas.GeoDataFrame and transformation array)
If inplace=False, returns the grid transformed and the transformation array.
The transformation array that can be used then to transform back simulation output rasters.
Examples
--------
>>> import pytools4dart as ptd
>>> from path import Path
>>> data_dir = Path(ptd.__file__) / 'data'
>>> voxfile = data_dir / 'forest.vox'
>>> vox = ptd.voxreader.voxel.from_vox(voxfile)
>>> vox.extent
(10.0, 20.0, 30.0, 40.0)
>>> vox.reduce_xy(inplace=True)
>>> vox.extent
(0.0, 0.0, 20.0, 20.0)
to_plots
def to_plots(
self,
pa_type='UL',
keep_columns=None,
reduce_xy=False,
pa_column='pad',
**kwargs
)
Convert to DART plots DataFrame to be included in simulation.
Parameters
----------
pa_type: str
If 'UL', pa_column is considered as a Plant Area Density (m2/m3)
If 'LAI', pa_column is considered as a Plant Area Index (m2/m2)
keep_columns: str or list of str
Columns from data to keep in plots DataFrame. If 'all',
reduce_xy: bool
If True, shift the grid minimum corner x,y=(0,0).
In that case, the transformation array is also returned.
pa_column: str
The plant or leaf area column name
kwargs: for retro-compatibility
Returns
-------
pandas.DataFrame | (pandas.DataFrame, list)
The plots DataFrame in DART plot file format.
If reduce_xy=True, the affine_transform parameters.
Examples
--------
>>> import pytools4dart as ptd
>>> ijk = [1, 2, 3]
>>> vox = ptd.voxreader.voxel().from_data(i=ijk, j=ijk, k=ijk, pad=ijk)
>>> vox.to_plots() # doctest: +ELLIPSIS +NORMALIZE_WHITESPACE
PLT_TYPE PT_1_X PT_1_Y ... PLT_HEI_MEA VEG_DENSITY_DEF VEG_UL
0 1 2.0 1.0 ... 1.0 1 1
1 1 3.0 2.0 ... 1.0 1 2
2 1 4.0 3.0 ... 1.0 1 3
<BLANKLINE>
[3 rows x 13 columns]
# See use case 3 and 6 for simulation cases
to_raster
def to_raster(
self,
raster_file,
crs=None,
use_transform=True,
aggregate_fun=None,
reproject=False
)
Convert to raster stack
Parameters
----------
raster_file: str
Raster file path.
crs: rasterio.crs
Coordinates reference system.
use_transform: bool
If True, the transformations stored in self.header are applied.
aggregate_fun: function
Function to aggregate the LAI/LAD values.
reproject: bool
If True and the transformation contains a rotation,
the raster is regridded on an x,y grid,
after applying the rotation that is not well supported by some softwares.
In that case the grid is aligned to x,y=(0,0),
the grid resolution is the same as the voxel grid,
and the resampling method is the nearest neighbour.
Returns
-------
str
Raster file path
Examples
--------
>>> import pytools4dart as ptd
>>> from path import Path
>>> import numpy as np
>>> voxfile = Path(ptd.__file__).parent / 'data' / 'forest.vox'
>>> vox = ptd.voxreader.voxel.from_vox(voxfile)
>>> vop = np.array([[0.453990499739275, 0.891006524188506, 0.0, -650363.659927172], [-0.891006524188506, 0.453990499739275, 0.0, -9491.23042430705], [0.0, 0.0, 1.0, 0.0], [0.0, 0.0, 0.0, 1.0]])
>>> ivop = np.linalg.inv(vop)
>>> ivop2D = (ivop[0,0],ivop[0,1],ivop[1,0],ivop[1,1],ivop[0,3],ivop[1,3])
>>> vox.affine_transform(ivop2D, inplace=True)
>>> raster_file = vox.to_raster('/tmp/test.tif', crs = '+init=epsg:2792')