#5209 closed defect (fixed)
Bug in 'cubic' resampling of gdalwarp
Reported by: | sprice | Owned by: | warmerdam |
---|---|---|---|
Priority: | normal | Milestone: | |
Component: | Algorithms | Version: | 1.10.0 |
Severity: | normal | Keywords: | warper |
Cc: | dron@… |
Description
I get artifacting with the attached file and the command: gdalwarp -tr 4 -4 -r cubic dem_p108r036.sm.tif dem_p108r036.sm-lg.tif
The coefficients seem to be wrong in the cubic interpolation. I done fixed them. To keep total multiplication operations down, I did some pre-multiplication of the 'dfDelta' arguments by 0.5.
I also removed the floor() operations because they are unneeded for all pixel indices greater than 0. Edge indices (less than and equal to zero) are already handled by GWKBilinearResample().
For correctness, I also added 0.5 to the interpolated value of shorts when casting.
And, just because, I pre-saved pointers so gratuitous pointer following isn't happening.
Attachments (2)
Change History (10)
by , 11 years ago
Attachment: | gdalwarpkernel.diff added |
---|
by , 11 years ago
Attachment: | dem_p108r036.sm.tif.gz added |
---|
comment:1 by , 11 years ago
Cc: | added |
---|
comment:3 by , 11 years ago
Component: | default → Algorithms |
---|---|
Keywords: | warper added |
comment:5 by , 11 years ago
Resolution: | → fixed |
---|---|
Status: | assigned → closed |
I've applied the floor and +0.5 changes. Other optimizations should really be held till there are more complete tests.
comment:6 by , 11 years ago
One fix to your patch:
- the line "*piValue = 32768;" should probably be "*piValue = -32768;"
Moving the dfDelta* and memory pointer invariants outside of the loop results in an overall 2% speedup in my tests. I'm running with Landsat 8 panchromatic and resampling to 3 m with a command like:
time gdalwarp -tr 3 -3 LC80210332013169LGN00_B8.TIF LC80210332013169LGN00_B8.spatch.resamp5.TIF
If one was resampling all of Landsat 8's panchromatic bands in this manner on my machine, the processor is spending an extra two days simply recomputing values. (for example)
comment:7 by , 11 years ago
Also, my tests were in Landsat 8 (short precision) and regular usage is in byte precision. That will have half as much data per calculation. Less time will be spent on I/O and correspondingly more time on the kernel operations. The optimization probably be slightly more significant with bytes.
comment:8 by , 10 years ago
The change "*piValue = 32768;" --> "*piValue = -32768;" was actually done in r26504
Hi Andrey,
SVN history shows that you are the original author of the bicubic interpolator (r5484, 10 years ago !). Seth's patch seems reasonable as it adopts the coefficients mentionned in http://en.wikipedia.org/wiki/Bicubic_interpolation#Bicubic_convolution_algorithm , which are also the ones mentionned in formula (18) in the Key's paper http://verona.fi-p.unam.mx/boris/practicas/CubConvInterp.pdf . Do you remember where the coefficients you proposed in the CubicConvolution() macro come from ?