Commit bd64ad38 authored by Gaetano Raffaele's avatar Gaetano Raffaele

ENH: possibility to freely specify output epsg code for co-registration (deals...

ENH: possibility to freely specify output epsg code for co-registration (deals with the L8 southern emisphere projection problem).
parent 83ffad1e
......@@ -22,7 +22,7 @@ def main(argv):
sys.exit("Platform not supported!")
try:
opts, args = getopt.getopt(argv, '', ['band-ref=', 'band-in=', 'step=', 'minstep=', 'minpoints=', 'maxpoints=', 'prec='])
opts, args = getopt.getopt(argv, '', ['band-ref=', 'band-in=', 'step=', 'minstep=', 'minpoints=', 'maxpoints=', 'prec=', 'outsrs='])
except getopt.GetoptError as err:
print str(err)
......@@ -35,6 +35,7 @@ def main(argv):
srs = osr.SpatialReference()
srs.ImportFromWkt(prj)
code = srs.GetAuthorityCode(None)
code_out = code
gsp = str(int(2 * round(max(abs(gt[1]), abs(gt[5])))))
vhrb = 3
......@@ -46,6 +47,7 @@ def main(argv):
enough = 20
prec = 3
maxgcp = -1
out_epsg = None
for opt, val in opts:
if opt == '--band-ref':
......@@ -62,6 +64,8 @@ def main(argv):
maxgcp = int(val)
elif opt == '--prec':
prec = float(val)
elif opt == '--outsrs':
code_out = val
tmphomo = randomword(16) + '.geom'
tmpfld = os.path.split(args[1])[0] + '/' + randomword(16) + '/'
......@@ -115,7 +119,7 @@ def main(argv):
ds.SetProjection('')
ds = None
cmd = ['otbcli_OrthoRectification', '-io.in \"' + args[1] + '?&geom=' + tmphomo + '\"', '-io.out', tmpfld + 'sits_band_coreg.tif', 'uint16', '-map', 'epsg', '-map.epsg.code', code, '-opt.gridspacing', gsp]
cmd = ['otbcli_OrthoRectification', '-io.in \"' + args[1] + '?&geom=' + tmphomo + '\"', '-io.out', tmpfld + 'sits_band_coreg.tif', 'uint16', '-map', 'epsg', '-map.epsg.code', code_out, '-opt.gridspacing', gsp]
cmdex = ' '.join(cmd)
subprocess.call(cmdex, shell=True)
......@@ -170,6 +174,6 @@ def main(argv):
if __name__ == '__main__':
if len(sys.argv) < 2:
sys.exit(
'Usage: python coregister.py [--band-ref <target ref band>] [--band-in <source band>] [--step <geobins initial step def. 256>] [--minstep <geobins minimum step def. 16>] [--minpoints <minimum # of homologous points def. 40>] [--maxpoints <maximum # of GCP points def. -1 unlimited> [--prec <colocalisation precision in pixels] <vhsr scene file> <stack to coregister>')
'Usage: python coregister.py [--band-ref <target ref band>] [--band-in <source band>] [--step <geobins initial step def. 256>] [--minstep <geobins minimum step def. 16>] [--minpoints <minimum # of homologous points def. 40>] [--maxpoints <maximum # of GCP points def. -1 unlimited> [--prec <colocalisation precision in pixels] [--outsrs <EPSG code for output>] <vhsr scene file> <stack to coregister>')
else:
main(sys.argv[1:])
......@@ -218,6 +218,7 @@ def main(argv):
# Define here the target SRS based on VHSR image
dst_srs = checkSRS(seg_src)
dst_srs_code = dst_srs.split(':')[1]
seg_fld = os.path.split(seg_src)[0] + '/../SEGMENTATION'
seg_output = seg_fld + '/segmentation.tif'
......@@ -288,7 +289,7 @@ def main(argv):
if not os.path.exists(ref_S2) or force_coreg:
cmd = ['otbcli_Superimpose', '-inr', refl, '-inm', ms_file, '-out', ref_S2, 'uint16']
subprocess.call(cmd, shell=sh)
cmd = ['python', 'coregister.py', '--band-ref', vhrbnd, '--band-in', bnd, '--maxpoints', str(mxp), ref_S2, refl]
cmd = ['python', 'coregister.py', '--band-ref', vhrbnd, '--band-in', bnd, '--maxpoints', str(mxp), '--outsrs', dst_srs_code, ref_S2, refl]
subprocess.call(cmd, shell=sh)
# Sentinel-2 THEIA
......
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment