/**************************************************************************** * * MODULE: r.out.vtk * * AUTHOR(S): Original author * Soeren Gebbert soerengebbert@gmx.de * 08 23 2005 Berlin * PURPOSE: Converts raster maps into the VTK-Ascii format * * COPYRIGHT: (C) 2005 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 #include #include #include #include #include "writeascii.h" #define MAIN #include "parameters.h" #include "globaldefs.h" /*local protos */ void CheckInputMaps(); /* ************************************************************************* */ /* Check the input maps **************************************************** */ /* ************************************************************************* */ void CheckInputMaps() { char *mapset = NULL; int i; /*Check for elevation map */ if (param.elevationmap->answer) { mapset = G_find_cell2(param.elevationmap->answer, ""); if (mapset == NULL) { G_fatal_error(_("Cell file [%s] not found\n"), param.elevationmap->answer); exit(EXIT_FAILURE); } } /*Check for normal input map */ if (param.input->answers != NULL) { for (i = 0; param.input->answers[i] != NULL; i++) { mapset = NULL; mapset = G_find_cell2(param.input->answers[i], ""); if (mapset == NULL) { G_fatal_error(_("Cell file [%s] not found\n"), param.input->answers[i]); exit(EXIT_FAILURE); } } } /*RGB raster maps */ if (param.rgbmaps->answers != NULL) { if (param.rgbmaps->answers[0] != NULL && param.rgbmaps->answers[1] != NULL && param.rgbmaps->answers[2] != NULL) { /*Loop over all input maps! */ for (i = 0; i < 3; i++) { mapset = NULL; mapset = G_find_cell2(param.rgbmaps->answers[i], ""); if (mapset == NULL) { G_fatal_error(_("rgb cell map [%s] not found\n"), param.rgbmaps->answers[i]); exit(EXIT_FAILURE); } } } else { G_fatal_error(_ ("Cannot create RGB data, please provide three maps [r,g,b]")); } } /*Vector raster maps */ if (param.vectmaps->answers != NULL) { if (param.vectmaps->answers[0] != NULL && param.vectmaps->answers[1] != NULL && param.vectmaps->answers[2] != NULL) { /*Loop over all input maps! */ for (i = 0; i < 3; i++) { mapset = NULL; mapset = G_find_cell2(param.vectmaps->answers[i], ""); if (mapset == NULL) { G_fatal_error(_("vector cell map [%s] not found\n"), param.vectmaps->answers[i]); exit(EXIT_FAILURE); } } } else { G_fatal_error(_ ("Cannot create vector data, please provide three maps [x,y,z]")); } } /*Give a warning if no output cell/point or rgb data is specified */ if (param.input->answers == NULL && param.rgbmaps->answers == NULL && param.vectmaps->answers == NULL) { G_warning ("No g3d data, RGB or vector maps are provided! Will only write the geometry."); } return; } /* ************************************************************************* */ /* MAIN ******************************************************************** */ /* ************************************************************************* */ int main(int argc, char *argv[]) { struct Cell_head region; struct Cell_head default_region; FILE *fp = NULL; struct GModule *module; int i = 0, polytype = 0; char *null_value, *mapset; int out_type; int fd; /*Normale maps ;) */ int rgbfd[3]; int vectfd[3]; int celltype[3] = { 0, 0, 0 }; int headertype; /* Initialize GRASS */ G_gisinit(argv[0]); module = G_define_module(); module->keywords = _("raster"); module->description = _("Converts raster maps into the VTK-Ascii format"); /* Get parameters from user */ SetParameters(); /* Have GRASS get inputs */ if (G_parser(argc, argv)) exit(EXIT_FAILURE); /*open the output */ if (param.output->answer) { fp = fopen(param.output->answer, "w"); if (fp == NULL) { perror(param.output->answer); G_usage(); exit(EXIT_FAILURE); } } else fp = stdout; /*Check the input maps */ CheckInputMaps(); /*Correct the coordinates, so the precision of VTK is not hurt :( */ if(param.coorcorr->answer){ /*Get the default region for coordiante correction*/ G_get_default_window(&default_region); /*Use the center of the current region as extent*/ y_extent = (default_region.north + default_region.south)/2; x_extent = (default_region.west + default_region.east)/2; } else { x_extent = 0; y_extent = 0; } /* Figure out the region from the map */ G_get_window(®ion); /*Set the null Value, maybe i have to check this? */ null_value = param.null_val->answer; /********************* WRITE ELEVATION *************************************/ if (param.elevationmap->answer) { /*If the elevation is set, write the correct Header */ if (param.usestruct->answer) { writeVTKStructuredElevationHeader(fp, region); } else { writeVTKPolygonalElevationHeader(fp, region); } G_debug(3, _("Open Raster file %s"), param.elevationmap->answer); mapset = G_find_cell2(param.elevationmap->answer, ""); out_type = G_raster_map_type(param.elevationmap->answer, mapset); /* open raster file */ fd = G_open_cell_old(param.elevationmap->answer, mapset); if (fd < 0) { G_fatal_error(_("Could not open map %s\n"), param.elevationmap->answer); exit(EXIT_FAILURE); } /*The write the Coordinates */ if (param.usestruct->answer) { writeVTKStructuredCoordinates(fd, fp, param.elevationmap->answer, region, out_type, null_value, atof(param.elevscale->answer)); } else { polytype = QUADS; /*The default */ if (param.usetriangle->answer) polytype = TRIANGLE_STRIPS; if (param.usevertices->answer) polytype = VERTICES; writeVTKPolygonalCoordinates(fd, fp, param.elevationmap->answer, region, out_type, null_value, atof(param.elevscale->answer), polytype); } G_close_cell(fd); } else { /*Should pointdata or celldata be written */ if (param.point->answer) headertype = 1; else headertype = 0; /*If no elevation is given, write the normal Header */ if (param.origin->answer) writeVTKNormalHeader(fp, region, atof(param.elevscale->answer) * (atof(param.elev->answer)), headertype); else writeVTKNormalHeader(fp, region, atof(param.elev->answer), headertype); } /******************** WRITE THE POINT OR CELL DATA HEADER ******************/ if (param.input->answers != NULL || param.rgbmaps->answers != NULL) { if (param.point->answer || param.elevationmap->answer) writeVTKPointDataHeader(fp, region); else writeVTKCellDataHeader(fp, region); } /********************** WRITE NORMAL DATA; CELL OR POINT *******************/ /*Loop over all input maps! */ if (param.input->answers != NULL) { for (i = 0; param.input->answers[i] != NULL; i++) { G_debug(3, _("Open Raster file %s"), param.input->answers[i]); mapset = NULL; mapset = G_find_cell2(param.input->answers[i], ""); out_type = G_raster_map_type(param.input->answers[i], mapset); /* open raster file */ fd = G_open_cell_old(param.input->answers[i], mapset); if (fd < 0) { G_fatal_error(_("Could not open map %s\n"), param.input->answers[i]); exit(EXIT_FAILURE); } /*Now write the data */ writeVTKData(fd, fp, param.input->answers[i], region, out_type, null_value); G_close_cell(fd); } } /********************** WRITE RGB IMAGE DATA; CELL OR POINT ****************/ if (param.rgbmaps->answers != NULL) { if (param.rgbmaps->answers[0] != NULL && param.rgbmaps->answers[1] != NULL && param.rgbmaps->answers[2] != NULL) { /*Loop over all three rgb input maps! */ for (i = 0; i < 3; i++) { G_debug(3, _("Open Raster file %s"), param.rgbmaps->answers[i]); mapset = NULL; mapset = G_find_cell2(param.rgbmaps->answers[i], ""); celltype[i] = G_raster_map_type(param.rgbmaps->answers[i], mapset); /* open raster file */ rgbfd[i] = G_open_cell_old(param.rgbmaps->answers[i], mapset); if (rgbfd[i] < 0) { G_fatal_error(_("Could not open map %s\n"), param.rgbmaps->answers[i]); exit(EXIT_FAILURE); } } /*Maps have to be from the same type */ if (celltype[0] == celltype[1] && celltype[0] == celltype[2]) { G_debug(3, _("Writing VTK ImageData\n")); out_type = celltype[0]; /*Now write the data */ writeVTKRGBImageData(rgbfd[0], rgbfd[1], rgbfd[2], fp, "RGB_Image", region, out_type); } else { G_warning(_ ("Wrong RGB maps. Maps should have the same type! RGB output not added!")); /*do nothing */ } /*Close the maps */ for (i = 0; i < 3; i++) G_close_cell(rgbfd[i]); } } /********************** WRITE VECTOR DATA; CELL OR POINT ****************/ if (param.vectmaps->answers != NULL) { if (param.vectmaps->answers[0] != NULL && param.vectmaps->answers[1] != NULL && param.vectmaps->answers[2] != NULL) { /*Loop over all three vect input maps! */ for (i = 0; i < 3; i++) { G_debug(3, _("Open Raster file %s"), param.vectmaps->answers[i]); mapset = NULL; mapset = G_find_cell2(param.vectmaps->answers[i], ""); celltype[i] = G_raster_map_type(param.vectmaps->answers[i], mapset); /* open raster file */ vectfd[i] = G_open_cell_old(param.vectmaps->answers[i], mapset); if (vectfd[i] < 0) { G_fatal_error(_("Could not open map %s\n"), param.vectmaps->answers[i]); exit(EXIT_FAILURE); } } /*Maps have to be from the same type */ if (celltype[0] == celltype[1] && celltype[0] == celltype[2]) { G_debug(3, _("Writing VTK Vector Data\n")); out_type = celltype[0]; /*Now write the data */ writeVTKVectorData(vectfd[0], vectfd[1], vectfd[2], fp, "Vector_Data", region, out_type); } else { G_warning(_ ("Wrong vector maps. Maps should have the same type! Vector output not added!")); /*do nothing */ } /*Close the maps */ for (i = 0; i < 3; i++) G_close_cell(vectfd[i]); } } if (param.output->answer && fp != NULL) if (fclose(fp)) { G_fatal_error(_("Error closing VTK-ASCII file")); exit(EXIT_FAILURE); } return 0; }