checking for not None nodata value before using the value. #BITBUCKET-3851.
This commit is contained in:
David Fisher
parent
732613cb68
commit
7153734a41
|
|
@ -399,10 +399,13 @@ def _crop_production_op(
|
|||
|
||||
"""
|
||||
result = numpy.empty(lulc_array.shape, dtype=numpy.float32)
|
||||
result[:] = _NODATA_YIELD
|
||||
valid_mask = ~numpy.isclose(lulc_array, landcover_nodata)
|
||||
if landcover_nodata is not None:
|
||||
result[:] = _NODATA_YIELD
|
||||
valid_mask = ~numpy.isclose(lulc_array, landcover_nodata)
|
||||
result[valid_mask] = 0.0
|
||||
else:
|
||||
result[:] = 0.0
|
||||
lulc_mask = lulc_array == crop_lucode
|
||||
result[valid_mask] = 0
|
||||
result[lulc_mask] = (
|
||||
yield_array[lulc_mask] * pixel_area_ha)
|
||||
return result
|
||||
|
|
@ -445,10 +448,13 @@ def _mask_observed_yield_op(
|
|||
|
||||
"""
|
||||
result = numpy.empty(lulc_array.shape, dtype=numpy.float32)
|
||||
result[:] = observed_yield_nodata
|
||||
valid_mask = ~numpy.isclose(lulc_array, landcover_nodata)
|
||||
if landcover_nodata is not None:
|
||||
result[:] = observed_yield_nodata
|
||||
valid_mask = ~numpy.isclose(lulc_array, landcover_nodata)
|
||||
result[valid_mask] = 0.0
|
||||
else:
|
||||
result[:] = 0.0
|
||||
lulc_mask = lulc_array == crop_lucode
|
||||
result[valid_mask] = 0.0
|
||||
result[lulc_mask] = (
|
||||
observed_yield_array[lulc_mask] * pixel_area_ha)
|
||||
return result
|
||||
|
|
@ -557,8 +563,11 @@ def tabulate_results(
|
|||
total_area = 0.0
|
||||
for _, band_values in pygeoprocessing.iterblocks(
|
||||
(landcover_raster_path, 1)):
|
||||
total_area += numpy.count_nonzero(
|
||||
~numpy.isclose(band_values, landcover_nodata))
|
||||
if landcover_nodata is not None:
|
||||
total_area += numpy.count_nonzero(
|
||||
~numpy.isclose(band_values, landcover_nodata))
|
||||
else:
|
||||
total_area += band_values.size
|
||||
result_table.write(
|
||||
'\n,total area (both crop and non-crop)\n,%f\n' % (
|
||||
total_area * pixel_area_ha))
|
||||
|
|
|
|||
|
|
@ -1,14 +1,12 @@
|
|||
"""InVEST Crop Production Percentile Model."""
|
||||
from __future__ import absolute_import
|
||||
import collections
|
||||
import re
|
||||
import os
|
||||
import logging
|
||||
|
||||
import numpy
|
||||
from osgeo import gdal
|
||||
from osgeo import osr
|
||||
from osgeo import ogr
|
||||
import pygeoprocessing
|
||||
import taskgraph
|
||||
|
||||
|
|
@ -594,10 +592,13 @@ def _mask_observed_yield_op(
|
|||
|
||||
"""
|
||||
result = numpy.empty(lulc_array.shape, dtype=numpy.float32)
|
||||
result[:] = observed_yield_nodata
|
||||
valid_mask = ~numpy.isclose(lulc_array, landcover_nodata)
|
||||
if landcover_nodata is not None:
|
||||
result[:] = observed_yield_nodata
|
||||
valid_mask = ~numpy.isclose(lulc_array, landcover_nodata)
|
||||
result[valid_mask] = 0.0
|
||||
else:
|
||||
result[:] = 0.0
|
||||
lulc_mask = lulc_array == crop_lucode
|
||||
result[valid_mask] = 0.0
|
||||
result[lulc_mask] = (
|
||||
observed_yield_array[lulc_mask] * pixel_area_ha)
|
||||
return result
|
||||
|
|
@ -691,8 +692,11 @@ def tabulate_regression_results(
|
|||
total_area = 0.0
|
||||
for _, band_values in pygeoprocessing.iterblocks(
|
||||
(landcover_raster_path, 1)):
|
||||
total_area += numpy.count_nonzero(
|
||||
~numpy.isclose(band_values, landcover_nodata))
|
||||
if landcover_nodata is not None:
|
||||
total_area += numpy.count_nonzero(
|
||||
~numpy.isclose(band_values, landcover_nodata))
|
||||
else:
|
||||
total_area += band_values.size
|
||||
result_table.write(
|
||||
'\n,total area (both crop and non-crop)\n,%f\n' % (
|
||||
total_area * pixel_area_ha))
|
||||
|
|
|
|||
|
|
@ -3,7 +3,6 @@ import unittest
|
|||
import tempfile
|
||||
import shutil
|
||||
import os
|
||||
import cProfile
|
||||
|
||||
import numpy
|
||||
from osgeo import gdal
|
||||
|
|
@ -33,7 +32,7 @@ class CropProductionTests(unittest.TestCase):
|
|||
def tearDown(self):
|
||||
"""Overriding tearDown function to remove temporary directory."""
|
||||
shutil.rmtree(self.workspace_dir)
|
||||
@unittest.skip("skip percentile")
|
||||
|
||||
def test_crop_production_percentile(self):
|
||||
"""Crop Production: test crop production percentile regression."""
|
||||
from natcap.invest import crop_production_percentile
|
||||
|
|
@ -50,10 +49,8 @@ class CropProductionTests(unittest.TestCase):
|
|||
'model_data_path': MODEL_DATA_PATH,
|
||||
'n_workers': '-1'
|
||||
}
|
||||
import time
|
||||
starttime = time.time()
|
||||
|
||||
crop_production_percentile.execute(args)
|
||||
print time.time() - starttime
|
||||
|
||||
agg_result_table_path = os.path.join(
|
||||
args['workspace_dir'], 'aggregate_results.csv')
|
||||
|
|
@ -77,6 +74,64 @@ class CropProductionTests(unittest.TestCase):
|
|||
pandas.testing.assert_frame_equal(
|
||||
expected_result_table, result_table, check_dtype=False)
|
||||
|
||||
def test_crop_production_percentile_no_nodata(self):
|
||||
"""Crop Production: test percentile model with undefined nodata raster.
|
||||
|
||||
Test with a landcover raster input that has no nodata value
|
||||
defined.
|
||||
"""
|
||||
from natcap.invest import crop_production_percentile
|
||||
|
||||
args = {
|
||||
'workspace_dir': self.workspace_dir,
|
||||
'results_suffix': '',
|
||||
'landcover_raster_path': os.path.join(
|
||||
SAMPLE_DATA_PATH, 'landcover.tif'),
|
||||
'landcover_to_crop_table_path': os.path.join(
|
||||
SAMPLE_DATA_PATH, 'landcover_to_crop_table.csv'),
|
||||
'model_data_path': MODEL_DATA_PATH,
|
||||
'n_workers': '-1'
|
||||
}
|
||||
|
||||
# Create a raster based on the test data geotransform, but smaller and
|
||||
# with no nodata value defined.
|
||||
base_lulc_info = pygeoprocessing.get_raster_info(
|
||||
args['landcover_raster_path'])
|
||||
base_geotransform = base_lulc_info['geotransform']
|
||||
origin_x = base_geotransform[0]
|
||||
origin_y = base_geotransform[3]
|
||||
|
||||
n = 9
|
||||
gtiff_driver = gdal.GetDriverByName('GTiff')
|
||||
raster_path = os.path.join(self.workspace_dir, 'small_raster.tif')
|
||||
new_raster = gtiff_driver.Create(
|
||||
raster_path, n, n, 1, gdal.GDT_Int32, options=[
|
||||
'TILED=YES', 'BIGTIFF=YES', 'COMPRESS=LZW',
|
||||
'BLOCKXSIZE=16', 'BLOCKYSIZE=16'])
|
||||
new_raster.SetProjection(base_lulc_info['projection'])
|
||||
new_raster.SetGeoTransform([origin_x, 1.0, 0.0, origin_y, 0.0, -1.0])
|
||||
new_band = new_raster.GetRasterBand(1)
|
||||
array = numpy.array(range(n*n), dtype=numpy.int32).reshape((n, n))
|
||||
array.fill(20) # 20 is present in the landcover_to_crop_table
|
||||
new_band.WriteArray(array)
|
||||
new_raster.FlushCache()
|
||||
new_band = None
|
||||
new_raster = None
|
||||
args['landcover_raster_path'] = raster_path
|
||||
|
||||
crop_production_percentile.execute(args)
|
||||
|
||||
result_table_path = os.path.join(
|
||||
args['workspace_dir'], 'result_table.csv')
|
||||
expected_result_table_path = os.path.join(
|
||||
TEST_DATA_PATH, 'expected_result_table_no_nodata.csv')
|
||||
expected_result_table = pandas.read_csv(
|
||||
expected_result_table_path)
|
||||
result_table = pandas.read_csv(
|
||||
result_table_path)
|
||||
pandas.testing.assert_frame_equal(
|
||||
expected_result_table, result_table, check_dtype=False)
|
||||
|
||||
def test_crop_production_percentile_bad_crop(self):
|
||||
"""Crop Production: test crop production with a bad crop name."""
|
||||
from natcap.invest import crop_production_percentile
|
||||
|
|
@ -133,7 +188,6 @@ class CropProductionTests(unittest.TestCase):
|
|||
with self.assertRaises(ValueError):
|
||||
crop_production_regression.execute(args)
|
||||
|
||||
@unittest.skip("skip percentile")
|
||||
def test_crop_production_regression(self):
|
||||
"""Crop Production: test crop production regression model."""
|
||||
from natcap.invest import crop_production_regression
|
||||
|
|
@ -155,12 +209,8 @@ class CropProductionTests(unittest.TestCase):
|
|||
'phosphorous_fertilization_rate': 8.4,
|
||||
'potassium_fertilization_rate': 14.2,
|
||||
}
|
||||
import time
|
||||
starttime = time.time()
|
||||
|
||||
crop_production_regression.execute(args)
|
||||
print time.time() - starttime
|
||||
# cProfile.runctx(
|
||||
# 'crop_production_regression.execute(args)', globals(), locals())
|
||||
|
||||
agg_result_table_path = os.path.join(
|
||||
args['workspace_dir'], 'aggregate_results.csv')
|
||||
|
|
@ -185,7 +235,7 @@ class CropProductionTests(unittest.TestCase):
|
|||
expected_result_table, result_table, check_dtype=False)
|
||||
|
||||
def test_crop_production_regression_no_nodata(self):
|
||||
"""Crop Production: test crop production regression model.
|
||||
"""Crop Production: test regression model with undefined nodata raster.
|
||||
|
||||
Test with a landcover raster input that has no nodata value
|
||||
defined.
|
||||
|
|
@ -207,9 +257,10 @@ class CropProductionTests(unittest.TestCase):
|
|||
'potassium_fertilization_rate': 14.2,
|
||||
}
|
||||
|
||||
# Create a raster based off the test data geotransform, but smaller and
|
||||
# Create a raster based on the test data geotransform, but smaller and
|
||||
# with no nodata value defined.
|
||||
base_lulc_info = pygeoprocessing.get_raster_info(args['landcover_raster_path'])
|
||||
base_lulc_info = pygeoprocessing.get_raster_info(
|
||||
args['landcover_raster_path'])
|
||||
base_geotransform = base_lulc_info['geotransform']
|
||||
origin_x = base_geotransform[0]
|
||||
origin_y = base_geotransform[3]
|
||||
|
|
@ -221,7 +272,6 @@ class CropProductionTests(unittest.TestCase):
|
|||
raster_path, n, n, 1, gdal.GDT_Int32, options=[
|
||||
'TILED=YES', 'BIGTIFF=YES', 'COMPRESS=LZW',
|
||||
'BLOCKXSIZE=16', 'BLOCKYSIZE=16'])
|
||||
# new_raster.SetProjection(srs.ExportToWkt())
|
||||
new_raster.SetProjection(base_lulc_info['projection'])
|
||||
new_raster.SetGeoTransform([origin_x, 1.0, 0.0, origin_y, 0.0, -1.0])
|
||||
new_band = new_raster.GetRasterBand(1)
|
||||
|
|
@ -233,28 +283,12 @@ class CropProductionTests(unittest.TestCase):
|
|||
new_raster = None
|
||||
args['landcover_raster_path'] = raster_path
|
||||
|
||||
import time
|
||||
starttime = time.time()
|
||||
crop_production_regression.execute(args)
|
||||
print time.time() - starttime
|
||||
# cProfile.runctx(
|
||||
# 'crop_production_regression.execute(args)', globals(), locals())
|
||||
|
||||
agg_result_table_path = os.path.join(
|
||||
args['workspace_dir'], 'aggregate_results.csv')
|
||||
expected_agg_result_table_path = os.path.join(
|
||||
TEST_DATA_PATH, 'expected_regression_aggregate_results.csv')
|
||||
expected_agg_result_table = pandas.read_csv(
|
||||
expected_agg_result_table_path)
|
||||
agg_result_table = pandas.read_csv(
|
||||
agg_result_table_path)
|
||||
pandas.testing.assert_frame_equal(
|
||||
expected_agg_result_table, agg_result_table, check_dtype=False)
|
||||
|
||||
result_table_path = os.path.join(
|
||||
args['workspace_dir'], 'result_table.csv')
|
||||
expected_result_table_path = os.path.join(
|
||||
TEST_DATA_PATH, 'expected_regression_result_table.csv')
|
||||
TEST_DATA_PATH, 'expected_regression_result_table_no_nodata.csv')
|
||||
expected_result_table = pandas.read_csv(
|
||||
expected_result_table_path)
|
||||
result_table = pandas.read_csv(
|
||||
|
|
|
|||
Loading…
Reference in New Issue