Predict on entire NAIP scene

src/annotate_naip_scenes.py

@@ -223,7 +223,9 @@ def save_road_annotation_for_naip_raster(counties, naip_file, naip):
             )
             road_geometries.append(road_geometry_transformed.buffer(road_buffer_meters))
 
-            road_geometries_for_mask.append(road_geometry_transformed.buffer(ROAD_BUFFER_METERS_MASK))
+            road_geometries_for_mask.append(
+                road_geometry_transformed.buffer(ROAD_BUFFER_METERS_MASK)
+            )
 
     road_values = rasterize(
         road_geometries,

src/constants.py

@@ -35,7 +35,7 @@
 ROAD_BUFFER_METERS = {"S1100": 5.0, "S1200": 5.0}
 
 # Note: roads are buffered by this amount in order to implement CDL developed masking logic
-ROAD_BUFFER_METERS_MASK = 25.0
+ROAD_BUFFER_METERS_MASK = 30.0
 
 HAS_BUILDINGS = "has_buildings"
 HAS_ROADS = "has_roads"

src/fit_model.py

@@ -75,7 +75,9 @@ def get_y_combined(y_cdl_recoded, y_road, y_road_for_mask, y_building, label_enc
     # (2) not roads or buildings, and (3) close to a road (without being covered by a road)
     mask = label_encoder.transform(["mask"])[0]
     developed = label_encoder.transform(["developed"])[0]
-    y_combined[np.where(np.logical_and(y_road_for_mask, y_combined == developed))] = mask
+    y_combined[
+        np.where(np.logical_and(y_road_for_mask, y_combined == developed))
+    ] = mask
 
     return y_combined
 
@@ -108,7 +110,9 @@ def get_annotated_scenes(naip_paths, label_encoder, cdl_mapping):
 
             y_road = road_annotation.read()
 
-        road_annotation_for_mask_path = os.path.join(ROAD_ANNOTATION_FOR_MASK_DIR, naip_file)
+        road_annotation_for_mask_path = os.path.join(
+            ROAD_ANNOTATION_FOR_MASK_DIR, naip_file
+        )
         with rasterio.open(road_annotation_for_mask_path) as road_annotation_for_mask:
 
             y_road_for_mask = road_annotation_for_mask.read()
@@ -117,7 +121,9 @@ def get_annotated_scenes(naip_paths, label_encoder, cdl_mapping):
         with rasterio.open(building_annotation_path) as building_annotation:
             y_building = building_annotation.read()
 
-        y_combined = get_y_combined(y_cdl_recoded, y_road, y_road_for_mask, y_building, label_encoder)
+        y_combined = get_y_combined(
+            y_cdl_recoded, y_road, y_road_for_mask, y_building, label_encoder
+        )
 
         # Note: swap NAIP and CDL shape from (band, height, width) to (width, height, band)
         annotated_scenes.append([np.swapaxes(X, 0, 2), np.swapaxes(y_combined, 0, 2)])

src/prediction.py

@@ -37,6 +37,69 @@ def get_colormap(label_encoder):
     }
 
 
+def get_pixel_predictions(model, n_pixel_classes, X_normalized, image_shape):
+
+    pixel_predictions = np.zeros(X_normalized.shape[:2] + (n_pixel_classes,))
+
+    output_names = get_output_names(model)
+    pixels_output_index = np.where(np.array(output_names) == PIXELS)[0][0]
+
+    n_predictions = np.zeros_like(pixel_predictions, dtype="uint8")
+
+    # Note: predicting on entire NAIP scene requires too much memory, so we cut the image into
+    # four parts (which may slightly overlap) and predict on each of them individually
+    prediction_width = image_shape[0] * (
+        ((X_normalized.shape[0] // image_shape[0]) // 2) + 1
+    )
+    prediction_height = image_shape[1] * (
+        ((X_normalized.shape[1] // image_shape[1]) // 2) + 1
+    )
+
+    # TODO Variables for slice indices, they're used in three places
+    X_normalized_top_left = X_normalized[
+        np.newaxis, :prediction_width, :prediction_height, :
+    ]
+    X_normalized_top_right = X_normalized[
+        np.newaxis, -prediction_width:, :prediction_height, :
+    ]
+    X_normalized_bottom_left = X_normalized[
+        np.newaxis, :prediction_width, -prediction_height:, :
+    ]
+    X_normalized_bottom_right = X_normalized[
+        np.newaxis, -prediction_width:, -prediction_height:, :
+    ]
+
+    model_predictions_top_left = model.predict(X_normalized_top_left)
+    model_predictions_top_right = model.predict(X_normalized_top_right)
+    model_predictions_bottom_left = model.predict(X_normalized_bottom_left)
+    model_predictions_bottom_right = model.predict(X_normalized_bottom_right)
+
+    pixel_predictions[
+        :prediction_width, :prediction_height, :
+    ] += model_predictions_top_left[pixels_output_index][0]
+    pixel_predictions[
+        -prediction_width:, :prediction_height, :
+    ] += model_predictions_top_right[pixels_output_index][0]
+    pixel_predictions[
+        :prediction_width, -prediction_height:, :
+    ] += model_predictions_bottom_left[pixels_output_index][0]
+    pixel_predictions[
+        -prediction_width:, -prediction_height:, :
+    ] += model_predictions_bottom_right[pixels_output_index][0]
+
+    n_predictions[:prediction_width, :prediction_height, :] += 1
+    n_predictions[-prediction_width:, :prediction_height, :] += 1
+    n_predictions[:prediction_width, -prediction_height:, :] += 1
+    n_predictions[-prediction_width:, -prediction_height:, :] += 1
+
+    has_predictions = np.where(n_predictions > 0)
+    pixel_predictions[has_predictions] = (
+        pixel_predictions[has_predictions] / n_predictions[has_predictions]
+    )
+
+    return pixel_predictions
+
+
 def predict_pixels_entire_scene(
     model, naip_path, X_mean_train, X_std_train, image_shape, label_encoder, colormap
 ):
@@ -53,27 +116,11 @@ def predict_pixels_entire_scene(
 
     # Predictions have shape (width, height, n_classes)
     n_pixel_classes = len(label_encoder.classes_)
-    pixel_predictions = np.zeros(X.shape[:2] + (n_pixel_classes,))
 
-    # TODO Predict on rest of scene
-    prediction_width = (image_shape[0] // 2) * (X_normalized.shape[0] // image_shape[0])
-    prediction_height = (image_shape[1] // 2) * (
-        X_normalized.shape[1] // image_shape[1]
+    pixel_predictions = get_pixel_predictions(
+        model, n_pixel_classes, X_normalized, image_shape
     )
 
-    X_normalized = X_normalized[
-        np.newaxis, :prediction_width, :prediction_height, :
-    ].copy()
-
-    model_predictions = model.predict(X_normalized)
-
-    output_names = get_output_names(model)
-    pixels_output_index = np.where(np.array(output_names) == PIXELS)[0][0]
-
-    pixel_predictions[:prediction_width, :prediction_height, :] = model_predictions[
-        pixels_output_index
-    ][0]
-
     profile["dtype"] = str(pixel_predictions.dtype)
     profile["count"] = n_pixel_classes
 
@@ -116,7 +163,7 @@ def load_X_mean_and_std_train(model_name):
     return np.load(infile_mean), np.load(infile_std)
 
 
-def main(model_name="./saved_models/cnn_land_cover_2019_11_10_02.h5"):
+def main(model_name="./saved_models/cnn_land_cover_2019_11_10_05.h5"):
 
     config = get_config(MODEL_CONFIG)
     label_encoder, _ = get_label_encoder_and_mapping()
@@ -127,7 +174,9 @@ def main(model_name="./saved_models/cnn_land_cover_2019_11_10_02.h5"):
     )
 
     # Note: avoid ValueError: Unknown loss function:masked_categorical_crossentropy when loading saved model
-    get_custom_objects().update({"masked_categorical_crossentropy": masked_categorical_crossentropy})
+    get_custom_objects().update(
+        {"masked_categorical_crossentropy": masked_categorical_crossentropy}
+    )
 
     model = load_model(model_name)