[GRASS-SVN] r51141 - grass/trunk/raster/r.proj
svn_grass at osgeo.org
svn_grass at osgeo.org
Thu Mar 22 06:28:16 EDT 2012
Author: mmetz
Date: 2012-03-22 03:28:16 -0700 (Thu, 22 Mar 2012)
New Revision: 51141
Modified:
grass/trunk/raster/r.proj/Makefile
grass/trunk/raster/r.proj/bilinear.c
grass/trunk/raster/r.proj/bilinear_f.c
grass/trunk/raster/r.proj/cubic.c
grass/trunk/raster/r.proj/cubic_f.c
grass/trunk/raster/r.proj/lanczos.c
grass/trunk/raster/r.proj/main.c
grass/trunk/raster/r.proj/nearest.c
grass/trunk/raster/r.proj/r.proj.h
grass/trunk/raster/r.proj/readcell.c
Log:
r.proj parallelization with openmp
Modified: grass/trunk/raster/r.proj/Makefile
===================================================================
--- grass/trunk/raster/r.proj/Makefile 2012-03-21 23:59:57 UTC (rev 51140)
+++ grass/trunk/raster/r.proj/Makefile 2012-03-22 10:28:16 UTC (rev 51141)
@@ -4,6 +4,9 @@
LIBES = $(GPROJLIB) $(RASTERLIB) $(GISLIB) $(MATHLIB)
DEPENDENCIES = $(GPROJDEP) $(RASTERDEP) $(GISDEP)
+
+EXTRA_LIBS = $(OMPLIB)
+EXTRA_CFLAGS = $(OMPCFLAGS)
EXTRA_INC = $(PROJINC)
include $(MODULE_TOPDIR)/include/Make/Module.make
Modified: grass/trunk/raster/r.proj/bilinear.c
===================================================================
--- grass/trunk/raster/r.proj/bilinear.c 2012-03-21 23:59:57 UTC (rev 51140)
+++ grass/trunk/raster/r.proj/bilinear.c 2012-03-22 10:28:16 UTC (rev 51141)
@@ -18,8 +18,8 @@
void p_bilinear(struct cache *ibuffer, /* input buffer */
void *obufptr, /* ptr in output buffer */
int cell_type, /* raster map type of obufptr */
- double *col_idx, /* column index */
- double *row_idx, /* row index */
+ double col_idx, /* column index */
+ double row_idx, /* row index */
struct Cell_head *cellhd /* information of output map */
)
{
@@ -31,8 +31,8 @@
FCELL c[2][2];
/* cut indices to integer */
- row = (int)floor(*row_idx - 0.5);
- col = (int)floor(*col_idx - 0.5);
+ row = (int)floor(row_idx - 0.5);
+ col = (int)floor(col_idx - 0.5);
/* check for out of bounds - if out of bounds set NULL value and return */
if (row < 0 || row + 1 >= cellhd->rows || col < 0 || col + 1 >= cellhd->cols) {
@@ -42,17 +42,17 @@
for (i = 0; i < 2; i++)
for (j = 0; j < 2; j++) {
- const FCELL *cellp = CPTR(ibuffer, row + i, col + j);
- if (Rast_is_f_null_value(cellp)) {
+ const FCELL cell = CVAL(ibuffer, row + i, col + j);
+ if (Rast_is_f_null_value(&cell)) {
Rast_set_null_value(obufptr, 1, cell_type);
return;
}
- c[i][j] = *cellp;
+ c[i][j] = cell;
}
/* do the interpolation */
- t = *col_idx - 0.5 - col;
- u = *row_idx - 0.5 - row;
+ t = col_idx - 0.5 - col;
+ u = row_idx - 0.5 - row;
result = Rast_interp_bilinear(t, u, c[0][0], c[0][1], c[1][0], c[1][1]);
Modified: grass/trunk/raster/r.proj/bilinear_f.c
===================================================================
--- grass/trunk/raster/r.proj/bilinear_f.c 2012-03-21 23:59:57 UTC (rev 51140)
+++ grass/trunk/raster/r.proj/bilinear_f.c 2012-03-22 10:28:16 UTC (rev 51141)
@@ -15,18 +15,18 @@
void p_bilinear_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 */
- double *row_idx, /* row index */
+ double col_idx, /* column index */
+ double row_idx, /* row index */
struct Cell_head *cellhd /* cell header of input layer */
)
{
/* start nearest neighbor to do some basic tests */
int row, col; /* row/col of nearest neighbor */
- FCELL *cellp, cell;
+ FCELL cell;
/* cut indices to integer */
- row = (int)floor(*row_idx);
- col = (int)floor(*col_idx);
+ 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) {
@@ -34,13 +34,12 @@
return;
}
- cellp = CPTR(ibuffer, row, col);
+ cell = CVAL(ibuffer, row, col);
/* if nearest is null, all the other interps will be null */
- if (Rast_is_f_null_value(cellp)) {
+ if (Rast_is_f_null_value(&cell)) {
Rast_set_null_value(obufptr, 1, cell_type);
return;
}
- cell = *cellp;
p_bilinear(ibuffer, obufptr, cell_type, col_idx, row_idx, cellhd);
/* fallback to nearest if bilinear is null */
Modified: grass/trunk/raster/r.proj/cubic.c
===================================================================
--- grass/trunk/raster/r.proj/cubic.c 2012-03-21 23:59:57 UTC (rev 51140)
+++ grass/trunk/raster/r.proj/cubic.c 2012-03-22 10:28:16 UTC (rev 51141)
@@ -21,8 +21,8 @@
void p_cubic(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) */
+ double col_idx, /* column index (decimal) */
+ double row_idx, /* row index (decimal) */
struct Cell_head *cellhd /* information of output map */
)
{
@@ -32,11 +32,11 @@
FCELL t, u; /* intermediate slope */
FCELL result; /* result of interpolation */
FCELL val[4]; /* buffer for temporary values */
- FCELL cell[4][4];
+ FCELL c[4][4];
/* cut indices to integer */
- row = (int)floor(*row_idx - 0.5);
- col = (int)floor(*col_idx - 0.5);
+ row = (int)floor(row_idx - 0.5);
+ col = (int)floor(col_idx - 0.5);
/* check for out of bounds of map - if out of bounds set NULL value */
if (row - 1 < 0 || row + 2 >= cellhd->rows ||
@@ -47,20 +47,20 @@
for (i = 0; i < 4; i++)
for (j = 0; j < 4; j++) {
- const FCELL *cellp = CPTR(ibuffer, row - 1 + i, col - 1 + j);
- if (Rast_is_f_null_value(cellp)) {
+ const FCELL cell = CVAL(ibuffer, row - 1 + i, col - 1 + j);
+ if (Rast_is_f_null_value(&cell)) {
Rast_set_null_value(obufptr, 1, cell_type);
return;
}
- cell[i][j] = *cellp;
+ c[i][j] = cell;
}
/* do the interpolation */
- t = *col_idx - 0.5 - col;
- u = *row_idx - 0.5 - row;
+ t = col_idx - 0.5 - col;
+ u = row_idx - 0.5 - row;
for (i = 0; i < 4; i++) {
- const FCELL *tmp = cell[i];
+ const FCELL *tmp = c[i];
val[i] = Rast_interp_cubic(t, tmp[0], tmp[1], tmp[2], tmp[3]);
}
Modified: grass/trunk/raster/r.proj/cubic_f.c
===================================================================
--- grass/trunk/raster/r.proj/cubic_f.c 2012-03-21 23:59:57 UTC (rev 51140)
+++ grass/trunk/raster/r.proj/cubic_f.c 2012-03-22 10:28:16 UTC (rev 51141)
@@ -16,18 +16,18 @@
void p_cubic_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 */
- double *row_idx, /* row index */
+ double col_idx, /* column index */
+ double row_idx, /* row index */
struct Cell_head *cellhd /* cell header of input layer */
)
{
/* start nearest neighbor to do some basic tests */
int row, col; /* row/col of nearest neighbor */
- FCELL *cellp, cell;
+ FCELL cell;
/* cut indices to integer */
- row = (int)floor(*row_idx);
- col = (int)floor(*col_idx);
+ 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) {
@@ -35,13 +35,12 @@
return;
}
- cellp = CPTR(ibuffer, row, col);
+ cell = CVAL(ibuffer, row, col);
/* if nearest is null, all the other interps will be null */
- if (Rast_is_f_null_value(cellp)) {
+ if (Rast_is_f_null_value(&cell)) {
Rast_set_null_value(obufptr, 1, cell_type);
return;
}
- cell = *cellp;
p_cubic(ibuffer, obufptr, cell_type, col_idx, row_idx, cellhd);
/* fallback to bilinear if cubic is null */
Modified: grass/trunk/raster/r.proj/lanczos.c
===================================================================
--- grass/trunk/raster/r.proj/lanczos.c 2012-03-21 23:59:57 UTC (rev 51140)
+++ grass/trunk/raster/r.proj/lanczos.c 2012-03-22 10:28:16 UTC (rev 51141)
@@ -20,8 +20,8 @@
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) */
+ double col_idx, /* column index (decimal) */
+ double row_idx, /* row index (decimal) */
struct Cell_head *cellhd /* information of output map */
)
{
@@ -30,11 +30,11 @@
int i, j, k;
double t, u; /* intermediate slope */
FCELL result; /* result of interpolation */
- DCELL cell[25];
+ DCELL c[25];
/* cut indices to integer */
- row = (int)floor(*row_idx);
- col = (int)floor(*col_idx);
+ 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 ||
@@ -46,20 +46,20 @@
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)) {
+ const FCELL cell = CVAL(ibuffer, row - 2 + i, col - 2 + j);
+ if (Rast_is_f_null_value(&cell)) {
Rast_set_null_value(obufptr, 1, cell_type);
return;
}
- cell[k++] = *cellp;
+ c[k++] = cell;
}
}
/* do the interpolation */
- t = *col_idx - 0.5 - col;
- u = *row_idx - 0.5 - row;
+ t = col_idx - 0.5 - col;
+ u = row_idx - 0.5 - row;
- result = Rast_interp_lanczos(t, u, cell);
+ result = Rast_interp_lanczos(t, u, c);
Rast_set_f_value(obufptr, result, cell_type);
}
@@ -67,18 +67,18 @@
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) */
+ 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;
+ FCELL cell;
/* cut indices to integer */
- row = (int)floor(*row_idx);
- col = (int)floor(*col_idx);
+ 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) {
@@ -86,13 +86,12 @@
return;
}
- cellp = CPTR(ibuffer, row, col);
+ cell = CVAL(ibuffer, row, col);
/* if nearest is null, all the other interps will be null */
- if (Rast_is_f_null_value(cellp)) {
+ if (Rast_is_f_null_value(&cell)) {
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 */
Modified: grass/trunk/raster/r.proj/main.c
===================================================================
--- grass/trunk/raster/r.proj/main.c 2012-03-21 23:59:57 UTC (rev 51140)
+++ grass/trunk/raster/r.proj/main.c 2012-03-22 10:28:16 UTC (rev 51141)
@@ -97,15 +97,13 @@
overwrite, /* Overwrite */
curr_proj; /* output projection (see gis.h) */
- void *obuffer, /* buffer that holds one output row */
- *obufptr; /* column ptr in output buffer */
+ void *obuffer; /* buffer that holds one output row */
+
struct cache *ibuffer; /* buffer that holds the input map */
func interpolate; /* interpolation routine */
- double xcoord1, xcoord2, /* temporary x coordinates */
- ycoord1, ycoord2, /* temporary y coordinates */
- col_idx, /* column index in input matrix */
- row_idx, /* row index in input matrix */
+ double xcoord2, /* temporary x coordinates */
+ ycoord2, /* temporary y coordinates */
onorth, osouth, /* save original border coords */
oeast, owest, inorth, isouth, ieast, iwest;
char north_str[30], south_str[30], east_str[30], west_str[30];
@@ -487,44 +485,51 @@
cell_size = Rast_cell_size(cell_type);
- xcoord1 = xcoord2 = outcellhd.west + (outcellhd.ew_res / 2);
- /**/ ycoord1 = ycoord2 = outcellhd.north - (outcellhd.ns_res / 2);
- /**/ G_important_message(_("Projecting..."));
+ xcoord2 = outcellhd.west + (outcellhd.ew_res / 2);
+ ycoord2 = outcellhd.north - (outcellhd.ns_res / 2);
+ G_important_message(_("Projecting..."));
G_percent(0, outcellhd.rows, 2);
for (row = 0; row < outcellhd.rows; row++) {
- obufptr = obuffer;
+ /* obufptr = obuffer */;
+ #pragma omp parallel for schedule (static)
+
for (col = 0; col < outcellhd.cols; col++) {
+ void *obufptr = (void *)((const unsigned char *)obuffer + col * cell_size);
+
+ double xcoord1 = xcoord2 + (col) * outcellhd.ew_res;
+ double ycoord1 = ycoord2;
+
/* project coordinates in output matrix to */
/* coordinates in input matrix */
if (pj_do_proj(&xcoord1, &ycoord1, &oproj, &iproj) < 0)
Rast_set_null_value(obufptr, 1, cell_type);
else {
/* convert to row/column indices of input matrix */
- col_idx = (xcoord1 - incellhd.west) / incellhd.ew_res;
- row_idx = (incellhd.north - ycoord1) / incellhd.ns_res;
+ /* column index in input matrix */
+ double col_idx = (xcoord1 - incellhd.west) / incellhd.ew_res;
+ /* row index in input matrix */
+ double row_idx = (incellhd.north - ycoord1) / incellhd.ns_res;
+
/* and resample data point */
interpolate(ibuffer, obufptr, cell_type,
- &col_idx, &row_idx, &incellhd);
+ col_idx, row_idx, &incellhd);
}
- obufptr = G_incr_void_ptr(obufptr, cell_size);
- xcoord2 += outcellhd.ew_res;
- xcoord1 = xcoord2;
- ycoord1 = ycoord2;
+ /* obufptr = G_incr_void_ptr(obufptr, cell_size); */
}
Rast_put_row(fdo, obuffer, cell_type);
- xcoord1 = xcoord2 = outcellhd.west + (outcellhd.ew_res / 2);
+ xcoord2 = outcellhd.west + (outcellhd.ew_res / 2);
ycoord2 -= outcellhd.ns_res;
- ycoord1 = ycoord2;
G_percent(row, outcellhd.rows - 1, 2);
}
Rast_close(fdo);
+ release_cache(ibuffer);
if (have_colors > 0) {
Rast_write_colors(mapname, G_mapset(), &colr);
Modified: grass/trunk/raster/r.proj/nearest.c
===================================================================
--- grass/trunk/raster/r.proj/nearest.c 2012-03-21 23:59:57 UTC (rev 51140)
+++ grass/trunk/raster/r.proj/nearest.c 2012-03-22 10:28:16 UTC (rev 51141)
@@ -12,17 +12,17 @@
void p_nearest(struct cache *ibuffer, /* input buffer */
void *obufptr, /* ptr in output buffer */
int cell_type, /* raster map type of obufptr */
- double *col_idx, /* column index in input matrix */
- double *row_idx, /* row index in input matrix */
+ double col_idx, /* column index in input matrix */
+ double row_idx, /* row index in input matrix */
struct Cell_head *cellhd /* cell header of input layer */
)
{
int row, col; /* row/col of nearest neighbor */
- FCELL *cellp;
+ FCELL cell;
/* cut indices to integer */
- row = (int)floor(*row_idx);
- col = (int)floor(*col_idx);
+ 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) {
@@ -30,12 +30,12 @@
return;
}
- cellp = CPTR(ibuffer, row, col);
+ cell = CVAL(ibuffer, row, col);
- if (Rast_is_f_null_value(cellp)) {
+ if (Rast_is_f_null_value(&cell)) {
Rast_set_null_value(obufptr, 1, cell_type);
return;
}
- Rast_set_f_value(obufptr, *cellp, cell_type);
+ Rast_set_f_value(obufptr, cell, cell_type);
}
Modified: grass/trunk/raster/r.proj/r.proj.h
===================================================================
--- grass/trunk/raster/r.proj/r.proj.h 2012-03-21 23:59:57 UTC (rev 51140)
+++ grass/trunk/raster/r.proj/r.proj.h 2012-03-22 10:28:16 UTC (rev 51141)
@@ -17,6 +17,7 @@
struct cache
{
int fd;
+ char *fname;
int stride;
int nblocks;
block **grid;
@@ -24,7 +25,7 @@
int *refs;
};
-typedef void (*func) (struct cache *, void *, int, double *, double *,
+typedef void (*func) (struct cache *, void *, int, double, double,
struct Cell_head *);
struct menu
@@ -38,27 +39,28 @@
struct pj_info *);
extern struct cache *readcell(int, const char *);
extern block *get_block(struct cache *, int);
+extern void release_cache(struct cache *);
/* declare resampling methods */
/* bilinear.c */
-extern void p_bilinear(struct cache *, void *, int, double *, double *,
+extern void p_bilinear(struct cache *, void *, int, double, double,
struct Cell_head *);
/* cubic.c */
-extern void p_cubic(struct cache *, void *, int, double *, double *,
+extern void p_cubic(struct cache *, void *, int, double, double,
struct Cell_head *);
/* nearest.c */
-extern void p_nearest(struct cache *, void *, int, double *, double *,
+extern void p_nearest(struct cache *, void *, int, double, double,
struct Cell_head *);
/* bilinear_f.c */
-extern void p_bilinear_f(struct cache *, void *, int, double *, double *,
+extern void p_bilinear_f(struct cache *, void *, int, double, double,
struct Cell_head *);
/* cubic_f.c */
-extern void p_cubic_f(struct cache *, void *, int, double *, double *,
+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 *,
+extern void p_lanczos(struct cache *, void *, int, double, double,
struct Cell_head *);
-extern void p_lanczos_f(struct cache *, void *, int, double *, double *,
+extern void p_lanczos_f(struct cache *, void *, int, double, double,
struct Cell_head *);
#if 1
@@ -66,7 +68,7 @@
#define BKIDX(c,y,x) ((y) * (c)->stride + (x))
#define BKPTR(c,y,x) ((c)->grid[BKIDX((c),(y),(x))])
#define BLOCK(c,y,x) (BKPTR((c),(y),(x)) ? BKPTR((c),(y),(x)) : get_block((c),BKIDX((c),(y),(x))))
-#define CPTR(c,y,x) (&(*BLOCK((c),HI((y)),HI((x))))[LO((y))][LO((x))])
+#define CVAL(c,y,x) ((*BLOCK((c),HI((y)),HI((x))))[LO((y))][LO((x))])
#else
Modified: grass/trunk/raster/r.proj/readcell.c
===================================================================
--- grass/trunk/raster/r.proj/readcell.c 2012-03-21 23:59:57 UTC (rev 51140)
+++ grass/trunk/raster/r.proj/readcell.c 2012-03-22 10:28:16 UTC (rev 51141)
@@ -22,7 +22,6 @@
int nrows;
int ncols;
int row;
- char *filename;
int nx, ny;
int nblocks;
int i;
@@ -46,20 +45,21 @@
c->nblocks = nblocks;
c->grid = (block **) G_calloc(nx * ny, sizeof(block *));
c->blocks = (block *) G_malloc(nblocks * sizeof(block));
- c->refs = (int *)G_malloc(nblocks * sizeof(int));
+ c->refs = (int *)G_calloc(nblocks, sizeof(int));
if (nblocks < nx * ny) {
/* Temporary file must be created in output location */
G__switch_env();
- filename = G_tempfile();
+ c->fname = G_tempfile();
G__switch_env();
- c->fd = open(filename, O_RDWR | O_CREAT | O_EXCL, 0600);
+ c->fd = open(c->fname, O_RDWR | O_CREAT | O_EXCL, 0600);
if (c->fd < 0)
G_fatal_error(_("Unable to open temporary file"));
- remove(filename);
}
- else
+ else {
c->fd = -1;
+ c->fname = NULL;
+ }
G_important_message(_("Allocating memory and reading input map..."));
G_percent(0, nrows, 5);
@@ -97,36 +97,57 @@
G_free(tmpbuf);
- if (c->fd < 0)
+ if (c->fd < 0) {
for (i = 0; i < c->nblocks; i++) {
c->grid[i] = &c->blocks[i];
c->refs[i] = i;
}
+ }
+ else {
+ close(c->fd);
+ c->fd = -1;
+ }
return c;
}
block *get_block(struct cache * c, int idx)
{
+ int fd;
int replace = rand() % c->nblocks;
block *p = &c->blocks[replace];
int ref = c->refs[replace];
off_t offset = (off_t) idx * sizeof(FCELL) << L2BSIZE;
- if (c->fd < 0)
+ if (c->fname == NULL)
G_fatal_error(_("Internal error: cache miss on fully-cached map"));
+ fd = open(c->fname, O_RDONLY, 0600);
+ if (fd < 0)
+ G_fatal_error(_("Unable to open temporary file"));
+
if (ref >= 0)
c->grid[ref] = NULL;
c->grid[idx] = p;
c->refs[replace] = idx;
- if (lseek(c->fd, offset, SEEK_SET) < 0)
+ if (lseek(fd, offset, SEEK_SET) < 0)
G_fatal_error(_("Error seeking on segment file"));
- if (read(c->fd, p, sizeof(block)) < 0)
- G_fatal_error(_("Error writing segment file"));
+ if (read(fd, p, sizeof(block)) < 0)
+ G_fatal_error(_("Error reading segment file"));
+
+ close(fd);
return p;
}
+
+void release_cache(struct cache *c)
+{
+ G_free(c->grid);
+ G_free(c->blocks);
+ G_free(c->refs);
+ remove(c->fname);
+ G_free(c);
+}
More information about the grass-commit
mailing list