import os
import structlog
from os.path import join, basename
import shutil
from lxml import etree
import textwrap
from .util import get_all_data_files, gdal, replace_ext
from .steps.georeference import set_gcps_from_gsc_footprint
logger = structlog.getLogger(__name__)
PROJECTIONS = {
"EUROPEAN": "EPSG:3035",
"ES_29": "EPSG:32629",
"ES_30": "EPSG:32630",
}
[docs]def rename_reference_dimap(
source_dir: os.PathLike,
target_dir: os.PathLike,
preprocess_config: dict,
search: str,
replace: str,
output_file_name: str = "imagery.dim",
):
filenames = get_all_data_files(source_dir, preprocess_config, ["*metadata.dim"])
# working in source_dir to be picked up later by move operation
target_filename = join(source_dir, output_file_name)
if len(filenames) != 0:
try:
# read file
with open(filenames[0]) as f:
# edit content
dimap = f.read()
replaced = dimap.replace(search, replace)
# write resulting file
with open(target_filename, "w") as w:
w.write(replaced)
os.unlink(filenames[0])
except Exception as e:
logger.warn(e)
try:
# check if .hdr file is created
# if not create it to source_dir to be picked up later by move operation
attempt_to_create_hdr(source_dir, source_dir, preprocess_config)
except Exception as e:
logger.warn(e)
try:
attempt_to_set_projection(source_dir, source_dir, preprocess_config)
except Exception as e:
logger.warn(e)
# copy everything else to target_dir
for input_filename in get_all_data_files(source_dir, preprocess_config, "**"):
target_filename = join(target_dir, basename(input_filename))
shutil.move(input_filename, target_filename)
[docs]def attempt_to_create_hdr(
source_dir, target_dir, preprocess_config, output: str = "imagery.hdr"
):
hdr_or_tif = get_all_data_files(source_dir, preprocess_config, ["*hdr", "*tif"])
if len(hdr_or_tif) == 0:
# attempt to recreate it from metadata.xml in source_dir
md_files = get_all_data_files(source_dir, preprocess_config, ["*metadata.xml"])
if len(md_files) >= 1:
hdr_data = {}
with open(md_files[0]) as g:
logger.debug("Parsing metadata file %s" % md_files[0])
tree = etree.parse(g)
root = tree.getroot()
for child in root:
if child.tag == "Image":
for child2 in child:
hdr_data[child2.tag] = child2.text
break
if len(hdr_data) > 0:
logger.debug("Used metadata for hdr creation %s" % hdr_data)
target_filename = join(target_dir, output)
with open(target_filename, "w") as hdr_output:
hdr_content = textwrap.dedent(
"""\
BYTEORDER M
LAYOUT BIL
NROWS {NROWS}
NCOLS {NCOLS}
NBANDS 4
NBITS 8
SKIPBYTES 0""".format(
NROWS=hdr_data["ROWS"],
NCOLS=hdr_data["COLUMNS"],
)
)
hdr_output.write(hdr_content)
logger.debug(
"Successfully written .hdr metadata file to %s" % target_filename
)
else:
logger.warn(
"Failed extraction of IMAGE xml metadata for .hdr creation."
)
else:
logger.warn(
".Hdr file not present but also metadata.xml not found. No way to fix"
".bil"
)
[docs]def attempt_to_set_projection(source_dir, target_dir, preprocess_config):
# check if there are no tifs and only bil + metadata.xml without dimap file
# search for GeoInformation><PROJECTION>European</PROJECTION></GeoInformation> to set
# crs for image from mapping by force
needed_files = get_all_data_files(
source_dir, preprocess_config, ["*.tif", "*.dim", "*.prj"]
)
if len(needed_files) > 0:
# expected files are there, return
return
# this means we have no dimap and a bil with possible georeferencing in metadata.xml
# (which we need to extract out)
md_files = get_all_data_files(source_dir, preprocess_config, ["*metadata.xml"])
if len(md_files) >= 1:
projection = None
with open(md_files[0]) as g:
logger.debug("Parsing metadata file %s" % md_files[0])
tree = etree.parse(g)
root = tree.getroot()
for child in root:
if child.tag == "GeoInformation":
for child2 in child:
if child2.tag == "PROJECTION":
projection = child2.text.upper()
break
break
if projection is not None:
logger.info("Found projection %s" % projection)
# get epsg code of projection by name
epsg_code = PROJECTIONS.get(projection, None)
input_filename = get_all_data_files(
source_dir, preprocess_config, ["*.bil"]
)[0]
target_filename = replace_ext(
join(target_dir, basename(input_filename)), ".tif"
)
# create a tif from bil and set its projection
ds = gdal.Open(input_filename)
ds = gdal.GetDriverByName("GTiff").CreateCopy(target_filename, ds, strict=0)
if not epsg_code:
logger.warning(
"Projection %s not in the mapping of projections" % projection
)
# use GSC bbox as a fallback if everything failed, risking bad orientation
set_gcps_from_gsc_footprint(ds, preprocess_config)
else:
ds.SetProjection(epsg_code)
del ds
del preprocess_config["stack_bands"]
else:
logger.warn("metadata.xml not found. No way to fix .bil without geotransform")