From dc6dda8c17a3b0bcbb3f7c7b73cf5a68aa487f1a Mon Sep 17 00:00:00 2001 From: Conor Doherty Date: Thu, 25 Jul 2019 15:41:29 -0700 Subject: [PATCH 1/5] general kcb test and bug fix Added point implementation of crop-specific ndvi -> kcb. The new test hit a bug in _kd_row_crop() - the +1 and **-1 were in the wrong order. This change broke some unit tests. --- openet/sims/model.py | 2 +- openet/sims/tests/test_b_model.py | 70 +++++++++++++++++++++++++++++++ 2 files changed, 71 insertions(+), 1 deletion(-) diff --git a/openet/sims/model.py b/openet/sims/model.py index da49fde..890dedc 100644 --- a/openet/sims/model.py +++ b/openet/sims/model.py @@ -428,7 +428,7 @@ def _kd_row_crop(self, fc): """ # First calculation is the fc.divide(0.7).lte(1) case return fc.multiply(self.m_l)\ - .min(fc.pow(fc.divide(0.7).multiply(self.h_max).pow(-1).add(1)))\ + .min(fc.pow(fc.divide(0.7).multiply(self.h_max).add(1).pow(-1)))\ .where( fc.divide(0.7).gt(1), fc.multiply(self.m_l).min(fc.pow(self.h_max.add(1).pow(-1))))\ diff --git a/openet/sims/tests/test_b_model.py b/openet/sims/tests/test_b_model.py index 374d538..7cc1b87 100644 --- a/openet/sims/tests/test_b_model.py +++ b/openet/sims/tests/test_b_model.py @@ -225,6 +225,76 @@ def test_Model_kc_rice_constant_value(ndvi, fc, expected, tol=0.0001): assert abs(output['kc'] - expected) <= tol +def ndvi_to_kc_point(ndvi, doy, crop_type): + crop_profile = data.cdl[crop_type] + + fc = min(max((1.26 * ndvi) - 0.18, 0), 1) + print(crop_profile['crop_class']) + if crop_profile['crop_class'] == 1: + h = crop_profile['h_max'] * min((fc / 0.7), 1) + fr = 1.0 + elif crop_profile['crop_class'] == 3 or crop_profile['crop_class'] == 2: + # Set fr based on doy + if doy < crop_profile['ls_start']: + fr = crop_profile['fr_mid'] + elif ls_start <= doy and doy <= ls_stop: + fr = fr_mid - ((doy - crop_profile['ls_start']) / + (crop_profile['ls_stop']) * (crop_profile['fr_mid'] - + crop_profile['fr_endi'])) + elif doy > crop_profile['ls_stop']: + fr = crop_profile['fr_end'] + + # Set h based on crop class + if crop_profile['crop_class'] == 3: + if fc > .5: + h = crop_profile['h_max'] + else: + h = crop_profile['h_max'] - 1 + elif crop_profile['crop_class'] == 2: + h = crop_profile['h_max'] + else: + return -1 + + kd = min(1, crop_profile['m_l'] * fc, fc ** (1/(1+h))) + kcb_full = fr * min(1 + (0.1 * crop_profile['h_max']), 1.2) + kc_min = 0.15 + kcb = kc_min + kd * (kcb_full - kc_min) + + # Crop class ceilings + if crop_profile['crop_class'] == 2: + kcb = min(kcb, 1.1) + elif crop_profile['crop_class'] == 3: + kcb = min(kcb, 1.2) + + return kcb + + +@pytest.mark.parametrize( + 'ndvi, doy, crop_type_num', + [ + # 1.26 * 0.8 - 0.18 = 0.828 + # ((0.828 ** 2) * -0.4771) + (1.4047 * 0.828) + 0.15 = 0.9859994736 + [0.8, 174, 1], + #[1.0, 200, 3], + #[.5, 200, 3], + #[.1, 200, 3], + [1.0, 200, 1], + [.5, 200, 1], + [.1, 200, 1], + [1.0, 200, 2], + [.5, 200, 2], + [.1, 200, 2], + ] +) +def test_Image_kc_constant_value(ndvi, doy, crop_type_num, tol=0.0001): + ndvi_img = ee.Image.constant(ndvi) + mod = default_model_obj(crop_type_kc_flag=True, + crop_type_source=crop_type_num) + output = utils.constant_image_value(mod.kc(ndvi_img)) + expected = ndvi_to_kc_point(ndvi, doy, crop_type_num) + assert abs(output['kc'] - expected) <= tol + + @pytest.mark.parametrize( 'fc, expected', [ From 71620cf02bb225ec5b62e3d06767ee8a06e05cc5 Mon Sep 17 00:00:00 2001 From: Charles Morton Date: Fri, 26 Jul 2019 07:25:55 -0700 Subject: [PATCH 2/5] Updated kd_row_crop test values Also some minor formatting (mainly leading zeros) --- openet/sims/tests/test_b_model.py | 32 ++++++++++++++++++------------- 1 file changed, 19 insertions(+), 13 deletions(-) diff --git a/openet/sims/tests/test_b_model.py b/openet/sims/tests/test_b_model.py index 7cc1b87..6afea7c 100644 --- a/openet/sims/tests/test_b_model.py +++ b/openet/sims/tests/test_b_model.py @@ -229,7 +229,7 @@ def ndvi_to_kc_point(ndvi, doy, crop_type): crop_profile = data.cdl[crop_type] fc = min(max((1.26 * ndvi) - 0.18, 0), 1) - print(crop_profile['crop_class']) + # print(crop_profile['crop_class']) if crop_profile['crop_class'] == 1: h = crop_profile['h_max'] * min((fc / 0.7), 1) fr = 1.0 @@ -246,7 +246,7 @@ def ndvi_to_kc_point(ndvi, doy, crop_type): # Set h based on crop class if crop_profile['crop_class'] == 3: - if fc > .5: + if fc > 0.5: h = crop_profile['h_max'] else: h = crop_profile['h_max'] - 1 @@ -255,7 +255,7 @@ def ndvi_to_kc_point(ndvi, doy, crop_type): else: return -1 - kd = min(1, crop_profile['m_l'] * fc, fc ** (1/(1+h))) + kd = min(1, crop_profile['m_l'] * fc, fc ** (1 / (1 + h))) kcb_full = fr * min(1 + (0.1 * crop_profile['h_max']), 1.2) kc_min = 0.15 kcb = kc_min + kd * (kcb_full - kc_min) @@ -276,14 +276,14 @@ def ndvi_to_kc_point(ndvi, doy, crop_type): # ((0.828 ** 2) * -0.4771) + (1.4047 * 0.828) + 0.15 = 0.9859994736 [0.8, 174, 1], #[1.0, 200, 3], - #[.5, 200, 3], - #[.1, 200, 3], + #[0.5, 200, 3], + #[0.1, 200, 3], [1.0, 200, 1], - [.5, 200, 1], - [.1, 200, 1], + [0.5, 200, 1], + [0.1, 200, 1], [1.0, 200, 2], - [.5, 200, 2], - [.1, 200, 2], + [0.5, 200, 2], + [0.1, 200, 2], ] ) def test_Image_kc_constant_value(ndvi, doy, crop_type_num, tol=0.0001): @@ -300,12 +300,18 @@ def test_Image_kc_constant_value(ndvi, doy, crop_type_num, tol=0.0001): [ # [-0.1, 0.0], [0.0, 0.0], - [0.45, 0.2418], - [0.6, 0.4454], - [0.7, 0.5857], + [0.1, 0.1668], + [0.2, 0.3591], + [0.3, 0.5229], + [0.4, 0.6521], + [0.45, 0.7051], # NDVI == 0.5 + [0.5, 0.7517], + [0.6, 0.8284], + [0.7, 0.8879], [0.8, 0.9283], + [0.9, 0.9655], [1.0, 1.0], - [1.1, 1.0], + [1.1, 1.0], # Check clamping ] ) def test_Model_kd_row_crop_constant_value(fc, expected, tol=0.0001): From 7a89ae973ff8991774e5edeea3d10863614fb292 Mon Sep 17 00:00:00 2001 From: Charles Morton Date: Fri, 26 Jul 2019 07:31:36 -0700 Subject: [PATCH 3/5] Moved kc functional test to bottom of model test script --- openet/sims/tests/test_b_model.py | 143 +++++++++++++++--------------- 1 file changed, 71 insertions(+), 72 deletions(-) diff --git a/openet/sims/tests/test_b_model.py b/openet/sims/tests/test_b_model.py index 6afea7c..555d7fc 100644 --- a/openet/sims/tests/test_b_model.py +++ b/openet/sims/tests/test_b_model.py @@ -225,76 +225,6 @@ def test_Model_kc_rice_constant_value(ndvi, fc, expected, tol=0.0001): assert abs(output['kc'] - expected) <= tol -def ndvi_to_kc_point(ndvi, doy, crop_type): - crop_profile = data.cdl[crop_type] - - fc = min(max((1.26 * ndvi) - 0.18, 0), 1) - # print(crop_profile['crop_class']) - if crop_profile['crop_class'] == 1: - h = crop_profile['h_max'] * min((fc / 0.7), 1) - fr = 1.0 - elif crop_profile['crop_class'] == 3 or crop_profile['crop_class'] == 2: - # Set fr based on doy - if doy < crop_profile['ls_start']: - fr = crop_profile['fr_mid'] - elif ls_start <= doy and doy <= ls_stop: - fr = fr_mid - ((doy - crop_profile['ls_start']) / - (crop_profile['ls_stop']) * (crop_profile['fr_mid'] - - crop_profile['fr_endi'])) - elif doy > crop_profile['ls_stop']: - fr = crop_profile['fr_end'] - - # Set h based on crop class - if crop_profile['crop_class'] == 3: - if fc > 0.5: - h = crop_profile['h_max'] - else: - h = crop_profile['h_max'] - 1 - elif crop_profile['crop_class'] == 2: - h = crop_profile['h_max'] - else: - return -1 - - kd = min(1, crop_profile['m_l'] * fc, fc ** (1 / (1 + h))) - kcb_full = fr * min(1 + (0.1 * crop_profile['h_max']), 1.2) - kc_min = 0.15 - kcb = kc_min + kd * (kcb_full - kc_min) - - # Crop class ceilings - if crop_profile['crop_class'] == 2: - kcb = min(kcb, 1.1) - elif crop_profile['crop_class'] == 3: - kcb = min(kcb, 1.2) - - return kcb - - -@pytest.mark.parametrize( - 'ndvi, doy, crop_type_num', - [ - # 1.26 * 0.8 - 0.18 = 0.828 - # ((0.828 ** 2) * -0.4771) + (1.4047 * 0.828) + 0.15 = 0.9859994736 - [0.8, 174, 1], - #[1.0, 200, 3], - #[0.5, 200, 3], - #[0.1, 200, 3], - [1.0, 200, 1], - [0.5, 200, 1], - [0.1, 200, 1], - [1.0, 200, 2], - [0.5, 200, 2], - [0.1, 200, 2], - ] -) -def test_Image_kc_constant_value(ndvi, doy, crop_type_num, tol=0.0001): - ndvi_img = ee.Image.constant(ndvi) - mod = default_model_obj(crop_type_kc_flag=True, - crop_type_source=crop_type_num) - output = utils.constant_image_value(mod.kc(ndvi_img)) - expected = ndvi_to_kc_point(ndvi, doy, crop_type_num) - assert abs(output['kc'] - expected) <= tol - - @pytest.mark.parametrize( 'fc, expected', [ @@ -386,9 +316,8 @@ def test_Model_kcb_constant_value(kd, doy, h_max, expected, tol=0.0001): assert abs(output['kcb'] - expected) <= tol - @pytest.mark.parametrize('crop_type', [1, 69, 66, 3]) -def test_Model_kc_constant_value(crop_type): +def test_Model_kc_crop_class_constant_value(crop_type): # Check that a number is returned for all crop classes output = utils.constant_image_value( default_model_obj(crop_type_source=crop_type).kc( @@ -467,3 +396,73 @@ def test_Model_kc_crop_class_2_clamping(): # # default_model_obj(crop_type_source=66, crop_type_kc_flag=False).kc( # # ndvi=ee.Image.constant(0.9))) # # assert output['kc'] == 1.2 + + +def ndvi_to_kc_point(ndvi, doy, crop_type): + crop_profile = data.cdl[crop_type] + + fc = min(max((1.26 * ndvi) - 0.18, 0), 1) + # print(crop_profile['crop_class']) + if crop_profile['crop_class'] == 1: + h = crop_profile['h_max'] * min((fc / 0.7), 1) + fr = 1.0 + elif crop_profile['crop_class'] == 3 or crop_profile['crop_class'] == 2: + # Set fr based on doy + if doy < crop_profile['ls_start']: + fr = crop_profile['fr_mid'] + elif ls_start <= doy and doy <= ls_stop: + fr = fr_mid - ((doy - crop_profile['ls_start']) / + (crop_profile['ls_stop']) * (crop_profile['fr_mid'] - + crop_profile['fr_endi'])) + elif doy > crop_profile['ls_stop']: + fr = crop_profile['fr_end'] + + # Set h based on crop class + if crop_profile['crop_class'] == 3: + if fc > 0.5: + h = crop_profile['h_max'] + else: + h = crop_profile['h_max'] - 1 + elif crop_profile['crop_class'] == 2: + h = crop_profile['h_max'] + else: + return -1 + + kd = min(1, crop_profile['m_l'] * fc, fc ** (1 / (1 + h))) + kcb_full = fr * min(1 + (0.1 * crop_profile['h_max']), 1.2) + kc_min = 0.15 + kcb = kc_min + kd * (kcb_full - kc_min) + + # Crop class ceilings + if crop_profile['crop_class'] == 2: + kcb = min(kcb, 1.1) + elif crop_profile['crop_class'] == 3: + kcb = min(kcb, 1.2) + + return kcb + + +@pytest.mark.parametrize( + 'ndvi, doy, crop_type_num', + [ + # 1.26 * 0.8 - 0.18 = 0.828 + # ((0.828 ** 2) * -0.4771) + (1.4047 * 0.828) + 0.15 = 0.9859994736 + [0.8, 174, 1], + #[1.0, 200, 3], + #[0.5, 200, 3], + #[0.1, 200, 3], + [1.0, 200, 1], + [0.5, 200, 1], + [0.1, 200, 1], + [1.0, 200, 2], + [0.5, 200, 2], + [0.1, 200, 2], + ] +) +def test_Image_kc_constant_value(ndvi, doy, crop_type_num, tol=0.0001): + ndvi_img = ee.Image.constant(ndvi) + mod = default_model_obj(crop_type_kc_flag=True, + crop_type_source=crop_type_num) + output = utils.constant_image_value(mod.kc(ndvi_img)) + expected = ndvi_to_kc_point(ndvi, doy, crop_type_num) + assert abs(output['kc'] - expected) <= tol From 2f80a70a25d6c93e24f52c74f9b4191381208b04 Mon Sep 17 00:00:00 2001 From: Conor Doherty Date: Fri, 26 Jul 2019 14:15:12 -0700 Subject: [PATCH 4/5] Fix test and add more cases --- openet/sims/model.py | 1 + openet/sims/tests/test_b_model.py | 32 ++++++++++++++++++++++++------- 2 files changed, 26 insertions(+), 7 deletions(-) diff --git a/openet/sims/model.py b/openet/sims/model.py index 890dedc..f719b14 100644 --- a/openet/sims/model.py +++ b/openet/sims/model.py @@ -402,6 +402,7 @@ def _kcb(self, kd, kc_min=0.15): return kd.multiply(kcb_full.subtract(kc_min)).add(kc_min)\ .rename(['kcb']) + def _kd_row_crop(self, fc): """Density coefficient for annual row crops (class 1) diff --git a/openet/sims/tests/test_b_model.py b/openet/sims/tests/test_b_model.py index 555d7fc..cc1ceb4 100644 --- a/openet/sims/tests/test_b_model.py +++ b/openet/sims/tests/test_b_model.py @@ -410,10 +410,10 @@ def ndvi_to_kc_point(ndvi, doy, crop_type): # Set fr based on doy if doy < crop_profile['ls_start']: fr = crop_profile['fr_mid'] - elif ls_start <= doy and doy <= ls_stop: - fr = fr_mid - ((doy - crop_profile['ls_start']) / - (crop_profile['ls_stop']) * (crop_profile['fr_mid'] - - crop_profile['fr_endi'])) + elif crop_profile['ls_start'] <= doy and doy <= crop_profile['ls_stop']: + fr = crop_profile['fr_mid'] - ((doy - crop_profile['ls_start'])\ + / (crop_profile['ls_stop'] - crop_profile['ls_start'])\ + * (crop_profile['fr_mid'] - crop_profile['fr_end'])) elif doy > crop_profile['ls_stop']: fr = crop_profile['fr_end'] @@ -455,14 +455,32 @@ def ndvi_to_kc_point(ndvi, doy, crop_type): [0.5, 200, 1], [0.1, 200, 1], [1.0, 200, 2], - [0.5, 200, 2], - [0.1, 200, 2], + [0.5, 250, 2], + [0.1, 300, 2], + [1.0, 200, 69], + [0.5, 200, 69], + [0.1, 200, 69], + [1.0, 250, 69], + [0.5, 250, 69], + [0.1, 250, 69], + [1.0, 300, 69], + [0.5, 300, 69], + [0.1, 300, 69], + [1.0, 200, 75], + [0.5, 200, 75], + [0.1, 200, 75], + [1.0, 300, 75], + [0.5, 300, 75], + [0.1, 300, 75], + [1.0, 301, 75], + [0.5, 301, 75], + [0.1, 301, 75], ] ) def test_Image_kc_constant_value(ndvi, doy, crop_type_num, tol=0.0001): ndvi_img = ee.Image.constant(ndvi) mod = default_model_obj(crop_type_kc_flag=True, - crop_type_source=crop_type_num) + crop_type_source=crop_type_num, doy=doy) output = utils.constant_image_value(mod.kc(ndvi_img)) expected = ndvi_to_kc_point(ndvi, doy, crop_type_num) assert abs(output['kc'] - expected) <= tol From a0492d29c9a9a6804fc63931cd8e8314334d97f1 Mon Sep 17 00:00:00 2001 From: Charles Morton Date: Mon, 29 Jul 2019 09:21:45 -0700 Subject: [PATCH 5/5] Bumped version number --- openet/sims/__init__.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/openet/sims/__init__.py b/openet/sims/__init__.py index 9074db4..9ea1df0 100644 --- a/openet/sims/__init__.py +++ b/openet/sims/__init__.py @@ -1,4 +1,4 @@ from .image import Image from .collection import Collection -__version__ = "0.0.5" +__version__ = "0.0.6"