1

Say a single-band 16-bit (UInt) raster with values representing elevation (meters) must be converted to a 32-bit RGBA PNG image (8-bits/band) for use in a web application. The image dimensions must be identical so there is no interpolation etc.

1. Can the original 16bit elevation value be back-calculated (algebraically) from the 32bit RGBA values if the image is rendered using color ramps?

2. How could the original image be converted into a 32-bit RGBA image (using R and G bands in QGIS for example), such that the original 16bit precision persists? How could the original elevation be back calculated algebraically from the utilized bands?

For context - certain web applications built on platforms such as OpenLayers3 have limitations on image format support, as do some browsers. For this reason, though not ideal, it can be necessary to back-calculate values from an alternative image format so that meaningful information can be conveyed to the user.

See this SE question for an example of the OpenLayers3 pixel query methodology.

  • 2
    32 bit RGBA isn't greyscale, it's Red, Green, Blue and Alpha - 4 bands of an unsigned byte each, that the values in the first 3 bands have identical values is superfluous, in rendering data has been lost by stretching a value into a lesser scale. I've seen this question in many forms and I will give you the same advice as all the others: you can't, don't even try. Contact the data custodian and request a clip of the DEM. – Michael Stimson Feb 11 at 4:45
  • 'Greyscale' is intended as an adjective implying that the image has no saturation. – CrystallineEntity Feb 12 at 0:02
  • I edited the question adding a followup and a note on the context of the question, since the response has largely focused on why one should not be doing this. I'm interested in the method of doing so. There are additional concerns such as gamma correction and palletting, which I am able to ignore for the purposes of this question. Thank you all for your responses. – CrystallineEntity Feb 12 at 0:11
  • 1
    You may need to re-think what you're doing. Data that has been back transformed this way will produce a result, but that result will not be "...meaningful information can be conveyed to the user". – jberrio Feb 12 at 1:18
1

You could do it if you have at least two points with known elevations, assuming the greyscale ramps was applied linearly. The two points could be for example the maximum and minimum values of the DEM, corresponding to the 0 and 255 greyscale values.

Using them, you need to calculate a linear equation to transform your data; that is an equation with the form

Y = mX + b

To calculate the equation, google how to do it based on two points (for example https://www.mathsisfun.com/algebra/line-equation-2points.html). Your two know points will be something like (0, min_elev) and (255, max_elev).

After you calculate m and b you can use this equation in QGIS' raster calculator, being able to back-calculate the elevations.

Be aware that even in the best of cases, some resolution will be lost due to rounding and decimal imprecisions. In the worst of cases, the grayscale was not applied linearly as you think and the equation will not be valid.

EDIT: (After considering Michael Stimson's comment) If you do it, you have no way to confirm if the transformation was applied correctly. You should not rely on the back-transformed data for any serious purpose. If anything, you can use the process just as a classroom exercise to learn how to use the raster calculator.

  • Am I to understand then, that in rendering the 16bit image into a 32bit image the range of possible elevation values has actually been reduced from 65,536, to 256? – CrystallineEntity Feb 11 at 2:23
  • 1
    I don't know what they were before but the RGB image now only has 256 possible values in a band. If all the bands are effectively a copy of each other, you could discard the duplicated bands and that leaves you with one band with those 256 possible values (with no decimals in between). – jberrio Feb 11 at 4:13
  • 1
    Even this may not work if the stretch that was used is standard deviations, custom or some other sort of algebraic stretch. There is no guarantee that the correlation is linear; If you want the DEM ask the organization that generated the raster rather than trying to estimate or interpolate which can be dangerous if the DEM is later relied upon. – Michael Stimson Feb 11 at 4:37
  • @Michael Stimson: Insofar as you trust QGIS to render linearly given the absolute minimum and maximum of the input raster, the image was rendered using a linear ramp. – CrystallineEntity Feb 12 at 0:16
  • The RGBA raster may not have been rendered using QGIS. In ArcGIS there are more options of rendering which are not linear and colour ramps that have multivariate ranges. Any resulting DEM interpolated would be highly suspect, an end user relying on the accuracy of such a DEM could make serious mistakes which can have legal consequences. This is why I say not to even attempt it. Have a go at this solution by all means but heavily disclaimer the resulting product. – Michael Stimson Feb 12 at 0:54
0

Disclosure this is an answer to my own question

Question 1. "Can the original 16bit elevation value be back-calculated (algebraically) from the 32bit RGBA values if the image is rendered using color ramps?"

Without processing the raster, I'm unsure, but my gut tells me it is possible although simple ramps do not seem sufficient.

Question 2. "How could the original image be converted into a 32-bit RGBA image (using R and G bands in QGIS for example), such that the original 16bit precision persists? How could the original elevation be back-calculated algebraically from the utilized bands?"

To do this using two bands is actually quite simple. Each 8bit band may contain a value from 0-255. Two bands of 8bit ==> one band of 16bit.

Deconstruct the 16bit raster

The built-in QGIS 3.4 raster calculator does not appear to have the functionality needed. GDAL(gdal_calc) has access to Numpy (Python library), and with it we can perform two raster calculations as follows:

  1. From the 16bit elevation raster input we divide by 255. We also must truncate the fractional portion of the resulting number using trunc() in the expression. This value will be encoded to one band of an 32bit RGBA.

    gdal_calc -calc "trunc(A/255)" -format GTiff -type Byte -A dem.tif -A_band 1 -outfile dem_div.tif

  2. We calculate the second band using the modulo function on the same 16bit input raster, which returns the integer remainder of a division. This value will be encoded in a second band of the 32bit RGBA.

    gdal_calc -calc "mod(A,255)" -format GTiff -type Byte -A dem.tif -A_band 1 -outfile dem_mod.tif

Construct the 32bit RGBA

Using GDAL or QGIS's merge raster command we assemble each 8 bit image into 2 different bands of an RGBA output raster.

Back-calculate from R & G bands

To back-calculate an original 16bit value we simply assemble the pieces with the JS or equivalent code:

(pixel[band0] * 255) + pixel[band1]

A few notes:

PNG is a lossless format and widely supported in web browsers, which makes this possible. This conversion will not hold true in a lossy format like JPEG, because the compression has changed the data.

Using raster data in GIS oriented web applications poses a real challenge primarily because of limited image format compatability. Newer webmap platforms are making improvements (Mapbox now allows for GeoTIFF), but limitations still exist in other systems. Consequently, it can be difficult to leverage the use of rasters as we would in primary GIS softwares.

This is a 16bit DEM rendered in greyscale 16bit DEM rendered in greyscale

This is an equivalent 32bit RGBA DEM rendered in red and green 32bit RGBA DEM rendered with red and green

Your Answer

By clicking "Post Your Answer", you acknowledge that you have read our updated terms of service, privacy policy and cookie policy, and that your continued use of the website is subject to these policies.

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