From b8c4fe98e4380f634718b89d4ab2569b750c30d1 Mon Sep 17 00:00:00 2001 From: bonaime Date: Wed, 11 Sep 2024 17:21:56 +0200 Subject: [PATCH 1/2] fix for reading gcp_file.txt --- dm/opendm/gcp.py | 28 ++++++++++++++-------------- dm/opendm/location.py | 28 ++++++++++++++-------------- 2 files changed, 28 insertions(+), 28 deletions(-) diff --git a/dm/opendm/gcp.py b/dm/opendm/gcp.py index 606f827..c8dd65b 100644 --- a/dm/opendm/gcp.py +++ b/dm/opendm/gcp.py @@ -14,21 +14,21 @@ def __init__(self, gcp_path): def read(self): if self.exists(): + # Read the gcp_file.txt with open(self.gcp_path, 'r') as f: - contents = f.read().decode('utf-8-sig').encode('utf-8').strip() - - lines = map(str.strip, contents.split('\n')) - if lines: - self.raw_srs = lines[0] # SRS - self.srs = location.parse_srs_header(self.raw_srs) - - for line in lines[1:]: - if line != "" and line[0] != "#": - parts = line.split() - if len(parts) >= 6: - self.entries.append(line) - else: - log.MM_WARNING("Malformed GCP line: %s" % line) + contents = [line.rstrip() for line in f] + + # Get the spatial reference system (SRS) header from the first line + self.raw_srs = contents[0].strip() + self.srs = location.parse_srs_header(self.raw_srs) + + for line in contents[1:]: + if line != "" and line[0] != "#": + parts = line.strip().split() + if len(parts) >= 6: + self.entries.append(line) + else: + log.MM_WARNING("Malformed GCP line: %s" % line) def iter_entries(self): for entry in self.entries: diff --git a/dm/opendm/location.py b/dm/opendm/location.py index 358c013..6782018 100644 --- a/dm/opendm/location.py +++ b/dm/opendm/location.py @@ -5,7 +5,7 @@ def extract_utm_coords(photos, images_path, output_coords_file): """ - Create a coordinate file containing the GPS positions of all cameras + Create a coordinate file containing the GPS positions of all cameras to be used later in the ODM toolchain for automatic georeferecing :param photos ([ODM_Photo]) list of photos :param images_path (str) path to dataset images @@ -14,7 +14,7 @@ def extract_utm_coords(photos, images_path, output_coords_file): """ if len(photos) == 0: raise Exception("No input images, cannot create coordinates file of GPS positions") - + utm_zone = None hemisphere = None coords = [] @@ -23,7 +23,7 @@ def extract_utm_coords(photos, images_path, output_coords_file): if photo.latitude is None or photo.longitude is None or photo.altitude is None: log.MM_ERROR("Failed parsing GPS position for %s, skipping" % photo.filename) continue - + if utm_zone is None: utm_zone, hemisphere = get_utm_zone_and_hemisphere_from(photo.longitude, photo.latitude) @@ -31,12 +31,12 @@ def extract_utm_coords(photos, images_path, output_coords_file): coord = convert_to_utm(photo.longitude, photo.latitude, photo.altitude, utm_zone, hemisphere) except: raise Exception("Failed to convert GPS position to UTM for %s" % photo.filename) - + coords.append(coord) if utm_zone is None: raise Exception("No images seem to have GPS information") - + # Calculate average dx = 0.0 dy = 0.0 @@ -54,7 +54,7 @@ def extract_utm_coords(photos, images_path, output_coords_file): f.write("%s %s\n" % (dx, dy)) for coord in coords: f.write("%s %s %s\n" % (coord[0] - dx, coord[1] - dy, coord[2])) - + def transform2(from_srs, to_srs, x, y): return transformer(from_srs, to_srs).TransformPoint(x, y, 0)[:2] @@ -73,14 +73,14 @@ def proj_srs_convert(srs): else: proj4 = srs.to_proj4() res.ImportFromProj4(proj4) - + return res def transformer(from_srs, to_srs): src = proj_srs_convert(from_srs) tgt = proj_srs_convert(to_srs) return osr.CoordinateTransformation(src, tgt) - + def get_utm_zone_and_hemisphere_from(lon, lat): """ Calculate the UTM zone and hemisphere that a longitude/latitude pair falls on @@ -106,7 +106,7 @@ def convert_to_utm(lon, lat, alt, utm_zone, hemisphere): p = Proj(proj='utm',zone=utm_zone,ellps='WGS84', preserve_units=True) else: p = Proj(proj='utm',zone=utm_zone,ellps='WGS84', preserve_units=True, south=True) - + x,y = p(lon, lat) return [x, y, alt] @@ -116,17 +116,17 @@ def parse_srs_header(header): :param header (str) line :return Proj object """ - log.MM_INFO('Parsing SRS header: %s' % header) + log.MM_INFO(f'Parsing spatial reference system (SRS) header: {header}') header = header.strip() ref = header.split(' ') try: if ref[0] == 'WGS84' and ref[1] == 'UTM': datum = ref[0] utm_pole = (ref[2][len(ref[2]) - 1]).upper() - utm_zone = int(ref[2][:len(ref[2]) - 1]) - + utm_zone = int(ref[2][:-1]) + proj_args = { - 'zone': utm_zone, + 'zone': utm_zone, 'datum': datum } @@ -140,7 +140,7 @@ def parse_srs_header(header): elif header.lower().startswith("epsg:"): srs = CRS.from_epsg(header.lower()[5:]) else: - log.MM_ERROR('Could not parse coordinates. Bad SRS supplied: %s' % header) + log.MM_ERROR(f'Could not parse coordinates. Bad SRS supplied: {header}') except RuntimeError as e: log.MM_ERROR('Uh oh! There seems to be a problem with your coordinates/GCP file.\n\n' 'The line: %s\n\n' From 180e05a9cbc0edb4ef3f4fd165379d8a179ffa3a Mon Sep 17 00:00:00 2001 From: Sylvain POULAIN Date: Fri, 13 Sep 2024 12:45:42 +0400 Subject: [PATCH 2/2] Update gcp.py --- dm/opendm/gcp.py | 79 +++++++++++++++++++++++++++++++++++++++--------- 1 file changed, 65 insertions(+), 14 deletions(-) diff --git a/dm/opendm/gcp.py b/dm/opendm/gcp.py index c8dd65b..49e83f5 100644 --- a/dm/opendm/gcp.py +++ b/dm/opendm/gcp.py @@ -14,26 +14,58 @@ def __init__(self, gcp_path): def read(self): if self.exists(): - # Read the gcp_file.txt with open(self.gcp_path, 'r') as f: - contents = [line.rstrip() for line in f] - - # Get the spatial reference system (SRS) header from the first line - self.raw_srs = contents[0].strip() - self.srs = location.parse_srs_header(self.raw_srs) - - for line in contents[1:]: - if line != "" and line[0] != "#": - parts = line.strip().split() - if len(parts) >= 6: - self.entries.append(line) - else: - log.MM_WARNING("Malformed GCP line: %s" % line) + contents = f.read().strip() + + # Strip eventual BOM characters + contents = contents.replace('\ufeff', '') + + lines = list(map(str.strip, contents.split('\n'))) + if lines: + self.raw_srs = lines[0] # SRS + self.srs = location.parse_srs_header(self.raw_srs) + + for line in contents[1:]: + if line != "" and line[0] != "#": + parts = line.strip().split() + if len(parts) >= 6: + self.entries.append(line) + else: + log.MM_WARNING("Malformed GCP line: %s" % line) def iter_entries(self): for entry in self.entries: yield self.parse_entry(entry) + def check_entries(self): + coords = {} + gcps = {} + errors = 0 + + for entry in self.iter_entries(): + k = entry.coords_key() + coords[k] = coords.get(k, 0) + 1 + if k not in gcps: + gcps[k] = [] + gcps[k].append(entry) + + for k in coords: + if coords[k] < 3: + description = "insufficient" if coords[k] < 2 else "not ideal" + for entry in gcps[k]: + log.ODM_WARNING(str(entry)) + log.ODM_WARNING("The number of images where the GCP %s has been tagged are %s" % (k, description)) + log.ODM_WARNING("You should tag at least %s more images" % (3 - coords[k])) + log.ODM_WARNING("=====================================") + errors += 1 + if len(coords) < 3: + log.ODM_WARNING("Low number of GCPs detected (%s). For best results use at least 5." % (3 - len(coords))) + log.ODM_WARNING("=====================================") + errors += 1 + + if errors > 0: + log.ODM_WARNING("Some issues detected with GCPs (but we're going to process this anyway)") + def parse_entry(self, entry): if entry: parts = entry.split() @@ -51,6 +83,25 @@ def entries_count(self): def exists(self): return bool(self.gcp_path and os.path.exists(self.gcp_path)) + def make_resized_copy(self, gcp_file_output, ratio): + """ + Creates a new resized GCP file from an existing GCP file. If one already exists, it will be removed. + :param gcp_file_output output path of new GCP file + :param ratio scale GCP coordinates by this value + :return path to new GCP file + """ + output = [self.raw_srs] + + for entry in self.iter_entries(): + entry.px *= ratio + entry.py *= ratio + output.append(str(entry)) + + with open(gcp_file_output, 'w') as f: + f.write('\n'.join(output) + '\n') + + return gcp_file_output + def wgs84_utm_zone(self): """ Finds the UTM zone where the first point of the GCP falls into