[gdal-dev] GDAL Cubic Interpolation

Previous Topic Next Topic
 
classic Classic list List threaded Threaded
5 messages Options
Reply | Threaded
Open this post in threaded view
|

[gdal-dev] GDAL Cubic Interpolation

Harvey Berger
Hi,

I subsampled a GLOBE image of elevation data by a subsampling factor of 12 to a resolution of 0.1 deg x 0.1 deg on integer multiples of 0.1 deg. The original GLOBE image is: Latitude from +90-1/240 deg to -90+1/240 deg in 1/120 deg steps, and Longitude from -180+1/240 deg to +180-1/240 deg in 1/120 deg steps. 

I understand that the cubic interpolation method in gdalwarp accounts for the subsampling factor, so 48x48 adjacent pixels influence the target pixel for a subsampling factor of 12. 

Does gdalwarp assume the original image is a sphere, where 48x48 adjacent pixels influence all target pixels, or a rectangle, where a limited number of pixels influence the target pixels at the top, bottom, left, and right edges of the subsampled image? If gdalwarp assumes the original image is a rectangle and uses a modified interpolation method at the four edges of the subsampled image, the modified method affects two rows and two columns at the top, bottom, left and right edges of the subsampled image. If gdalwarp assumes the original image is a rectangle, is there a command option so gdalwarp considers spherical wrapping when interpolating? Alternatively, I assume I could augment the original GLOBE image with the appropriate additional rows and columns on all four edges. If I did that, how do I tell gdalwarp the extent of the augmented input image?

Kind regards,

Harvey

_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|

Re: GDAL Cubic Interpolation

Even Rouault-2

Harvey,

 

>

> I subsampled a GLOBE image of elevation data by a subsampling factor of 12

> to a resolution of 0.1 deg x 0.1 deg on integer multiples of 0.1 deg. The

> original GLOBE image is: Latitude from +90-1/240 deg to -90+1/240 deg in

> 1/120 deg steps, and Longitude from -180+1/240 deg to +180-1/240 deg in

> 1/120 deg steps.

>

> I understand that the cubic interpolation method in gdalwarp accounts for

> the subsampling factor, so 48x48 adjacent pixels influence the target pixel

> for a subsampling factor of 12.

>

> Does gdalwarp assume the original image is a sphere, where 48x48 adjacent

> pixels influence all target pixels, or a rectangle, where a limited number

> of pixels influence the target pixels at the top, bottom, left, and right

> edges of the subsampled image?

 

A rectangle, with fallback to bilinear resampling at edges

 

> If gdalwarp assumes the original image is a

> rectangle and uses a modified interpolation method at the four edges of the

> subsampled image, the modified method affects two rows and two columns at

> the top, bottom, left and right edges of the subsampled image. If gdalwarp

> assumes the original image is a rectangle, is there a command option so

> gdalwarp considers spherical wrapping when interpolating?

 

Nothing directly straightforward. Your below plan is the a reasonable solution.

 

> Alternatively, I

> assume I could augment the original GLOBE image with the appropriate

> additional rows and columns on all four edges. If I did that, how do I tell

> gdalwarp the extent of the augmented input image?

 

You'll have to recreate a properly referenced input image, for example with

 

gdal_translate -ullr upperleftx upperlefty lowerleftx lowerlefty

 

Even

 

--

Spatialys - Geospatial professional services

http://www.spatialys.com


_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|

Re: GDAL Cubic Interpolation

Joaquim Luis
Alternatively,  use GMT's grdsample with -fg option to force the knowledge (if it's not already in the file) that the Earth is round.

http://gmt.soest.hawaii.edu/doc/5.4.2/grdsample.html



Harvey,

 

>

> I subsampled a GLOBE image of elevation data by a subsampling factor of 12

> to a resolution of 0.1 deg x 0.1 deg on integer multiples of 0.1 deg. The

> original GLOBE image is: Latitude from +90-1/240 deg to -90+1/240 deg in

> 1/120 deg steps, and Longitude from -180+1/240 deg to +180-1/240 deg in

> 1/120 deg steps.

>

> I understand that the cubic interpolation method in gdalwarp accounts for

> the subsampling factor, so 48x48 adjacent pixels influence the target pixel

> for a subsampling factor of 12.

>

> Does gdalwarp assume the original image is a sphere, where 48x48 adjacent

> pixels influence all target pixels, or a rectangle, where a limited number

> of pixels influence the target pixels at the top, bottom, left, and right

> edges of the subsampled image?

 

A rectangle, with fallback to bilinear resampling at edges

 

> If gdalwarp assumes the original image is a

> rectangle and uses a modified interpolation method at the four edges of the

> subsampled image, the modified method affects two rows and two columns at

> the top, bottom, left and right edges of the subsampled image. If gdalwarp

> assumes the original image is a rectangle, is there a command option so

> gdalwarp considers spherical wrapping when interpolating?

 

Nothing directly straightforward. Your below plan is the a reasonable solution.

 

> Alternatively, I

> assume I could augment the original GLOBE image with the appropriate

> additional rows and columns on all four edges. If I did that, how do I tell

> gdalwarp the extent of the augmented input image?

 

You'll have to recreate a properly referenced input image, for example with

 

gdal_translate -ullr upperleftx upperlefty lowerleftx lowerlefty

 

Even

 

--

Spatialys - Geospatial professional services

http://www.spatialys.com





_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|

Re: GDAL Cubic Interpolation

Harvey Berger
Thanks Joaquim and Even,

Per my email I augmented the original GLOBE image and created a properly referenced image using gdal_translate as follows: I did an imread from QGIS/GDAL into Octave, added the appropriate 24 rows and columns on each side, and did an imwrite back to QGIS/GDAL to continue processing; however, I'm having a problem with negative elevations probably due to data type confusion.

1. GDAL reports the data type of the original GLOBE tiff file as int16. (file = output1.tiff)
2. After using imread in Octave, Octave reports the data type as uint16.
3. I retained the uint16 data type in Octave and added the additional rows and columns.
4. I did an imwrite to create the augmented file. (file=output2.tiff)
5. GDAL reports the data type of the augmented files as uint16.
6. I performed the following GDAL command: gdal_translate -a_ullr ulx uly lrx lry src_dataset dst_dataset. (file=output3.tiff)
7. I then used gdalwarp to create the subsampled and cubic interpolated file (file=sample.tiff); GDAL report the data type as int16.

The problem is with locations with negative elevations that result in erroneous results when I subsequent subsample using gdal_warp. Here are the elevations reported by GDAL for one example pixel in the Caspian Sea:
output1.tiff (int16): -28 m
output2.tiff (uint16): 65508 m
output3.tiff (uint16): 65508 m

Unfortunately, Octave will not perform an imwrite with a int16 data type. Is it possible to fix this within QGIS/GDAL or is there another solution?

Thanks again for your help.

Kind regards,

Harvey

On Mon, Oct 9, 2017 at 9:32 AM, Joaquim Luis <[hidden email]> wrote:
Alternatively,  use GMT's grdsample with -fg option to force the knowledge (if it's not already in the file) that the Earth is round.




Harvey,

 

>

> I subsampled a GLOBE image of elevation data by a subsampling factor of 12

> to a resolution of 0.1 deg x 0.1 deg on integer multiples of 0.1 deg. The

> original GLOBE image is: Latitude from +90-1/240 deg to -90+1/240 deg in

> 1/120 deg steps, and Longitude from -180+1/240 deg to +180-1/240 deg in

> 1/120 deg steps.

>

> I understand that the cubic interpolation method in gdalwarp accounts for

> the subsampling factor, so 48x48 adjacent pixels influence the target pixel

> for a subsampling factor of 12.

>

> Does gdalwarp assume the original image is a sphere, where 48x48 adjacent

> pixels influence all target pixels, or a rectangle, where a limited number

> of pixels influence the target pixels at the top, bottom, left, and right

> edges of the subsampled image?

 

A rectangle, with fallback to bilinear resampling at edges

 

> If gdalwarp assumes the original image is a

> rectangle and uses a modified interpolation method at the four edges of the

> subsampled image, the modified method affects two rows and two columns at

> the top, bottom, left and right edges of the subsampled image. If gdalwarp

> assumes the original image is a rectangle, is there a command option so

> gdalwarp considers spherical wrapping when interpolating?

 

Nothing directly straightforward. Your below plan is the a reasonable solution.

 

> Alternatively, I

> assume I could augment the original GLOBE image with the appropriate

> additional rows and columns on all four edges. If I did that, how do I tell

> gdalwarp the extent of the augmented input image?

 

You'll have to recreate a properly referenced input image, for example with

 

gdal_translate -ullr upperleftx upperlefty lowerleftx lowerlefty

 

Even

 

--

Spatialys - Geospatial professional services

http://www.spatialys.com





_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev


_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|

Re: GDAL Cubic Interpolation

Harvey Berger
To follow-up on my latest email, I confirmed that GDAL is interpreting the pixels in the augmented image file from Octave as uint16 rather than int16, and this results in incorrect elevation values after cubic interpolation and subsampling (e.g., a pixel elevation of -28 m is misinterpreted as 65,508 m). Octave has a typecast function that changes a data type without changing the data. Does GDAL have an equivalent function .... or does some field need to be changed in the augmented tiff image file?

Kind regards,

Harvey

On Wed, Oct 11, 2017 at 10:06 AM, Harvey Berger <[hidden email]> wrote:
Thanks Joaquim and Even,

Per my email I augmented the original GLOBE image and created a properly referenced image using gdal_translate as follows: I did an imread from QGIS/GDAL into Octave, added the appropriate 24 rows and columns on each side, and did an imwrite back to QGIS/GDAL to continue processing; however, I'm having a problem with negative elevations probably due to data type confusion.

1. GDAL reports the data type of the original GLOBE tiff file as int16. (file = output1.tiff)
2. After using imread in Octave, Octave reports the data type as uint16.
3. I retained the uint16 data type in Octave and added the additional rows and columns.
4. I did an imwrite to create the augmented file. (file=output2.tiff)
5. GDAL reports the data type of the augmented files as uint16.
6. I performed the following GDAL command: gdal_translate -a_ullr ulx uly lrx lry src_dataset dst_dataset. (file=output3.tiff)
7. I then used gdalwarp to create the subsampled and cubic interpolated file (file=sample.tiff); GDAL report the data type as int16.

The problem is with locations with negative elevations that result in erroneous results when I subsequent subsample using gdal_warp. Here are the elevations reported by GDAL for one example pixel in the Caspian Sea:
output1.tiff (int16): -28 m
output2.tiff (uint16): 65508 m
output3.tiff (uint16): 65508 m

Unfortunately, Octave will not perform an imwrite with a int16 data type. Is it possible to fix this within QGIS/GDAL or is there another solution?

Thanks again for your help.

Kind regards,

Harvey

On Mon, Oct 9, 2017 at 9:32 AM, Joaquim Luis <[hidden email]> wrote:
Alternatively,  use GMT's grdsample with -fg option to force the knowledge (if it's not already in the file) that the Earth is round.




Harvey,

 

>

> I subsampled a GLOBE image of elevation data by a subsampling factor of 12

> to a resolution of 0.1 deg x 0.1 deg on integer multiples of 0.1 deg. The

> original GLOBE image is: Latitude from +90-1/240 deg to -90+1/240 deg in

> 1/120 deg steps, and Longitude from -180+1/240 deg to +180-1/240 deg in

> 1/120 deg steps.

>

> I understand that the cubic interpolation method in gdalwarp accounts for

> the subsampling factor, so 48x48 adjacent pixels influence the target pixel

> for a subsampling factor of 12.

>

> Does gdalwarp assume the original image is a sphere, where 48x48 adjacent

> pixels influence all target pixels, or a rectangle, where a limited number

> of pixels influence the target pixels at the top, bottom, left, and right

> edges of the subsampled image?

 

A rectangle, with fallback to bilinear resampling at edges

 

> If gdalwarp assumes the original image is a

> rectangle and uses a modified interpolation method at the four edges of the

> subsampled image, the modified method affects two rows and two columns at

> the top, bottom, left and right edges of the subsampled image. If gdalwarp

> assumes the original image is a rectangle, is there a command option so

> gdalwarp considers spherical wrapping when interpolating?

 

Nothing directly straightforward. Your below plan is the a reasonable solution.

 

> Alternatively, I

> assume I could augment the original GLOBE image with the appropriate

> additional rows and columns on all four edges. If I did that, how do I tell

> gdalwarp the extent of the augmented input image?

 

You'll have to recreate a properly referenced input image, for example with

 

gdal_translate -ullr upperleftx upperlefty lowerleftx lowerlefty

 

Even

 

--

Spatialys - Geospatial professional services

http://www.spatialys.com





_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev



_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev