Westerly geo predictions offset from imagery

[PatBall1]

Why do some predictions at the west of a orthomosaic become offset when reprojected into geo coordinates?

[zeppehpt1]

Hello there,

also having a similar situation, but my predictions become offset at the north when reprojected into geo coordinates.

[PatBall1]

Hi @zeppehpt1 I am away in the field until next week. I will have a look at correcting this bug as soon as I have an opportunity.

[zeppehpt1]

Thanks for letting me know.

I also did some investigations and found a way to fix it for my inference example. Maybe it also works for your predictions.

1. I passed the folder with the tiff tiles to the project_to_geojson function and iterated over the jsons and tiff tiles with zip.
2. Opened the individual tiles with rasterio to get the right reference/orientation for the georeferencing step.
3. Finally, I used rasterio to transform the coordinates into the correct geographic space:

   raster_transform = data.transform
   x_coord,y_coord = rasterio.transform.xy(transform=raster_transform,
                           rows=y_coord,
                           cols=x_coord)

However, the disadvantage of this approach is that it runs longer.

[PatBall1]

@zeppehpt1 thanks for sharing your investigations. I am aiming to fix this issue today. May I please ask, which CRS you were using when this problem arose?

[zeppehpt1]

@PatBall1 Sure, I used EPSG: 32632.

[PatBall1]

@zeppehpt1 would you mind sharing the code you described above? I am going through and integrating some changes and it would be very helpful. Thanks.

[zeppehpt1]

def project_to_geojson(data_dir, output_fold=None, pred_fold=None):  # noqa:N803
    """Projects json predictions back in geographic space.

    Takes a json and changes it to a geojson so it can overlay with orthomosaic. Another copy is produced to overlay
    with PNGs.
    """
    Path(output_fold).mkdir(parents=True, exist_ok=True)
    entries = os.listdir(pred_fold)
    entries.sort()
    print("count json",len(entries))

    data_files = glob.glob(data_dir + '*.tif')
    data_files.sort()
    print("count tiffs",len(data_files))

    for file,raster_tile in tqdm(zip(entries, data_files), total=len(entries)):
        if ".json" in file:
            data = rasterio.open(raster_tile)

            # scale to deal with the resolutio
            scalingx = data.transform[0]
            scalingy = -data.transform[4]

            # create a dictionary for each file to store data used multiple times
            img_dict = {}
            img_dict["filename"] = file

            file_mins = file.replace(".json", "")
            file_mins_split = file_mins.split("_")
            minx = int(file_mins_split[-5])
            miny = int(file_mins_split[-4])
            tile_height = int(file_mins_split[-3])
            buffer = int(file_mins_split[-2])
            height = (tile_height + 2 * buffer) / scalingx
            epsg = file_mins_split[-1]
            # create a geofile for each tile --> the EPSG value should be done
            # automatically
            geofile = {
                "type": "FeatureCollection",
                "crs": {
                    "type": "name",
                    "properties": {
                        "name": "urn:ogc:def:crs:EPSG::" + epsg
                    },
                },
                "features": [],
            }
            # update the image dictionary to store all information cleanly
            img_dict.update({"minx": minx, "miny": miny, "height": height, "buffer": buffer})
            #print("Img dict:", img_dict)

            # load the json file we need to convert into a geojson
            with open(pred_fold + "/" + img_dict["filename"]) as prediction_file:
                datajson = json.load(prediction_file)
            # print("data_json:",datajson)
            # json file is formated as a list of segmentation polygons so cycle through each one
            for crown_data in datajson:
                # just a check that the crown image is correct
                if str(minx) + "_" + str(miny) in crown_data["image_id"]:
                    crown = crown_data["segmentation"]
                    confidence_score = crown_data["score"]

                    # changing the coords from RLE format so can be read as numbers, here the numbers are
                    # integers so a bit of info on position is lost
                    mask_of_coords = mask_util.decode(crown)
                    crown_coords = polygon_from_mask(mask_of_coords)
                    if crown_coords == 0:
                        continue
                    moved_coords = []

                    # coords from json are in a list of [x1, y1, x2, y2,... ] so convert them to [[x1, y1], ...]
                    # format and at the same time rescale them so they are in the correct position for QGIS
                    for c in range(0, len(crown_coords), 2):
                        x_coord = crown_coords[c]
                        y_coord = crown_coords[c + 1]

                        raster_transform = data.transform
                        x_coord,y_coord = rasterio.transform.xy(transform=raster_transform,
                            rows=y_coord,
                            cols=x_coord)

                        moved_coords.append([x_coord, y_coord])

                    geofile["features"].append({
                        "type": "Feature",
                        "properties": {
                            "Confidence_score": confidence_score
                        },
                        "geometry": {
                            "type": "Polygon",
                            "coordinates": [moved_coords],
                        },
                    })

            # Check final form is correct - compare to a known geojson file if error appears.
            # print("geofile",geofile)

            output_geo_file = os.path.join(output_fold, img_dict["filename"].replace(".json", ".geojson"))
            # print("output location:", output_geo_file)
            with open(output_geo_file, "w") as dest:
                json.dump(geofile, dest)

This is the 'project_to_geojson' function that I used.

[PatBall1]

@zeppehpt1 I have incorporated your suggestions into the project_to_geojson function. Thanks for your help on that and please let me know if you have any other issues/suggestions.