0

I am trying to convert a pixel value above 0.02 to be nodata. However, the output file comes out with no projection. When I assign projection, QGIS loads it in the wrong place.

gdal_calc.py -A "000-IMG_0037_4_gcp.tif" --outfile="000-IMG_0037_4_gcp_noglint8.tif" --calc="((A<=0.02)*(A>0))*A+((A>0.02)*0)" --NoDataValue=0

The gdalinfo of the original file:

Files: /drive/000-IMG_0037_4_gcp.tif
Size is 1048, 1422
Coordinate System is `'
GCP Projection = 
PROJCS["WGS 84 / UTM zone 16N",
    GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0,
        AUTHORITY["EPSG","8901"]],
    UNIT["degree",0.0174532925199433,
        AUTHORITY["EPSG","9122"]],
    AUTHORITY["EPSG","4326"]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",-87],
    PARAMETER["scale_factor",0.9996],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
    AUTHORITY["EPSG","9001"]],
    AXIS["Easting",EAST],
    AXIS["Northing",NORTH],
    AUTHORITY["EPSG","32616"]]
GCP[  0]: Id=1, Info=
      (0,0) -> (280191.450981,3352631.04297,0)
GCP[  1]: Id=2, Info=
      (1048.9984,0) -> (280298.687984,3352631.04297,0)
GCP[  2]: Id=3, Info=
      (1048.9984,1422.4381) -> (280298.687984,3352484.8363,0)
GCP[  3]: Id=4, Info=
      (0,0) -> (280191.450981,3352631.04297,0)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0, 1422.0)
Upper Right ( 1048.0,    0.0)
Lower Right ( 1048.0, 1422.0)
Center      (  524.0,  711.0)
Band 1 Block=1048x1 Type=Float32, ColorInterp=Gray

The output file gdalinfo:

Driver: GTiff/GeoTIFF
Files: /drive/000-IMG_0037_4_gcp_noglint8.tif
Size is 1048, 1422
Coordinate System is `'
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0, 1422.0)
Upper Right ( 1048.0,    0.0)
Lower Right ( 1048.0, 1422.0)
Center      (  524.0,  711.0)
Band 1 Block=1048x1 Type=Float32, ColorInterp=Gray
  NoData Value=0

Is there a way to make sure that the output image comes with the same srs of the input raster?

1

Your original file isn't actually georeferenced, it uses GCPs. I'm not sure why gdal_calc ignores the GCPs, but you can workaround using gdalwarp to georeference the original file.

The following should work:

gdalwarp 000-IMG_0037_4_gcp.tif 000-IMG_0037_4_warped.tif
gdal_calc.py -A "000-IMG_0037_4_warped.tif" --outfile="000-IMG_0037_4_warped_noglint8.tif" --calc="((A<=0.02)*(A>0))*A+((A>0.02)*0)" --NoDataValue=0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.