New g3d module for raster3d map retiling implemented with tests. Added more g3d std. parameter.

git-svn-id: https://svn.osgeo.org/grass/grass/trunk@46742 15284696-431f-4ddb-bdfa-cd5b030d7da7

include/gis.h

@@ -189,6 +189,10 @@ static const char *GRASS_copyright __attribute__ ((unused))
     G_OPT_R3_OUTPUT,		/* new output raster3d map */
     G_OPT_R3_MAP,		/* old input raster3d map */
     G_OPT_R3_MAPS,		/* old input raster3d maps */
+    G_OPT_R3_TYPE,              /* Type (FCELL or DCELL) of a new created raster3d map */
+    G_OPT_R3_PRECISION,         /* The precision of the new generated raster3d map */
+    G_OPT_R3_TILE_DIMENSION,    /* The tile dimension of a new generated raster3d map */
+    G_OPT_R3_COMPRESSION,       /* The kind of compression of a new created raster3d map */
 
     G_OPT_V_INPUT,		/* old input vector map */
     G_OPT_V_INPUTS,		/* old input vector maps */

lib/g3d/g3dparam.c

@@ -42,43 +42,13 @@ void G3d_setStandard3dInputParams()
 {
     param = G3d_malloc(sizeof(G3d_paramType));
 
-    param->type = G_define_option();
-    param->type->key = "type";
-    param->type->type = TYPE_STRING;
-    param->type->required = NO;
-    param->type->multiple = NO;
-    param->type->answer = "default";
-    param->type->options = "default,double,float";
-    param->type->description = _("Data type used in the output file");
-
-    param->precision = G_define_option();
-    param->precision->key = "precision";
-    param->precision->type = TYPE_STRING;
-    param->precision->required = NO;
-    param->precision->multiple = NO;
-    param->precision->answer = "default";
-    param->precision->description =
-	_("Precision used in the output file (default, max, or 0 to 52)");
-
-    param->compression = G_define_option();
-    param->compression->key = "compression";
-    param->compression->type = TYPE_STRING;
-    param->compression->required = NO;
-    param->compression->multiple = NO;
-    param->compression->answer = "default";
-    param->compression->options = "default,rle,none";
-    param->compression->description =
-	_("The compression method used in the output file");
-
-    param->dimension = G_define_option();
-    param->dimension->key = "tiledimension";
-    param->dimension->type = TYPE_STRING;
-    param->dimension->required = NO;
-    param->dimension->multiple = NO;
-    param->dimension->key_desc = "XxYxZ";
-    param->dimension->answer = "default";
-    param->dimension->description =
-	_("The dimensions of the tiles used in the output file");
+    param->type = G_define_standard_option(G_OPT_R3_TYPE);
+
+    param->precision = G_define_standard_option(G_OPT_R3_PRECISION);
+
+    param->compression = G_define_standard_option(G_OPT_R3_COMPRESSION);
+
+    param->dimension = G_define_standard_option(G_OPT_R3_TILE_DIMENSION);
 }
 
 /*----------------------------------------------------------------------------*/

lib/g3d/retile.c

@@ -37,10 +37,10 @@ retileNocache(void *map, const char *nameOut, int tileX, int tileY, int tileZ)
 
     G3d_getNofTilesMap(map2, &nx, &ny, &nz);
 
-    for (z = 0; z < nz; z++)
+    for (z = 0; z < nz; z++) {
+        G_percent(z, nz, 1);
 	for (y = 0; y < ny; y++)
 	    for (x = 0; x < nx; x++) {
-
 		G3d_getBlock(map, x * tileX, y * tileY, z * tileZ,
 			     tileX, tileY, tileZ, data, typeIntern);
 		if (!G3d_writeTile
@@ -49,7 +49,10 @@ retileNocache(void *map, const char *nameOut, int tileX, int tileY, int tileZ)
 		    G3d_fatalError
 			("G3d_retileNocache: error in G3d_writeTile");
 	    }
-
+    }
+    
+    G_percent(1, 1, 1);
+        
     G3d_freeTiles(data);
     G3d_closeCell(map2);
 }
@@ -116,6 +119,7 @@ G3d_retile(void *map, const char *nameOut, int tileX, int tileY, int tileZ)
     G3d_getCoordsMap(map, &rows, &cols, &depths);
 
     for (z = 0; z < depths; z++) {
+        G_percent(z, depths, 1);
 	G3d_coord2tileCoord(map2, x, y, z, &xTile, &yTile, &zTile,
 			    &xOffs, &yOffs, &zOffs);
 	if (zTile > prev) {
@@ -133,6 +137,7 @@ G3d_retile(void *map, const char *nameOut, int tileX, int tileY, int tileZ)
 	    }
     }
 
+    G_percent(1, 1, 1);
     if (!G3d_flushAllTiles(map2))
 	G3d_fatalError("G3d_retile: error in G3d_flushAllTiles");
     if (!G3d_closeCell(map2))

lib/gis/parser_standard_options.c

@@ -302,7 +302,45 @@ struct Option *G_define_standard_option(int opt)
 	Opt->gisprompt = "old,grid3,3d-raster";
 	Opt->description = _("Name of 3D raster map(s)");
 	break;
-
+    case G_OPT_R3_TYPE:
+        Opt->key = "type";
+        Opt->type = TYPE_STRING;
+        Opt->required = NO;
+        Opt->multiple = NO;
+        Opt->answer = "default";
+        Opt->options = "default,double,float";
+        Opt->description = _("Data type used in the output file");
+	break;
+    case G_OPT_R3_PRECISION:
+        Opt->key = "precision";
+        Opt->type = TYPE_STRING;
+        Opt->required = NO;
+        Opt->multiple = NO;
+        Opt->answer = "default";
+        Opt->description =
+            _("Precision used in the output file (default, max, or 0 to 52)");
+	break;
+    case G_OPT_R3_COMPRESSION:
+        Opt->key = "compression";
+        Opt->type = TYPE_STRING;
+        Opt->required = NO;
+        Opt->multiple = NO;
+        Opt->answer = "default";
+        Opt->options = "default,rle,none";
+        Opt->description =
+            _("The compression method used in the output file");
+	break;
+    case G_OPT_R3_TILE_DIMENSION:
+        Opt->key = "tiledimension";
+        Opt->type = TYPE_STRING;
+        Opt->required = NO;
+        Opt->multiple = NO;
+        Opt->key_desc = "XxYxZ";
+        Opt->answer = "default";
+        Opt->description =
+            _("The dimensions of the tiles used in the output file (XxYxZ)");
+	break;
+        
 	/*vector maps */
     case G_OPT_V_INPUT:
 	Opt->key = "input";

raster3d/Makefile

@@ -10,6 +10,7 @@ SUBDIRS = \
 	r3.out.ascii \
 	r3.out.v5d \
 	r3.out.vtk \
+	r3.retile \
 	r3.stats \
 	r3.timestamp \
 	r3.to.rast \

raster3d/r3.info/r3.info.main.c

@@ -66,6 +66,7 @@ int main(int argc, char *argv[])
     struct Flag *tflag;
     struct Flag *timestampflag;
     struct Flag *gflag;
+    struct Flag *iflag;
     struct Flag *hflag;
     int data_type;
 
@@ -104,6 +105,10 @@ int main(int argc, char *argv[])
     gflag->key = 'g';
     gflag->description = _("Print 3D raster map region only");
 
+    iflag = G_define_flag();
+    iflag->key = 'i';
+    iflag->description = _("Print 3D raster tile number and size information nly");
+    
     hflag = G_define_flag();
     hflag->key = 'h';
     hflag->description = _("Print 3D raster history instead of info");
@@ -148,7 +153,7 @@ int main(int argc, char *argv[])
 
     /*Show the info if no flag is set */
     if (!rflag->answer && !sflag->answer && !tflag->answer && !gflag->answer
-	&& !timestampflag->answer && !hflag->answer) {
+	&& !timestampflag->answer && !hflag->answer && !iflag->answer) {
 	divider('+');
 
 	if (G_asprintf(&line, "Layer:    %-29.29s  Date: %s", name,
@@ -254,7 +259,7 @@ int main(int argc, char *argv[])
 	    else
 		G_fatal_error(_("Cannot allocate memory for string"));
 
-	    if (G_asprintf(&line, "  Number of tiles:      %ld",
+	    if (G_asprintf(&line, "  Number of tiles:      %d",
                 ((G3D_Map* )g3map)->nTiles) > 0)
 		printline(line);
 	    else
@@ -437,6 +442,14 @@ int main(int argc, char *argv[])
 		    data_type == DCELL_TYPE ? "DCELL" :
 		    "??");
 
+	}			/*Resolution */
+	else if (iflag->answer) {
+	    fprintf(out, "tilenumx=%d\n", ((G3D_Map* )g3map)->nx);
+	    fprintf(out, "tilenumy=%d\n", ((G3D_Map* )g3map)->ny);
+	    fprintf(out, "tilenumz=%d\n", ((G3D_Map* )g3map)->nz);
+	    fprintf(out, "tiledimx=%d\n", ((G3D_Map* )g3map)->tileX);
+	    fprintf(out, "tiledimy=%d\n", ((G3D_Map* )g3map)->tileY);
+	    fprintf(out, "tiledimz=%d\n", ((G3D_Map* )g3map)->tileZ);
 	}			/*History output */
 	else if (hflag->answer) {
 	    if (hist_ok) {

raster3d/r3.retile/Makefile

@@ -0,0 +1,10 @@
+MODULE_TOPDIR = ../..
+
+PGM=r3.retile
+
+LIBES = $(G3DLIB) $(GISLIB)
+DEPENDENCIES = $(G3DDEP) $(GISDEP)
+
+include $(MODULE_TOPDIR)/include/Make/Module.make
+
+default: cmd

raster3d/r3.retile/main.c

@@ -0,0 +1,158 @@
+
+/****************************************************************************
+ *
+ * MODULE:       r3.retile
+ *   	    	
+ * AUTHOR(S):    Original author 
+ *               Soeren Gebbert soerengebbert <at> googlemail <dot> co
+ * 
+ * PURPOSE:      Retiles an existing G3D map with user defined x, y and z tile size
+ *
+ * COPYRIGHT:    (C) 2011 by the GRASS Development Team
+ *
+ *               This program is free software under the GNU General Public
+ *   	    	License (>=v2). Read the file COPYING that comes with GRASS
+ *   	    	for details.
+ *
+ *****************************************************************************/
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <grass/gis.h>
+#include <grass/raster.h>
+#include <grass/G3d.h>
+#include <grass/glocale.h>
+
+/*- Parameters and global variables -----------------------------------------*/
+typedef struct {
+    struct Option *input, *output, *tiling;
+    struct Flag *cache;
+} paramType;
+
+paramType param; /*Parameters */
+
+/*- prototypes --------------------------------------------------------------*/
+void fatal_error(void *map, int *fd, int depths, char *errorMsg); /*Simple Error message */
+void set_params(); /*Fill the paramType structure */
+void g3d_to_raster(void *map, G3D_Region region, int *fd); /*Write the raster */
+int open_output_map(const char *name, int res_type); /*opens the outputmap */
+void close_output_map(int fd); /*close the map */
+
+
+
+/* ************************************************************************* */
+/* Error handling ********************************************************** */
+
+/* ************************************************************************* */
+void fatal_error(void *map, int *fd, int depths, char *errorMsg)
+{
+    int i;
+
+    /* Close files and exit */
+    if (map != NULL) {
+        if (!G3d_closeCell(map))
+            G3d_fatalError(_("Unable to close the 3d raster map"));
+    }
+
+    if (fd != NULL) {
+        for (i = 0; i < depths; i++)
+            Rast_unopen(fd[i]);
+    }
+
+    G3d_fatalError(errorMsg);
+    exit(EXIT_FAILURE);
+
+}
+
+/* ************************************************************************* */
+/* Set up the arguments we are expecting ********************************** */
+
+/* ************************************************************************* */
+void set_params()
+{
+    param.input = G_define_standard_option(G_OPT_R3_INPUT);
+
+    param.output = G_define_standard_option(G_OPT_R3_OUTPUT);
+    param.output->description = _("Name of the retiled 3D raster map");
+    
+    param.tiling = G_define_standard_option(G_OPT_R3_TILE_DIMENSION);
+    
+    param.cache = G_define_flag();
+    param.cache->key = 'c';
+    param.cache->description = "Disable tile caching";
+}
+
+
+/* ************************************************************************* */
+/* Main function, open the G3D map and create the raster maps ************** */
+
+/* ************************************************************************* */
+int main(int argc, char *argv[])
+{
+    struct GModule *module;
+    G3D_Map *map = NULL;
+    int tileX, tileY, tileZ;
+    char *mapset;
+
+    /* Initialize GRASS */
+    G_gisinit(argv[0]);
+
+    module = G_define_module();
+    G_add_keyword(_("raster3d"));
+    G_add_keyword(_("retile"));
+    G_add_keyword(_("voxel"));
+    module->description = _("Retiles an existing G3D map with user defined x, y and z tile size");
+
+    /* Get parameters from user */
+    set_params();
+
+    /* Have GRASS get inputs */
+    if (G_parser(argc, argv))
+        exit(EXIT_FAILURE);
+
+    G_debug(3, _("Open 3d raster map <%s>"), param.input->answer);
+
+    mapset = G_find_grid3(param.input->answer, "");
+    
+    if (mapset == NULL)
+        G3d_fatalError(_("3d raster map <%s> not found"),
+                       param.input->answer);
+
+    /*Set the defaults */
+    G3d_initDefaults();
+    
+    if(!param.cache->answer)
+        map = G3d_openCellOld(param.input->answer, mapset, G3D_DEFAULT_WINDOW, 
+                          G3D_TILE_SAME_AS_FILE, G3D_USE_CACHE_DEFAULT);
+    else
+        map = G3d_openCellOld(param.input->answer, mapset, G3D_DEFAULT_WINDOW, 
+                          G3D_TILE_SAME_AS_FILE, G3D_NO_CACHE);
+
+    if (map == NULL)
+        G3d_fatalError(_("Unable to open 3d raster map <%s>"),
+                       param.input->answer);
+
+    /* Get the tile dimension */
+    G3d_getTileDimension(&tileX, &tileY, &tileZ);    
+    if (strcmp(param.tiling->answer, "default") != 0) {
+	if (sscanf(param.tiling->answer, "%dx%dx%d",
+		   &tileX, &tileY, &tileZ) != 3) {
+	    G3d_fatalError(_("G3d_getStandard3dParams: tile dimension value invalid"));
+	}
+    }
+    
+    if(!param.cache->answer)
+        G_message("Retile map with tile cache enabled");
+    else
+        G_message("Retile map without tile caching");
+    
+    G3d_retile(map, param.output->answer, tileX, tileY, tileZ);
+
+    /* Close files and exit */
+    if (!G3d_closeCell(map))
+        fatal_error(map, NULL, 0, _("Error closing 3d raster map"));
+
+    map = NULL;
+
+    return (EXIT_SUCCESS);
+}

raster3d/r3.retile/r3.retile.html

@@ -0,0 +1,16 @@
+<H2>DESCRIPTION</H2>
+
+Writes a copy of an existing G3D map with a user defined 
+number of tiles in x, y and z direction. 
+The precision and the type of the original G3D map are used for the new retiled G3D map. 
+
+<H2>SEE ALSO</H2>
+
+<EM><A HREF="r3.info.html">r3.cross.rast</A></EM><br>
+<EM><A HREF="r3.out.ascii.html">r3.out.ascii</A></EM><br>
+<EM><A HREF="g.region.html">g.region</A></EM><br>
+
+<H2>AUTHOR</H2>
+Soeren Gebbert
+
+<p><i>Last changed: $Date: 2010-09-24 14:05:56 +0200 (Fr, 24. Sep 2010) $</i>

raster3d/r3.retile/test.r3.retile.sh

@@ -0,0 +1,41 @@
+# The module r3.retile is the object of ths test case
+
+
+# We need to set a specific region in the
+# @preprocess step of this test. We generate
+# voxel data with r3.mapcalc. The region setting 
+# should work for UTM and LL test locations
+g.region s=0 n=47 w=0 e=35 b=0 t=23 res=10 res3=1 -p3
+# Now create a voxel map with value = col + row + depth. 
+r3.mapcalc --o expr="volume = col() + row() + depth()"
+# Add null value information
+r3.mapcalc --o expr="volume_null = if(row() == 1 || row() == 5, null(), volume)"
+
+# We @test the retiling with different tile sizes, with and without tile caching.
+# The validation is based on @raster3d maps and r3.info generated @files. 
+# We use @precision=0 for data export
+r3.retile input=volume_null output=test_retile_map_1  tiledim=1x1x1
+r3.info -i map=test_retile_map_1 > test_retile_map_1_info.txt
+r3.retile input=volume_null output=test_retile_map_2  tiledim=1x1x1 -c
+r3.info -i map=test_retile_map_2 > test_retile_map_2_info.txt
+r3.retile input=volume_null output=test_retile_map_3  tiledim=8x8x4 
+r3.info -i map=test_retile_map_3 > test_retile_map_3_info.txt
+r3.retile input=volume_null output=test_retile_map_4  tiledim=8x8x4 -c
+r3.info -i map=test_retile_map_4 > test_retile_map_4_info.txt
+r3.retile input=volume_null output=test_retile_map_5  tiledim=16x16x16 
+r3.info -i map=test_retile_map_5 > test_retile_map_5_info.txt
+r3.retile input=volume_null output=test_retile_map_6  tiledim=16x16x16 -c
+r3.info -i map=test_retile_map_6 > test_retile_map_6_info.txt
+r3.retile input=volume_null output=test_retile_map_7  tiledim=35x47x23 
+r3.info -i map=test_retile_map_7 > test_retile_map_7_info.txt
+r3.retile input=volume_null output=test_retile_map_8  tiledim=35x47x23 -c
+r3.info -i map=test_retile_map_8 > test_retile_map_8_info.txt
+r3.retile input=volume_null output=test_retile_map_9  tiledim=34x46x22 
+r3.info -i map=test_retile_map_9 > test_retile_map_9_info.txt
+r3.retile input=volume_null output=test_retile_map_10 tiledim=34x46x22 -c
+r3.info -i map=test_retile_map_10 > test_retile_map_10_info.txt
+
+# Create the validation files
+# for map in `g.mlist type=rast3d pattern=test_retile_map_*` ; do
+#  r3.out.ascii input=${map} output=${map}.ref dp=0
+# done

raster3d/r3.retile/test_retile_map_10_info.ref

@@ -0,0 +1,6 @@
+tilenumx=2
+tilenumy=2
+tilenumz=2
+tiledimx=34
+tiledimy=46
+tiledimz=22

raster3d/r3.retile/test_retile_map_1_info.ref

@@ -0,0 +1,6 @@
+tilenumx=35
+tilenumy=47
+tilenumz=23
+tiledimx=1
+tiledimy=1
+tiledimz=1

raster3d/r3.retile/test_retile_map_2_info.ref

@@ -0,0 +1,6 @@
+tilenumx=35
+tilenumy=47
+tilenumz=23
+tiledimx=1
+tiledimy=1
+tiledimz=1

raster3d/r3.retile/test_retile_map_3_info.ref

@@ -0,0 +1,6 @@
+tilenumx=5
+tilenumy=6
+tilenumz=6
+tiledimx=8
+tiledimy=8
+tiledimz=4

raster3d/r3.retile/test_retile_map_4_info.ref

@@ -0,0 +1,6 @@
+tilenumx=5
+tilenumy=6
+tilenumz=6
+tiledimx=8
+tiledimy=8
+tiledimz=4

raster3d/r3.retile/test_retile_map_5_info.ref

@@ -0,0 +1,6 @@
+tilenumx=3
+tilenumy=3
+tilenumz=2
+tiledimx=16
+tiledimy=16
+tiledimz=16

raster3d/r3.retile/test_retile_map_6_info.ref

@@ -0,0 +1,6 @@
+tilenumx=3
+tilenumy=3
+tilenumz=2
+tiledimx=16
+tiledimy=16
+tiledimz=16

raster3d/r3.retile/test_retile_map_7_info.ref

@@ -0,0 +1,6 @@
+tilenumx=1
+tilenumy=1
+tilenumz=1
+tiledimx=35
+tiledimy=47
+tiledimz=23

raster3d/r3.retile/test_retile_map_8_info.ref

@@ -0,0 +1,6 @@
+tilenumx=1
+tilenumy=1
+tilenumz=1
+tiledimx=35
+tiledimy=47
+tiledimz=23

raster3d/r3.retile/test_retile_map_9_info.ref

@@ -0,0 +1,6 @@
+tilenumx=2
+tilenumy=2
+tilenumz=2
+tiledimx=34
+tiledimy=46
+tiledimz=22