I am trying to rectify a tif of of sonar images collected from a moving ship. The tifs are in the form of long ribbons. I have ground control points (GCPs) that match coordinates on the tif (row and column of a particular pixel) to the coordinates of the ship track. Every fifth row of pixels in the image has one GCP
Examples of the GCP points:
(1904.00P ,23086.00L) -> (-67.1203008E,42.0984950N,0.00)
(0.00P ,23081.00L) -> (-67.1203325E,42.0993950N,0.00)
(1904.00P ,23076.00L) -> (-67.1202966E,42.0984951N,0.00)
(0.00P ,23066.00L) -> (-67.1202984E,42.0993966N,0.00)
The tifs themselves have no GCS and are completely rectangular. However, the track to which I am trying to georeference them is curved, occasionally turns 90 degrees, and in some cases will double back on itself.
What i've been trying to do is the following (adapted from this stack overflow thread) :
# open tif ds = gdal.Open(sonar_tif, gdal.GA_Update) # Get Source Projections (PCS) ds_prj = ds.GetProjection() if ds_prj: print("Source Projection:",ds_prj,"\n") else: print("Source Projection: None\n") # new instance of SpatialReference class from osr module drs = osr.SpatialReference() print("drs is",drs) # Set Geographic Coordinate Reference System (GCS) srs = osr.SpatialReference() srs.ImportFromEPSG(4326) # WGS84 GCS if srs.IsProjected: print ("Target Projection:",srs.GetAttrValue('projcs')) print ("Target GCS:",srs.GetAttrValue('geogcs')) # Attach GCP points to tif ds.SetGCPs(gcp_list, ds_prj) ds = None # closes dataset, unlocks files d_out = None
If I understand correctly, this should perform the same action as gdal_translate, and I should have embedded GCP points and a geotiff. However, the output I get is a big black square.
Am I doing something fundamentally incorrect?
Is there a coordinate system mismatch that I might be unaware of?
Do I have to run gdal_warp before I can see a fully georeferenced image?