[GRASS-SVN] r44062 - grass/trunk/raster/r.proj
svn_grass at osgeo.org
svn_grass at osgeo.org
Thu Oct 28 03:16:24 EDT 2010
Author: mmetz
Date: 2010-10-28 00:16:24 -0700 (Thu, 28 Oct 2010)
New Revision: 44062
Added:
grass/trunk/raster/r.proj/lanczos.c
Modified:
grass/trunk/raster/r.proj/main.c
grass/trunk/raster/r.proj/r.proj.h
grass/trunk/raster/r.proj/r.proj.html
Log:
add lanczos and lanczos fallback to r.proj
Added: grass/trunk/raster/r.proj/lanczos.c
===================================================================
--- grass/trunk/raster/r.proj/lanczos.c (rev 0)
+++ grass/trunk/raster/r.proj/lanczos.c 2010-10-28 07:16:24 UTC (rev 44062)
@@ -0,0 +1,109 @@
+/*
+ * Name
+ * lanczos.c -- use lanczos interpolation for given row, col
+ *
+ * Description
+ * lanczos returns the value in the buffer that is the result of
+ * lanczos interpolation for the given row, column indices.
+ * If the given row or column is outside the bounds of the input map,
+ * the corresponding point in the output map is set to NULL.
+ *
+ * If single surrounding points in the interpolation matrix are
+ * NULL they where filled with their neighbor
+ */
+
+#include <grass/gis.h>
+#include <grass/raster.h>
+#include <math.h>
+#include "r.proj.h"
+
+void p_lanczos(struct cache *ibuffer, /* input buffer */
+ void *obufptr, /* ptr in output buffer */
+ int cell_type, /* raster map type of obufptr */
+ double *col_idx, /* column index (decimal) */
+ double *row_idx, /* row index (decimal) */
+ struct Cell_head *cellhd /* information of output map */
+ )
+{
+ int row; /* row indices for interp */
+ int col; /* column indices for interp */
+ int i, j, k;
+ double t, u; /* intermediate slope */
+ FCELL result; /* result of interpolation */
+ DCELL cell[25];
+
+ /* cut indices to integer */
+ row = (int)floor(*row_idx);
+ col = (int)floor(*col_idx);
+
+ /* check for out of bounds of map - if out of bounds set NULL value */
+ if (row - 2 < 0 || row + 2 >= cellhd->rows ||
+ col - 2 < 0 || col + 2 >= cellhd->cols) {
+ Rast_set_null_value(obufptr, 1, cell_type);
+ return;
+ }
+
+ k = 0;
+ for (i = 0; i < 5; i++) {
+ for (j = 0; j < 5; j++) {
+ const FCELL *cellp = CPTR(ibuffer, row - 2 + i, col - 2 + j);
+ if (Rast_is_f_null_value(cellp)) {
+ Rast_set_null_value(obufptr, 1, cell_type);
+ return;
+ }
+ cell[k++] = *cellp;
+ }
+ }
+
+ /* do the interpolation */
+ t = *col_idx - 0.5 - col;
+ u = *row_idx - 0.5 - row;
+
+ result = Rast_interp_lanczos(u, t, cell);
+
+ Rast_set_f_value(obufptr, result, cell_type);
+}
+
+void p_lanczos_f(struct cache *ibuffer, /* input buffer */
+ void *obufptr, /* ptr in output buffer */
+ int cell_type, /* raster map type of obufptr */
+ double *col_idx, /* column index (decimal) */
+ double *row_idx, /* row index (decimal) */
+ struct Cell_head *cellhd /* information of output map */
+ )
+{
+ int row; /* row indices for interp */
+ int col; /* column indices for interp */
+ FCELL *cellp, cell;
+
+ /* cut indices to integer */
+ row = (int)floor(*row_idx);
+ col = (int)floor(*col_idx);
+
+ /* check for out of bounds - if out of bounds set NULL value */
+ if (row < 0 || row >= cellhd->rows || col < 0 || col >= cellhd->cols) {
+ Rast_set_null_value(obufptr, 1, cell_type);
+ return;
+ }
+
+ cellp = CPTR(ibuffer, row, col);
+ /* if nearest is null, all the other interps will be null */
+ if (Rast_is_f_null_value(cellp)) {
+ Rast_set_null_value(obufptr, 1, cell_type);
+ return;
+ }
+ cell = *cellp;
+
+ p_lanczos(ibuffer, obufptr, cell_type, col_idx, row_idx, cellhd);
+ /* fallback to bicubic if lanczos is null */
+ if (Rast_is_f_null_value(obufptr)) {
+ p_cubic(ibuffer, obufptr, cell_type, col_idx, row_idx, cellhd);
+ /* fallback to bilinear if cubic is null */
+ if (Rast_is_f_null_value(obufptr)) {
+ p_bilinear(ibuffer, obufptr, cell_type, col_idx, row_idx, cellhd);
+ /* fallback to nearest if bilinear is null */
+ if (Rast_is_f_null_value(obufptr))
+ Rast_set_f_value(obufptr, cell, cell_type);
+ }
+ }
+}
Modified: grass/trunk/raster/r.proj/main.c
===================================================================
--- grass/trunk/raster/r.proj/main.c 2010-10-28 07:11:32 UTC (rev 44061)
+++ grass/trunk/raster/r.proj/main.c 2010-10-28 07:16:24 UTC (rev 44062)
@@ -49,6 +49,7 @@
* Glynn Clements 2006: Use G_interp_* functions, modified
* version of r.proj which uses a tile cache instead of loading the
* entire map into memory.
+* Markus Metz 2010: lanczos and lanczos fallback interpolation methods
*****************************************************************************/
@@ -68,8 +69,10 @@
{p_nearest, "nearest", "nearest neighbor"},
{p_bilinear, "bilinear", "bilinear"},
{p_cubic, "cubic", "cubic convolution"},
+ {p_lanczos, "lanczos", "lanczos filter"},
{p_bilinear_f, "bilinear_f", "bilinear with fallback"},
{p_cubic_f, "cubic_f", "cubic convolution with fallback"},
+ {p_lanczos_f, "lanczos_f", "lanczos filter with fallback"},
{NULL, NULL, NULL}
};
Modified: grass/trunk/raster/r.proj/r.proj.h
===================================================================
--- grass/trunk/raster/r.proj/r.proj.h 2010-10-28 07:11:32 UTC (rev 44061)
+++ grass/trunk/raster/r.proj/r.proj.h 2010-10-28 07:16:24 UTC (rev 44062)
@@ -55,6 +55,11 @@
/* cubic_f.c */
extern void p_cubic_f(struct cache *, void *, int, double *, double *,
struct Cell_head *);
+/* lanczos.c */
+extern void p_lanczos(struct cache *, void *, int, double *, double *,
+ struct Cell_head *);
+extern void p_lanczos_f(struct cache *, void *, int, double *, double *,
+ struct Cell_head *);
#if 1
Modified: grass/trunk/raster/r.proj/r.proj.html
===================================================================
--- grass/trunk/raster/r.proj/r.proj.html 2010-10-28 07:11:32 UTC (rev 44061)
+++ grass/trunk/raster/r.proj/r.proj.html 2010-10-28 07:16:24 UTC (rev 44062)
@@ -94,8 +94,8 @@
map from a different location, projects it and write it out to the current
location.
<br>
-The projected data is resampled with one of three different methods:
-nearest neighbor, bilinear and cubic convolution.
+The projected data is resampled with one of four different methods:
+nearest neighbor, bilinear, cubic convolution or lanczos.
<p>
The <em>method=nearest</em> method, which performs a nearest neighbor assignment,
is the fastest of the three resampling methods. It is primarily used for
@@ -104,19 +104,24 @@
value of the cell based on a weighted distance average of the 4 surrounding
cells in the input map. The <em>method=cubic</em> method determines the new value of
the cell based on a weighted distance average of the 16 surrounding cells in
-the input map.
+the input map. The <em>method=lanzcos</em> method determines the new value of
+the cell based on a weighted distance average of the 25 surrounding cells in
+the input map. Compared to cubic, lanczos puts a higher weight on cells close
+to the center and a lower weight on cells away from the center, resulting in
+slightly better contrast.
<p>
-The bilinear and cubic interpolation methods are most appropriate for
-continuous data and cause some smoothing. Both options should not be used
+The bilinear, cubic and lanczos interpolation methods are most appropriate for
+continuous data and cause some smoothing. The amount of smoothing decreases from
+bilinear to cubic to lanczos. These options should not be used
with categorical data, since the cell values will be altered.
<p>
-In the bilinear and cubic methods, if any of the surrounding cells used to
+In the bilinear, cubic and lanczos methods, if any of the surrounding cells used to
interpolate the new cell value are null, the resulting cell will be null, even if
the nearest cell is not null. This will cause some thinning along null borders,
-such as the coasts of land areas in a DEM. The bilinear_f and cubic_f
+such as the coasts of land areas in a DEM. The bilinear_f, cubic_f and lanczos_f
interpolation methods can be used if thinning along null edges is not desired.
These methods "fall back" to simpler interpolation methods along null borders.
-That is, from cubic to bilinear to nearest.
+That is, from lanczos to cubic to bilinear to nearest.
<p>
If nearest neighbor assignment is used, the output map has the same raster
format as the input map. If any of the interpolations is used, the
@@ -155,7 +160,7 @@
<tt>g.region -a res=5 -p</tt>. Note that this is just a rough guide.
<p>
-A more involvedi, but more accurate, way to do this is to generate a vector
+A more involved, but more accurate, way to do this is to generate a vector
"box" map of the region in the source location using
<em><a href="v.in.region.html">v.in.region</a></em>.
This "box" map is then reprojected into the target location with
More information about the grass-commit
mailing list