Skip to content

HTTPS clone URL

Subversion checkout URL

You can clone with
or
.
Download ZIP
Browse files

initial import of points2grid

  • Loading branch information...
commit 6364eb918664a7add8a341fb9632aef92d5dc4b2 0 parents
@hobu hobu authored
12 CHANGELOG
@@ -0,0 +1,12 @@
+Version 1.0.1:
+* Updated LICENSE to a more general 4-clause BSD license
+
+Version 1.0:
+
+* Input point clouds in ASCII and LAS formats
+* Outputs in Arc ASCII and GRD ASCII formats
+* DEM products including min, max, mean, and idw values
+* Point count values for every grid cell
+* Configurable search radius and grid resolution
+* Null filling, using a configurable moving window, to fill cells that have
+ null values
77 CoreInterp.h
@@ -0,0 +1,77 @@
+/*
+*
+COPYRIGHT AND LICENSE
+
+Copyright (c) 2011 The Regents of the University of California.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or
+without modification, are permitted provided that the following
+conditions are met:
+
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above
+copyright notice, this list of conditions and the following
+disclaimer in the documentation and/or other materials provided
+with the distribution.
+
+3. All advertising materials mentioning features or use of this
+software must display the following acknowledgement: This product
+includes software developed by the San Diego Supercomputer Center.
+
+4. Neither the names of the Centers nor the names of the contributors
+may be used to endorse or promote products derived from this
+software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS''
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS
+OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
+USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
+AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
+ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+POSSIBILITY OF SUCH DAMAGE.
+*
+*
+* Based on the notes by Prof. Ramon Arrowsmith(ramon.arrowsmith@asu.edu)
+* Authors: Han S Kim (hskim@cs.ucsd.edu), Sriram Krishnan (sriram@sdsc.edu)
+*
+*/
+
+#ifndef _CORE_INTERP_H_
+#define _CORE_INTERP_H_
+
+class CoreInterp
+{
+ public:
+ CoreInterp() {};
+ virtual ~CoreInterp() {};
+
+ virtual int init() = 0;
+ virtual int update(double data_x, double data_y, double data_z) = 0;
+ virtual int finish(char *outputName, int outputFormat, unsigned int outputType) = 0;
+
+ protected:
+ double GRID_DIST_X;
+ double GRID_DIST_Y;
+
+ int GRID_SIZE_X; // total size of a grid
+ int GRID_SIZE_Y; //
+
+ // for outputting
+ double min_x;
+ double max_x;
+ double min_y;
+ double max_y;
+
+ // for DEM filling
+ int window_size;
+};
+
+#endif
64 Global.h
@@ -0,0 +1,64 @@
+/*
+*
+COPYRIGHT AND LICENSE
+
+Copyright (c) 2011 The Regents of the University of California.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or
+without modification, are permitted provided that the following
+conditions are met:
+
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above
+copyright notice, this list of conditions and the following
+disclaimer in the documentation and/or other materials provided
+with the distribution.
+
+3. All advertising materials mentioning features or use of this
+software must display the following acknowledgement: This product
+includes software developed by the San Diego Supercomputer Center.
+
+4. Neither the names of the Centers nor the names of the contributors
+may be used to endorse or promote products derived from this
+software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS''
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS
+OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
+USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
+AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
+ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+POSSIBILITY OF SUCH DAMAGE.
+*
+*/
+
+#ifndef _GLOBAL_H_
+#define _GLOBAL_H_
+
+static const unsigned int OUTPUT_TYPE_MIN = 0x00000001;
+static const unsigned int OUTPUT_TYPE_MAX = 0x00000010;
+static const unsigned int OUTPUT_TYPE_MEAN = 0x00000100;
+static const unsigned int OUTPUT_TYPE_IDW = 0x00001000;
+static const unsigned int OUTPUT_TYPE_DEN = 0x00010000;
+static const unsigned int OUTPUT_TYPE_ALL = 0x00011111;
+
+enum OUTPUT_FORMAT {
+ OUTPUT_FORMAT_ALL = 0,
+ OUTPUT_FORMAT_ARC_ASCII,
+ OUTPUT_FORMAT_GRID_ASCII,
+};
+
+enum INPUT_FORMAT {
+ INPUT_ASCII = 0,
+ INPUT_LAS = 1,
+};
+
+#endif
177 GridFile.cpp
@@ -0,0 +1,177 @@
+/*
+*
+COPYRIGHT AND LICENSE
+
+Copyright (c) 2011 The Regents of the University of California.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or
+without modification, are permitted provided that the following
+conditions are met:
+
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above
+copyright notice, this list of conditions and the following
+disclaimer in the documentation and/or other materials provided
+with the distribution.
+
+3. All advertising materials mentioning features or use of this
+software must display the following acknowledgement: This product
+includes software developed by the San Diego Supercomputer Center.
+
+4. Neither the names of the Centers nor the names of the contributors
+may be used to endorse or promote products derived from this
+software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS''
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS
+OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
+USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
+AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
+ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+POSSIBILITY OF SUCH DAMAGE.
+*
+*
+* Based on the notes by Prof. Ramon Arrowsmith(ramon.arrowsmith@asu.edu)
+* Authors: Han S Kim (hskim@cs.ucsd.edu), Sriram Krishnan (sriram@sdsc.edu)
+*
+*/
+
+
+#include "GridFile.h"
+#include <float.h>
+#include <sys/mman.h>
+#include <fcntl.h>
+#include <unistd.h>
+#include <errno.h>
+#include <sys/types.h>
+#include <sys/stat.h>
+#include <string.h>
+
+GridFile::GridFile(int id, char *_fname, int _size_x, int _size_y)
+{
+ FILE *fp;
+ size_t n = 0;
+
+ ID = id;
+ strncpy(fname, _fname, sizeof(fname));
+ size_x = _size_x;
+ size_y = _size_y;
+ inMemory = false;
+
+ // for each file, you have to initialize every points
+ // then map.initialized = true;
+ if((fp = fopen64(fname, "w+")) == NULL)
+ {
+ fprintf(stderr, "%s fopen error %d(%s) \n", fname, errno, strerror(errno));
+ }
+
+ GridPoint init_value = {DBL_MAX, -DBL_MAX, 0, 0, 0, 0, 0, 0};
+
+ for(int i = 0; i < size_x * size_y; i++)
+ n += fwrite(&init_value, sizeof(GridPoint), 1, fp);
+
+ cout << id << ". file size: " << n * sizeof(GridPoint) << endl;
+ fclose(fp);
+
+ /*
+ if((filedes = open(fname, O_RDWR)) < 0)
+ {
+ fprintf(stderr, "%s open error %d(%s)\n", fname, errno, strerror(errno));
+ }
+ */
+
+ //cout << "file initialized: " << ID << endl;
+}
+
+GridFile::~GridFile()
+{
+ if(inMemory == true){
+ munmap(interp, sizeof(GridPoint) * size_x * size_y);
+ inMemory = false;
+ interp = NULL;
+ }
+
+ if(filedes >= 0)
+ close(filedes);
+
+ unlink(fname);
+
+ //cout << "file closed: " << ID << endl;
+}
+
+int GridFile::getId()
+{
+ return ID;
+}
+
+// memory map
+int GridFile::map()
+{
+ //filedes = fileno(fp);
+ if((filedes = open(fname, O_RDWR)) < 0)
+ {
+ fprintf(stderr, "%s open error %d(%s)\n", fname, errno, strerror(errno));
+ }
+
+ if((interp = (GridPoint *)mmap(0, sizeof(GridPoint) * size_x * size_y, PROT_READ | PROT_WRITE, MAP_SHARED, filedes, 0)) == MAP_FAILED)
+ {
+ fprintf(stderr, "mmap error %d(%s) \n", errno, strerror(errno));
+ return -1;
+ }
+ inMemory = true;
+
+ return 0;
+
+ //cout << "map size: " << sizeof(GridPoint) * size_x * size_y << endl;
+ //cout << "sizeof(GridPoint): " << sizeof(GridPoint) << endl;
+ //cout << "mmaped: " << ID << endl;
+}
+
+int GridFile::unmap()
+{
+ int rc;
+
+ // have to delete previous information??
+ // have to overwrite
+
+ // we can track the changes but that scheme requires memory usage.
+ // But our main goal is to fully utilize memory.
+
+ if(interp != NULL)
+ {
+ rc = munmap(interp, sizeof(GridPoint) * size_x * size_y);
+
+ if(rc < 0)
+ {
+ fprintf(stderr, "munmap error %d(%s)\n", errno, strerror(errno));
+ return -1;
+ }
+
+ inMemory = false;
+ interp = NULL;
+ close(filedes);
+ filedes = -1;
+
+ //cout << "unmapped: " << ID << endl;
+ }
+
+ return 0;
+}
+
+bool GridFile::isInMemory()
+{
+ return inMemory;
+}
+
+unsigned int GridFile::getMemSize()
+{
+ return size_x * size_y * sizeof(GridPoint);
+}
81 GridFile.h
@@ -0,0 +1,81 @@
+/*
+*
+COPYRIGHT AND LICENSE
+
+Copyright (c) 2011 The Regents of the University of California.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or
+without modification, are permitted provided that the following
+conditions are met:
+
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above
+copyright notice, this list of conditions and the following
+disclaimer in the documentation and/or other materials provided
+with the distribution.
+
+3. All advertising materials mentioning features or use of this
+software must display the following acknowledgement: This product
+includes software developed by the San Diego Supercomputer Center.
+
+4. Neither the names of the Centers nor the names of the contributors
+may be used to endorse or promote products derived from this
+software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS''
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS
+OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
+USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
+AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
+ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+POSSIBILITY OF SUCH DAMAGE.
+*
+*
+* Based on the notes by Prof. Ramon Arrowsmith(ramon.arrowsmith@asu.edu)
+* Authors: Han S Kim (hskim@cs.ucsd.edu), Sriram Krishnan (sriram@sdsc.edu)
+*
+*/
+
+#ifndef _GRID_FILE_H_
+#define _GRID_FILE_H_
+
+using namespace std;
+
+#include <iostream>
+#include <fstream>
+#include "GridPoint.h"
+
+class GridFile
+{
+ public:
+ GridFile(int id, char *_fname, int _size_x, int _size_y);
+ ~GridFile();
+
+ int getId();
+ int map();
+ int unmap();
+ bool isInMemory();
+ unsigned int getMemSize();
+
+ GridPoint *interp;
+
+ private:
+ //ofstream fout;
+
+ int filedes;
+ int ID;
+ int size_x;
+ int size_y;
+ bool inMemory;
+ char fname[1024];
+};
+
+#endif
131 GridMap.cpp
@@ -0,0 +1,131 @@
+/*
+*
+COPYRIGHT AND LICENSE
+
+Copyright (c) 2011 The Regents of the University of California.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or
+without modification, are permitted provided that the following
+conditions are met:
+
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above
+copyright notice, this list of conditions and the following
+disclaimer in the documentation and/or other materials provided
+with the distribution.
+
+3. All advertising materials mentioning features or use of this
+software must display the following acknowledgement: This product
+includes software developed by the San Diego Supercomputer Center.
+
+4. Neither the names of the Centers nor the names of the contributors
+may be used to endorse or promote products derived from this
+software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS''
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS
+OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
+USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
+AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
+ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+POSSIBILITY OF SUCH DAMAGE.
+*
+*
+* Based on the notes by Prof. Ramon Arrowsmith(ramon.arrowsmith@asu.edu)
+* Authors: Han S Kim (hskim@cs.ucsd.edu), Sriram Krishnan (sriram@sdsc.edu)
+*
+*/
+
+
+#include "GridMap.h"
+
+
+GridMap::GridMap(int _id, int _size_x, int _lower_bound, int _upper_bound, int _overlap_lower_bound, int _overlap_upper_bound, bool _initialized, char *fname)
+{
+ id = _id;
+ lowerBound = _lower_bound;
+ upperBound = _upper_bound;
+ overlapLowerBound = _overlap_lower_bound;
+ overlapUpperBound = _overlap_upper_bound;
+ initialized = _initialized;
+
+ gridFile = new GridFile(id, fname, _size_x, _overlap_upper_bound - _overlap_lower_bound + 1);
+}
+GridMap::~GridMap()
+{
+ if(gridFile != NULL)
+ delete gridFile;
+}
+
+int GridMap::getLowerBound()
+{
+ return lowerBound;
+}
+int GridMap::getUpperBound()
+{
+ return upperBound;
+}
+int GridMap::getOverlapLowerBound()
+{
+ return overlapLowerBound;
+}
+int GridMap::getOverlapUpperBound()
+{
+ return overlapUpperBound;
+}
+
+bool GridMap::isInitialized()
+{
+ return initialized;
+}
+int GridMap::getId()
+{
+ return id;
+}
+GridFile *GridMap::getGridFile()
+{
+ return gridFile;
+}
+
+void GridMap::setLowerBound(int _lower_bound)
+{
+ lowerBound = _lower_bound;
+}
+void GridMap::setUpperBound(int _upper_bound)
+{
+ upperBound = _upper_bound;
+}
+void GridMap::setOverlapLowerBound(int _overlap_lower_bound)
+{
+ overlapLowerBound = _overlap_lower_bound;
+}
+void GridMap::setOverlapUpperBound(int _overlap_upper_bound)
+{
+ overlapUpperBound = _overlap_upper_bound;
+}
+void GridMap::setInitialized(bool _initialized)
+{
+ initialized = _initialized;
+}
+void GridMap::setId(int _id)
+{
+ id = _id;
+}
+/*
+void GridMap::setGridFile(string fname)
+{
+ if(gridFile != NULL)
+ delete gridFile;
+
+ gridFile = new GridFile(id, fname, _size_x, _overlap_upper_bound - _overlap_lower_bound + 1);
+}
+*/
+
100 GridMap.h
@@ -0,0 +1,100 @@
+/*
+*
+COPYRIGHT AND LICENSE
+
+Copyright (c) 2011 The Regents of the University of California.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or
+without modification, are permitted provided that the following
+conditions are met:
+
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above
+copyright notice, this list of conditions and the following
+disclaimer in the documentation and/or other materials provided
+with the distribution.
+
+3. All advertising materials mentioning features or use of this
+software must display the following acknowledgement: This product
+includes software developed by the San Diego Supercomputer Center.
+
+4. Neither the names of the Centers nor the names of the contributors
+may be used to endorse or promote products derived from this
+software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS''
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS
+OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
+USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
+AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
+ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+POSSIBILITY OF SUCH DAMAGE.
+*
+*
+* Based on the notes by Prof. Ramon Arrowsmith(ramon.arrowsmith@asu.edu)
+* Authors: Han S Kim (hskim@cs.ucsd.edu), Sriram Krishnan (sriram@sdsc.edu)
+*
+*/
+
+#ifndef _GRID_MAP_H_
+#define _GRID_MAP_H_
+
+using namespace std;
+
+#include <iostream>
+#include <fstream>
+#include <stdlib.h>
+#include <string>
+#include "GridFile.h"
+
+class GridMap
+{
+ public:
+ GridMap(int _id,
+ int _size_x,
+ int _lower_bound,
+ int _upper_bound,
+ int _overlap_lower_bound,
+ int _overlap_uuper_bound,
+ bool _initialized,
+ char * fname);
+ ~GridMap();
+
+ public:
+ int getId();
+ int getLowerBound();
+ int getUpperBound();
+ int getOverlapLowerBound();
+ int getOverlapUpperBound();
+ GridFile *getGridFile();
+
+ bool isInitialized();
+
+ void setId(int _id);
+ void setLowerBound(int _lower_bound);
+ void setUpperBound(int _upper_bound);
+ void setOverlapLowerBound(int _overlap_lower_bound);
+ void setOverlapUpperBound(int _overlap_upper_bound);
+ void setInitialized(bool _initialized);
+ //void setGridFile(string fname);
+
+ private:
+ int lowerBound;
+ int upperBound;
+ int overlapLowerBound;
+ int overlapUpperBound;
+
+ bool initialized;
+ int id;
+ GridFile * gridFile;
+};
+
+#endif
62 GridPoint.h
@@ -0,0 +1,62 @@
+/*
+*
+COPYRIGHT AND LICENSE
+
+Copyright (c) 2011 The Regents of the University of California.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or
+without modification, are permitted provided that the following
+conditions are met:
+
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above
+copyright notice, this list of conditions and the following
+disclaimer in the documentation and/or other materials provided
+with the distribution.
+
+3. All advertising materials mentioning features or use of this
+software must display the following acknowledgement: This product
+includes software developed by the San Diego Supercomputer Center.
+
+4. Neither the names of the Centers nor the names of the contributors
+may be used to endorse or promote products derived from this
+software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS''
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS
+OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
+USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
+AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
+ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+POSSIBILITY OF SUCH DAMAGE.
+*
+*
+* Based on the notes by Prof. Ramon Arrowsmith(ramon.arrowsmith@asu.edu)
+* Authors: Han S Kim (hskim@cs.ucsd.edu), Sriram Krishnan (sriram@sdsc.edu)
+*
+*/
+
+#ifndef _GRID_POINT_H_
+#define _GRID_POINT_H_
+
+typedef struct _GridPoint
+{
+ double Zmin;
+ double Zmax;
+ double Zmean;
+ unsigned int count;
+ double Zidw;
+ double sum;
+ int empty;
+ int filled;
+} GridPoint;
+
+#endif
34 INSTALL
@@ -0,0 +1,34 @@
+PREREQUISITES
+-------------
+Points2Grid requires two libraries to be installed -
+
+1) curl: Please install libcurl from http://curl.haxx.se/. Points2Grid has
+been tested to work with version 7.20.0 of libcurl.
+
+2) liblas: Please install liblas from http://liblas.org. Points2Grid has
+been tested to work with version 1.2.1 of liblas.
+
+CONFIGURATION
+-------------
+Please adjust the MEM_LIMIT variable in the Interpolation.h file to reflect
+the memory limitations of your system. As a rule of thumb, set the
+MEM_LIMIT to be (available memory in bytes)/55. For jobs that generate
+grids cells greater than the MEM_LIMIT, points2grid goes "out-of-core",
+resulting in significantly worse performance.
+
+COMPILATION
+-----------
+
+If you have installed curl and liblas at non-standard locations, please set
+your CPLUS_INCLUDE_PATH and LIBRARY_PATH to respectively point to the
+locations of the include files and libraries of the installed software
+pre-requisites.
+
+If you are on a Mac (OSX), edit the Makefile to uncomment the Mac-specific
+CFLAGS setting. If you are on Solaris or Linux, there is no need to update
+the Makefile.
+
+To compile points2grid, simply use the Makefile by typing "make". The
+Makefile has been confirmed to work with gcc3 and gcc4 on Linux, Solaris
+and OSX.
+
699 InCoreInterp.cpp
@@ -0,0 +1,699 @@
+/*
+*
+COPYRIGHT AND LICENSE
+
+Copyright (c) 2011 The Regents of the University of California.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or
+without modification, are permitted provided that the following
+conditions are met:
+
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above
+copyright notice, this list of conditions and the following
+disclaimer in the documentation and/or other materials provided
+with the distribution.
+
+3. All advertising materials mentioning features or use of this
+software must display the following acknowledgement: This product
+includes software developed by the San Diego Supercomputer Center.
+
+4. Neither the names of the Centers nor the names of the contributors
+may be used to endorse or promote products derived from this
+software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS''
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS
+OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
+USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
+AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
+ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+POSSIBILITY OF SUCH DAMAGE.
+*
+*
+* Based on the notes by Prof. Ramon Arrowsmith(ramon.arrowsmith@asu.edu)
+* Authors: Han S Kim (hskim@cs.ucsd.edu), Sriram Krishnan (sriram@sdsc.edu)
+*
+*/
+
+#include "InCoreInterp.h"
+
+#include <float.h>
+#include <math.h>
+#include "Interpolation.h"
+#include "Global.h"
+#include <time.h>
+#include <sys/times.h>
+#include <stdio.h>
+#include <string.h>
+
+InCoreInterp::InCoreInterp(double dist_x, double dist_y,
+ int size_x, int size_y,
+ double r_sqr,
+ double _min_x, double _max_x,
+ double _min_y, double _max_y,
+ int _window_size)
+{
+ GRID_DIST_X = dist_x;
+ GRID_DIST_Y = dist_y;
+
+ GRID_SIZE_X = size_x;
+ GRID_SIZE_Y = size_y;
+
+ radius_sqr = r_sqr;
+
+ min_x = _min_x;
+ max_x = _max_x;
+ min_y = _min_y;
+ max_y = _max_y;
+
+ window_size = _window_size;
+
+ cout << "InCoreInterp created successfully" << endl;
+}
+
+InCoreInterp::~InCoreInterp()
+{
+ int i;
+ for(i = 0; i < GRID_SIZE_X; i++)
+ delete interp[i];
+ delete interp;
+}
+
+int InCoreInterp::init()
+{
+ int i, j;
+
+ interp = (GridPoint**)malloc(sizeof(GridPoint *) * GRID_SIZE_X);
+ //interp = new GridPoint*[GRID_SIZE_X];
+ if(interp == NULL)
+ {
+ cout << "InCoreInterp::init() new allocate error" << endl;
+ return -1;
+ }
+
+ for(i = 0; i < GRID_SIZE_X; i++)
+ {
+ interp[i] = (GridPoint *)malloc(sizeof(GridPoint) * GRID_SIZE_Y);
+ //interp[i] = new GridPoint[GRID_SIZE_Y];
+ if(interp[i] == NULL)
+ {
+ cout << "InCoreInterp::init() new allocate error" << endl;
+ return -1;
+ }
+ }
+
+ for(i = 0; i < GRID_SIZE_X; i++)
+ for(j = 0; j < GRID_SIZE_Y; j++)
+ {
+ interp[i][j].Zmin = DBL_MAX;
+ interp[i][j].Zmax = -DBL_MAX;
+ interp[i][j].Zmean = 0;
+ interp[i][j].count = 0;
+ interp[i][j].Zidw = 0;
+ interp[i][j].sum = 0;
+ interp[i][j].empty = 0;
+ interp[i][j].filled = 0;
+ }
+
+ cout << "InCoreInterp::init() done" << endl;
+
+ return 0;
+}
+
+int InCoreInterp::update(double data_x, double data_y, double data_z)
+{
+ double x;
+ double y;
+
+ int lower_grid_x;
+ int lower_grid_y;
+
+ //cout << GRID_DIST_X << " " << GRID_DIST_Y;
+ lower_grid_x = (int)floor((double)data_x/GRID_DIST_X);
+ lower_grid_y = (int)floor((double)data_y/GRID_DIST_Y);
+
+ if(lower_grid_x > GRID_SIZE_X || lower_grid_y > GRID_SIZE_Y)
+ {
+ cout << "larger at (" << lower_grid_x << "," << lower_grid_y << "): ("<< data_x << ", " << data_y << ")" << endl;
+ return 0;
+ }
+
+ //printf("lower_grid_x: %d, grid_y: %d, arrX: %.2f, arrY: %.2f\n", lower_grid_x, lower_grid_y, arrX[i], arrY[i]);
+ x = (data_x - (lower_grid_x) * GRID_DIST_X);
+ y = (data_y - (lower_grid_y) * GRID_DIST_Y);
+
+ //cout << "(" << data_x << " " << data_y << ")=(" << lower_grid_x << ", " << lower_grid_y << "): ";
+ //printf("(%f %f) = (%d, %d): ", data_x, data_y, lower_grid_x, lower_grid_y);
+
+ //if(lower_grid_y == 30 && data_y > GRID_DIST_Y * lower_grid_y)
+ //printf("(%f %f) = (%d, %d)\n", data_x, data_y, lower_grid_x, lower_grid_y);
+
+ update_first_quadrant(data_z, lower_grid_x+1, lower_grid_y+1, GRID_DIST_X - x, GRID_DIST_Y - y);
+ update_second_quadrant(data_z, lower_grid_x, lower_grid_y+1, x, GRID_DIST_Y - y);
+ update_third_quadrant(data_z, lower_grid_x, lower_grid_y, x, y);
+ update_fourth_quadrant(data_z, lower_grid_x+1, lower_grid_y, GRID_DIST_X - x, y);
+
+ //cout << "test" << endl;
+ return 0;
+}
+
+int InCoreInterp::finish(char *outputName, int outputFormat, unsigned int outputType)
+{
+ int rc;
+ int i,j;
+
+ //struct tms tbuf;
+ clock_t t0, t1;
+
+ for(i = 0; i < GRID_SIZE_X; i++)
+ for(j = 0; j < GRID_SIZE_Y; j++)
+ {
+ if(interp[i][j].Zmin == DBL_MAX){
+ // interp[i][j].Zmin = NAN;
+ interp[i][j].Zmin = 0;
+ }
+
+ if(interp[i][j].Zmax == -DBL_MAX){
+ //interp[i][j].Zmax = NAN;
+ interp[i][j].Zmax = 0;
+ }
+
+ if(interp[i][j].count != 0) {
+ interp[i][j].Zmean /= interp[i][j].count;
+ interp[i][j].empty = 1;
+ } else{
+ //interp[i][j].Zmean = NAN;
+ interp[i][j].Zmean = 0;
+ }
+
+ if(interp[i][j].sum != 0 && interp[i][j].sum != -1)
+ interp[i][j].Zidw /= interp[i][j].sum;
+ else if (interp[i][j].sum == -1){
+ // do nothing
+ } else{
+ //interp[i][j].Zidw = NAN;
+ interp[i][j].Zidw = 0;
+ }
+ }
+
+ // Sriram's edit: Fill zeros using the window size parameter
+ if (window_size != 0) {
+ int window_dist = window_size / 2;
+ for (int i = 0; i < GRID_SIZE_X; i++)
+ for (int j = 0; j < GRID_SIZE_Y; j++)
+ {
+ if (interp[i][j].empty == 0) {
+ double new_sum=0.0;
+ for (int p = i - window_dist; p <= i + window_dist; p++) {
+ for (int q = j - window_dist; q <= j + window_dist; q++) {
+ if ((p >= 0) && (p < GRID_SIZE_X) && (q >=0) && (q < GRID_SIZE_Y)) {
+ if ((p == i) && (q == j))
+ continue;
+
+ if (interp[p][q].empty != 0) {
+ double distance = max(fabs(p-i), fabs(q-j));
+ interp[i][j].Zmean += interp[p][q].Zmean/(pow(distance,Interpolation::WEIGHTER));
+ interp[i][j].Zidw += interp[p][q].Zidw/(pow(distance,Interpolation::WEIGHTER));
+ interp[i][j].Zmin += interp[p][q].Zmin/(pow(distance,Interpolation::WEIGHTER));
+ interp[i][j].Zmax += interp[p][q].Zmax/(pow(distance,Interpolation::WEIGHTER));
+
+ new_sum += 1/(pow(distance,Interpolation::WEIGHTER));
+ }
+ }
+ }
+ }
+ if (new_sum > 0) {
+ interp[i][j].Zmean /= new_sum;
+ interp[i][j].Zidw /= new_sum;
+ interp[i][j].Zmin /= new_sum;
+ interp[i][j].Zmax /= new_sum;
+ interp[i][j].filled = 1;
+ }
+ }
+ }
+ }
+
+ t0 = clock();
+ //t0 = times(&tbuf);
+
+ if((rc = outputFile(outputName, outputFormat, outputType)) < 0)
+ {
+ cout << "InCoreInterp::finish outputFile error" << endl;
+ return -1;
+ }
+
+ t1 = clock();
+ //t1 = times(&tbuf);
+
+ printf("Output Execution time: %10.2f\n", (double)(t1 - t0)/ CLOCKS_PER_SEC);
+
+
+ return 0;
+}
+
+//////////////////////////////////////////////////////
+// Private Methods
+//////////////////////////////////////////////////////
+
+void InCoreInterp::update_first_quadrant(double data_z, int base_x, int base_y, double x, double y)
+{
+ int i;
+ int j;
+ //double temp;
+
+ //printf("radius: %f ", radius_sqrt);
+
+ for(i = base_x; i < GRID_SIZE_X; i++)
+ {
+ for(j = base_y; j < GRID_SIZE_Y; j++)
+ {
+ /*
+ temp = ( ((i - base_x)*GRID_DIST + x) * ((i - base_x)*GRID_DIST + x) +
+ ((j - base_y)*GRID_DIST + y) * ((j - base_y)*GRID_DIST + y)) ;
+ printf("%f ", temp);
+ */
+
+ double distance = ((i - base_x)*GRID_DIST_X + x) * ((i - base_x)*GRID_DIST_X + x) +
+ ((j - base_y)*GRID_DIST_Y + y) * ((j - base_y)*GRID_DIST_Y + y) ;
+
+ if(distance <= radius_sqr)
+ {
+ //printf("(%d %d) ", i, j);
+ //interp[i][j]++;
+
+ // update GridPoint
+ updateGridPoint(i, j, data_z, distance);
+
+ } else if(j == base_y) {
+ //printf("return ");
+ return;
+ } else {
+ //printf("break ");
+ break;
+ }
+ }
+ }
+
+ //cout << "test2" << endl;
+}
+
+
+void InCoreInterp::update_second_quadrant(double data_z, int base_x, int base_y, double x, double y)
+{
+ int i;
+ int j;
+
+ for(i = base_x; i >= 0; i--)
+ {
+ for(j = base_y; j < GRID_SIZE_Y; j++)
+ {
+ double distance = ((base_x - i)*GRID_DIST_X + x) * ((base_x - i)*GRID_DIST_X + x) +
+ ((j - base_y)*GRID_DIST_Y + y) * ((j - base_y)*GRID_DIST_Y + y);
+
+ if(distance <= radius_sqr)
+ {
+ //printf("(%d %d) ", i, j);
+ //interp[i][j]++;
+
+ updateGridPoint(i, j, data_z, sqrt(distance));
+
+ } else if(j == base_y){
+ return;
+ } else {
+ break;
+ }
+ }
+ }
+
+}
+
+
+void InCoreInterp::update_third_quadrant(double data_z, int base_x, int base_y, double x, double y)
+{
+ int i;
+ int j;
+
+ for(i = base_x; i >= 0; i--)
+ {
+ for(j = base_y; j >= 0; j--)
+ {
+ double distance = ((base_x - i)*GRID_DIST_X + x) * ((base_x - i)*GRID_DIST_X + x) +
+ ((base_y - j)*GRID_DIST_Y + y) * ((base_y - j)*GRID_DIST_Y + y);
+
+ if(distance <= radius_sqr)
+ {
+ //if(j == 30)
+ //printf("(%d %d)\n", i, j);
+ //interp[i][j]++;
+ updateGridPoint(i, j, data_z, sqrt(distance));
+ } else if(j == base_y){
+ return;
+ } else {
+ break;
+ }
+ }
+ }
+}
+
+void InCoreInterp::update_fourth_quadrant(double data_z, int base_x, int base_y, double x, double y)
+{
+ int i, j;
+
+ for(i = base_x; i < GRID_SIZE_X; i++)
+ {
+ for(j = base_y; j >= 0; j--)
+ {
+ double distance = ((i - base_x)*GRID_DIST_X + x) * ((i - base_x)*GRID_DIST_X + x) +
+ ((base_y - j)*GRID_DIST_Y + y) * ((base_y - j)*GRID_DIST_Y + y);
+
+ if(distance <= radius_sqr)
+ {
+ //printf("(%d %d) ", i, j);
+ //interp[i][j]++;
+ updateGridPoint(i, j, data_z, sqrt(distance));
+ } else if (j == base_y) {
+ return ;
+ } else {
+ break;
+ }
+ }
+ }
+}
+
+void InCoreInterp::updateGridPoint(int x, int y, double data_z, double distance)
+{
+ if(interp[x][y].Zmin > data_z)
+ interp[x][y].Zmin = data_z;
+ if(interp[x][y].Zmax < data_z)
+ interp[x][y].Zmax = data_z;
+
+ interp[x][y].Zmean += data_z;
+ interp[x][y].count++;
+
+ double dist = pow(distance, Interpolation::WEIGHTER);
+
+ if(interp[x][y].sum != -1) {
+ if(dist != 0) {
+ interp[x][y].Zidw += data_z/dist;
+ interp[x][y].sum += 1/dist;
+ } else {
+ interp[x][y].Zidw = data_z;
+ interp[x][y].sum = -1;
+ }
+ } else {
+ // do nothing
+ }
+}
+
+void InCoreInterp::printArray()
+{
+ int i, j;
+
+ for(i = 0; i < GRID_SIZE_X; i++)
+ {
+ for(j = 1; j < GRID_SIZE_Y; j++)
+ {
+ cout << interp[i][j].Zmin << ", " << interp[i][j].Zmax << ", ";
+ cout << interp[i][j].Zmean << ", " << interp[i][j].Zidw << endl;
+ //printf("%.20f ", interp[i][j].Zmax);
+ }
+ //printf("\n");
+ }
+ cout << endl;
+}
+
+int InCoreInterp::outputFile(char *outputName, int outputFormat, unsigned int outputType)
+{
+ int i,j,k;
+
+ FILE **arcFiles;
+ char arcFileName[1024];
+
+ FILE **gridFiles;
+ char gridFileName[1024];
+
+ char *ext[5] = {".min", ".max", ".mean", ".idw", ".den"};
+ unsigned int type[5] = {OUTPUT_TYPE_MIN, OUTPUT_TYPE_MAX, OUTPUT_TYPE_MEAN, OUTPUT_TYPE_IDW, OUTPUT_TYPE_DEN};
+ int numTypes = 5;
+
+
+
+ // open ArcGIS files
+ if(outputFormat == OUTPUT_FORMAT_ARC_ASCII || outputFormat == OUTPUT_FORMAT_ALL)
+ {
+ if((arcFiles = (FILE **)malloc(sizeof(FILE *) * numTypes)) == NULL)
+ {
+ cout << "Arc File open error: " << endl;
+ return -1;
+ }
+
+ for(i = 0; i < numTypes; i++)
+ {
+ if(outputType & type[i])
+ {
+ strncpy(arcFileName, outputName, sizeof(arcFileName));
+ strncat(arcFileName, ext[i], strlen(ext[i]));
+ strncat(arcFileName, ".asc", strlen(".asc"));
+
+ if((arcFiles[i] = fopen64(arcFileName, "w+")) == NULL)
+ {
+ cout << "File open error: " << arcFileName << endl;
+ return -1;
+ }
+ }else{
+ arcFiles[i] = NULL;
+ }
+ }
+ } else {
+ arcFiles = NULL;
+ }
+
+ // open Grid ASCII files
+ if(outputFormat == OUTPUT_FORMAT_GRID_ASCII || outputFormat == OUTPUT_FORMAT_ALL)
+ {
+ if((gridFiles = (FILE **)malloc(sizeof(FILE *) * numTypes)) == NULL)
+ {
+ cout << "File array allocation error" << endl;
+ return -1;
+ }
+
+ for(i = 0; i < numTypes; i++)
+ {
+ if(outputType & type[i])
+ {
+ strncpy(gridFileName, outputName, sizeof(arcFileName));
+ strncat(gridFileName, ext[i], strlen(ext[i]));
+ strncat(gridFileName, ".grid", strlen(".grid"));
+
+ if((gridFiles[i] = fopen64(gridFileName, "w+")) == NULL)
+ {
+ cout << "File open error: " << gridFileName << endl;
+ return -1;
+ }
+ }else{
+ gridFiles[i] = NULL;
+ }
+ }
+ } else {
+ gridFiles = NULL;
+ }
+
+
+ //struct tms tbuf;
+ //clock_t t0, t1;
+
+ //t0 = times(&tbuf);
+
+ // print ArcGIS headers
+ if(arcFiles != NULL)
+ {
+ for(i = 0; i < numTypes; i++)
+ {
+ if(arcFiles[i] != NULL)
+ {
+ fprintf(arcFiles[i], "ncols %d\n", GRID_SIZE_X);
+ fprintf(arcFiles[i], "nrows %d\n", GRID_SIZE_Y);
+ fprintf(arcFiles[i], "xllcorner %f\n", min_x);
+ fprintf(arcFiles[i], "yllcorner %f\n", min_y);
+ fprintf(arcFiles[i], "cellsize %f\n", GRID_DIST_X);
+ fprintf(arcFiles[i], "NODATA_value -9999\n");
+ }
+ }
+ }
+
+ // print Grid headers
+ if(gridFiles != NULL)
+ {
+ for(i = 0; i < numTypes; i++)
+ {
+ if(gridFiles[i] != NULL)
+ {
+ fprintf(gridFiles[i], "north: %f\n", max_y);
+ fprintf(gridFiles[i], "south: %f\n", min_y);
+ fprintf(gridFiles[i], "east: %f\n", max_x);
+ fprintf(gridFiles[i], "west: %f\n", min_x);
+ fprintf(gridFiles[i], "rows: %d\n", GRID_SIZE_Y);
+ fprintf(gridFiles[i], "cols: %d\n", GRID_SIZE_X);
+ }
+ }
+ }
+
+ // print data
+ for(i = GRID_SIZE_Y - 1; i >= 0; i--)
+ {
+ for(j = 0; j < GRID_SIZE_X; j++)
+ {
+ if(arcFiles != NULL)
+ {
+ // Zmin
+ if(arcFiles[0] != NULL)
+ {
+ if(interp[j][i].empty == 0 &&
+ interp[j][i].filled == 0)
+ fprintf(arcFiles[0], "-9999 ");
+ else
+ fprintf(arcFiles[0], "%f ", interp[j][i].Zmin);
+ }
+
+ // Zmax
+ if(arcFiles[1] != NULL)
+ {
+ if(interp[j][i].empty == 0 &&
+ interp[j][i].filled == 0)
+ fprintf(arcFiles[1], "-9999 ");
+ else
+ fprintf(arcFiles[1], "%f ", interp[j][i].Zmax);
+ }
+
+ // Zmean
+ if(arcFiles[2] != NULL)
+ {
+ if(interp[j][i].empty == 0 &&
+ interp[j][i].filled == 0)
+ fprintf(arcFiles[2], "-9999 ");
+ else
+ fprintf(arcFiles[2], "%f ", interp[j][i].Zmean);
+ }
+
+ // Zidw
+ if(arcFiles[3] != NULL)
+ {
+ if(interp[j][i].empty == 0 &&
+ interp[j][i].filled == 0)
+ fprintf(arcFiles[3], "-9999 ");
+ else
+ fprintf(arcFiles[3], "%f ", interp[j][i].Zidw);
+ }
+
+ // count
+ if(arcFiles[4] != NULL)
+ {
+ if(interp[j][i].empty == 0 &&
+ interp[j][i].filled == 0)
+ fprintf(arcFiles[4], "-9999 ");
+ else
+ fprintf(arcFiles[4], "%d ", interp[j][i].count);
+ }
+ }
+
+ if(gridFiles != NULL)
+ {
+ // Zmin
+ if(gridFiles[0] != NULL)
+ {
+ if(interp[j][i].empty == 0 &&
+ interp[j][i].filled == 0)
+ fprintf(gridFiles[0], "-9999 ");
+ else
+ fprintf(gridFiles[0], "%f ", interp[j][i].Zmin);
+ }
+
+ // Zmax
+ if(gridFiles[1] != NULL)
+ {
+ if(interp[j][i].empty == 0 &&
+ interp[j][i].filled == 0)
+ fprintf(gridFiles[1], "-9999 ");
+ else
+ fprintf(gridFiles[1], "%f ", interp[j][i].Zmax);
+ }
+
+ // Zmean
+ if(gridFiles[2] != NULL)
+ {
+ if(interp[j][i].empty == 0 &&
+ interp[j][i].filled == 0)
+ fprintf(gridFiles[2], "-9999 ");
+ else
+ fprintf(gridFiles[2], "%f ", interp[j][i].Zmean);
+ }
+
+ // Zidw
+ if(gridFiles[3] != NULL)
+ {
+ if(interp[j][i].empty == 0 &&
+ interp[j][i].filled == 0)
+ fprintf(gridFiles[3], "-9999 ");
+ else
+ fprintf(gridFiles[3], "%f ", interp[j][i].Zidw);
+ }
+
+ // count
+ if(gridFiles[4] != NULL)
+ {
+ if(interp[j][i].empty == 0 &&
+ interp[j][i].filled == 0)
+ fprintf(gridFiles[4], "-9999 ");
+ else
+ fprintf(gridFiles[4], "%d ", interp[j][i].count);
+ }
+ }
+ }
+ if(arcFiles != NULL)
+ for(k = 0; k < numTypes; k++)
+ {
+ if(arcFiles[k] != NULL)
+ fprintf(arcFiles[k], "\n");
+ }
+ if(gridFiles != NULL)
+ for(k = 0; k < numTypes; k++)
+ {
+ if(gridFiles[k] != NULL)
+ fprintf(gridFiles[k], "\n");
+ }
+ }
+
+ // close files
+ if(gridFiles != NULL)
+ {
+ for(i = 0; i < numTypes; i++)
+ {
+ if(gridFiles[i] != NULL)
+ fclose(gridFiles[i]);
+ }
+ }
+
+ if(arcFiles != NULL)
+ {
+ for(i = 0; i < numTypes; i++)
+ {
+ if(arcFiles[i] != NULL)
+ fclose(arcFiles[i]);
+ }
+ }
+
+ //t1 = times(&tbuf);
+ //printf("elapsed time: %10.10f\n", (t1 - t0)/(double) sysconf(_SC_CLK_TCK));
+
+ return 0;
+}
+
87 InCoreInterp.h
@@ -0,0 +1,87 @@
+/*
+*
+COPYRIGHT AND LICENSE
+
+Copyright (c) 2011 The Regents of the University of California.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or
+without modification, are permitted provided that the following
+conditions are met:
+
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above
+copyright notice, this list of conditions and the following
+disclaimer in the documentation and/or other materials provided
+with the distribution.
+
+3. All advertising materials mentioning features or use of this
+software must display the following acknowledgement: This product
+includes software developed by the San Diego Supercomputer Center.
+
+4. Neither the names of the Centers nor the names of the contributors
+may be used to endorse or promote products derived from this
+software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS''
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS
+OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
+USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
+AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
+ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+POSSIBILITY OF SUCH DAMAGE.
+*
+*
+* Based on the notes by Prof. Ramon Arrowsmith(ramon.arrowsmith@asu.edu)
+* Authors: Han S Kim (hskim@cs.ucsd.edu), Sriram Krishnan (sriram@sdsc.edu)
+*
+*/
+
+#ifndef _IN_CORE_INTERP_H_
+#define _IN_CORE_INTERP_H_
+
+using namespace std;
+
+#include <iostream>
+#include "GridPoint.h"
+#include "CoreInterp.h"
+
+class InCoreInterp : public CoreInterp
+{
+ public:
+ InCoreInterp() {};
+ InCoreInterp(double dist_x, double dist_y,
+ int size_x, int size_y,
+ double r_sqr,
+ double _min_x, double _max_x,
+ double _min_y, double _max_y,
+ int _window_size);
+ ~InCoreInterp();
+
+ virtual int init();
+ virtual int update(double data_x, double data_y, double data_z);
+ virtual int finish(char *outputName, int outputFormat, unsigned int outputType);
+
+ private:
+ GridPoint **interp;
+ double radius_sqr;
+
+ private:
+ void update_first_quadrant(double data_z, int base_x, int base_y, double x, double y);
+ void update_second_quadrant(double data_z, int base_x, int base_y, double x, double y);
+ void update_third_quadrant(double data_z, int base_x, int base_y, double x, double y);
+ void update_fourth_quadrant(double data_z, int base_x, int base_y, double x, double y);
+
+ void updateGridPoint(int x, int y, double data_z, double distance);
+ void printArray();
+ int outputFile(char *outputName, int outputFormat, unsigned int outputType);
+};
+
+#endif
395 Interpolation.cpp
@@ -0,0 +1,395 @@
+/*
+*
+COPYRIGHT AND LICENSE
+
+Copyright (c) 2011 The Regents of the University of California.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or
+without modification, are permitted provided that the following
+conditions are met:
+
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above
+copyright notice, this list of conditions and the following
+disclaimer in the documentation and/or other materials provided
+with the distribution.
+
+3. All advertising materials mentioning features or use of this
+software must display the following acknowledgement: This product
+includes software developed by the San Diego Supercomputer Center.
+
+4. Neither the names of the Centers nor the names of the contributors
+may be used to endorse or promote products derived from this
+software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS''
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS
+OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
+USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
+AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
+ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+POSSIBILITY OF SUCH DAMAGE.
+*
+*
+* Based on the notes by Prof. Ramon Arrowsmith(ramon.arrowsmith@asu.edu)
+* Authors: Han S Kim (hskim@cs.ucsd.edu), Sriram Krishnan (sriram@sdsc.edu)
+*
+*/
+
+#include "Interpolation.h"
+#include "Global.h"
+#include <string.h>
+#include <math.h>
+#include <float.h>
+#include <stddef.h>
+#include <stdlib.h>
+#include <time.h>
+#include <sys/times.h>
+#include <stdio.h>
+#include <liblas/laspoint.hpp>
+#include <liblas/lasreader.hpp>
+#include <fstream> // std::ifstream
+#include <iostream> // std::cout
+#include <string.h>
+
+/////////////////////////////////////////////////////////////
+// Public Methods
+/////////////////////////////////////////////////////////////
+
+Interpolation::Interpolation(double x_dist, double y_dist, double radius, int _window_size) : GRID_DIST_X (x_dist), GRID_DIST_Y(y_dist)
+{
+ data_count = 0;
+ radius_sqr = radius * radius;
+ window_size = _window_size;
+
+ min_x = DBL_MAX;
+ min_y = DBL_MAX;
+
+ max_x = -DBL_MAX;
+ max_y = -DBL_MAX;
+}
+
+Interpolation::~Interpolation()
+{
+ delete interp;
+}
+
+int Interpolation::init(char *inputName, int inputFormat)
+{
+ //unsigned int i;
+
+ //struct tms tbuf;
+ clock_t t0, t1;
+
+
+ //////////////////////////////////////////////////////////////////////
+ // MIN/MAX SEARCHING
+ // TODO: the size of data, min, max of each coordinate
+ // are required to implement streaming processing....
+ //
+ // This code can be eliminated if database can provide these values
+ //////////////////////////////////////////////////////////////////////
+
+ //t0 = times(&tbuf);
+ t0 = clock();
+
+ if(inputName == NULL)
+ {
+ cout << "Wrong Input File Name" << endl;
+ return -1;
+ }
+
+ printf("inputName: '%s'\n", inputName);
+
+ if (inputFormat == INPUT_ASCII) {
+ FILE *fp;
+ char line[1024];
+ double data_x, data_y;
+ //double data_z;
+
+ if((fp = fopen64(inputName, "r")) == NULL)
+ {
+ cout << "file open error" << endl;
+ return -1;
+ }
+
+ // throw the first line away - it contains the header
+ fgets(line, sizeof(line), fp);
+
+ // read the data points to find min and max values
+ while(fgets(line, sizeof(line), fp) != NULL)
+ {
+ data_x = atof(strtok(line, ",\n"));
+ if(min_x > data_x) min_x = data_x;
+ if(max_x < data_x) max_x = data_x;
+
+ data_y = atof(strtok(NULL, ",\n"));
+ if(min_y > data_y) min_y = data_y;
+ if(max_y < data_y) max_y = data_y;
+
+ data_count++;
+
+
+ /*
+ // tokenizing
+ string strLine(line);
+
+ // first token
+ string::size_type pos = strLine.find_first_of(",",0);
+ string::size_type lastPos = strLine.find_first_not_of(",",0);
+
+ string strX = strLine.substr(lastPos, pos - lastPos);
+
+ // second token
+ lastPos = strLine.find_first_not_of(",", pos);
+ pos = strLine.find_first_of(",", lastPos);
+
+ string strY = strLine.substr(lastPos, pos - lastPos);
+
+ // third token
+ lastPos = strLine.find_first_not_of(",", pos);
+ pos = strLine.find_first_of("\n", lastPos);
+
+ string strZ = strLine.substr(lastPos, pos - lastPos);
+
+ // data conversion
+ arrX[data_count] = atof(strX.c_str());
+ if(min_x > arrX[data_count]) min_x = arrX[data_count];
+ if(max_x < arrX[data_count]) max_x = arrX[data_count];
+
+ arrY[data_count] = atof(strY.c_str());
+ if(min_y > arrY[data_count]) min_y = arrY[data_count];
+ if(max_y < arrY[data_count]) max_y = arrY[data_count];
+
+ arrZ[data_count++] = atof(strZ.c_str());
+ */
+ }
+
+ fclose(fp);
+ } else { // las input
+
+ try {
+ std::ifstream ifs;
+ ifs.open(inputName, std::ios::in | std::ios::binary);
+
+ // create a las file reader
+ liblas::LASReader reader(ifs);
+
+ /// get header information
+ liblas::LASHeader const& header = reader.GetHeader();
+
+ min_x = header.GetMinX();
+ min_y = header.GetMinY();
+ max_x = header.GetMaxX();
+ max_y = header.GetMaxY();
+ data_count = header.GetPointRecordsCount();
+
+ ifs.close();
+ } catch (std::runtime_error &e) {
+ cout << "error while reading LAS file: verify that the input is valid LAS file" << endl;
+ return -1;
+ }
+ }
+
+ t1 = clock();
+ //t1 = times(&tbuf);
+ //printf("Min/Max searching time: %10.2f\n", (t1 - t0)/(double) sysconf(_SC_CLK_TCK));
+ printf("Min/Max searching time: %10.2f\n", (double)(t1 - t0)/CLOCKS_PER_SEC);
+
+
+ //t0 = times(&tbuf);
+
+ //////////////////////////////////////////////////////////////////////
+ // Intialization Step excluding min/max searching
+ //////////////////////////////////////////////////////////////////////
+ /*
+ for(i = 0; i < data_count; i++)
+ {
+ arrX[i] -= min_x;
+ arrY[i] -= min_y;
+ //printf("%f,%f,%f\n", arrX[i] , arrY[i] ,arrZ[i]);
+ }
+ */
+
+ cout << "min_x: " << min_x << ", max_x: " << max_x << ", min_y: " << min_y << ", max_y: " << max_y << endl;
+
+ GRID_SIZE_X = (int)(ceil((max_x - min_x)/GRID_DIST_X)) + 1;
+ GRID_SIZE_Y = (int)(ceil((max_y - min_y)/GRID_DIST_Y)) + 1;
+
+ cout << "GRID_SIZE_X " << GRID_SIZE_X << endl;
+ cout << "GRID_SIZE_Y " << GRID_SIZE_Y << endl;
+
+ // if the size is too big to fit in memory,
+ // then
+ // construct out-of-core structure
+ if(GRID_SIZE_X * GRID_SIZE_Y > MEM_LIMIT)
+ {
+ core_mode = OUTCORE;
+ cout << "Using out of core interp code" << endl;;
+
+ interp = new OutCoreInterp(GRID_DIST_X, GRID_DIST_Y, GRID_SIZE_X, GRID_SIZE_Y, radius_sqr, min_x, max_x, min_y, max_y, window_size);
+ if(interp == NULL)
+ {
+ cout << "OutCoreInterp construction error" << endl;
+ return -1;
+ }
+
+ cout << "Interpolation uses out-of-core algorithm" << endl;
+
+ // else
+ // do normal
+ } else {
+ core_mode = INCORE;
+ cout << "Using incore interp code" << endl;
+
+ interp = new InCoreInterp(GRID_DIST_X, GRID_DIST_Y, GRID_SIZE_X, GRID_SIZE_Y, radius_sqr, min_x, max_x, min_y, max_y, window_size);
+
+ cout << "Interpolation uses in-core algorithm" << endl;
+ }
+
+ if(interp->init() < 0)
+ {
+ cout << "inter->init() error" << endl;
+ return -1;
+ }
+
+ cout << "Interpolation::init() done successfully" << endl;
+
+ //t1 = times(&tbuf);
+ //printf("Initializing time: %10.2f\n", (t1 - t0)/(double) sysconf(_SC_CLK_TCK));
+
+ return 0;
+}
+
+int Interpolation::interpolation(char *inputName,
+ char *outputName,
+ int inputFormat,
+ int outputFormat,
+ unsigned int outputType)
+{
+ int rc;
+ //unsigned int i;
+ double data_x, data_y;
+ double data_z;
+
+ //struct tms tbuf;
+ //clock_t t0, t1;
+
+ printf("Interpolation Starts\n");
+
+ //t0 = times(&tbuf);
+
+ //cout << "data_count: " << data_count << endl;
+
+ /*
+ if((rc = interp->init()) < 0)
+ {
+ cout << "inter->init() error" << endl;
+ return -1;
+ }
+ */
+
+ if (inputFormat == INPUT_ASCII) {
+ FILE *fp;
+ char line[1024];
+
+ if((fp = fopen64(inputName, "r")) == NULL)
+ {
+ printf("file open error\n");
+ return -1;
+ }
+
+ // throw the first line away - it contains the header
+ fgets(line, sizeof(line), fp);
+
+ // read every point and generate DEM
+ while(fgets(line, sizeof(line), fp) != NULL)
+ {
+ data_x = atof(strtok(line, ",\n"));
+ data_y = atof(strtok(NULL, ",\n"));
+ data_z = atof(strtok(NULL, ",\n"));
+
+ data_x -= min_x;
+ data_y -= min_y;
+
+ //if((rc = interp->update(arrX[i], arrY[i], arrZ[i])) < 0)
+ if((rc = interp->update(data_x, data_y, data_z)) < 0)
+ {
+ cout << "interp->update() error while processing " << endl;
+ return -1;
+ }
+ }
+
+ fclose(fp);
+ } else { // input format is LAS
+
+ // instantiate a reader for the LAS file
+ std::ifstream ifs;
+ ifs.open(inputName, std::ios::in | std::ios::binary);
+
+ liblas::LASReader reader(ifs);
+
+ // process every point in the LAS file, and generate DEM
+ while (reader.ReadNextPoint()) {
+ liblas::LASPoint const& p = reader.GetPoint();
+
+ data_x = p.GetX();
+ data_y = p.GetY();
+ data_z = p.GetZ();
+
+ data_x -= min_x;
+ data_y -= min_y;
+
+ //if((rc = interp->update(arrX[i], arrY[i], arrZ[i])) < 0)
+ if((rc = interp->update(data_x, data_y, data_z)) < 0)
+ {
+ cout << "interp->update() error while processing " << endl;
+ return -1;
+ }
+ }
+
+ ifs.close();
+ }
+
+ if((rc = interp->finish(outputName, outputFormat, outputType)) < 0)
+ {
+ cout << "interp->finish() error" << endl;
+ return -1;
+ }
+
+ cout << "Interpolation::interpolation() done successfully" << endl;
+
+ //t1 = times(&tbuf);
+ //printf("Interpolation + Output time: %10.2f\n", (t1 - t0)/(double) sysconf(_SC_CLK_TCK));
+
+ return 0;
+}
+
+void Interpolation::setRadius(double r)
+{
+ radius_sqr = r * r;
+}
+
+unsigned int Interpolation::getDataCount()
+{
+ return data_count;
+}
+
+unsigned int Interpolation::getGridSizeX()
+{
+ return GRID_SIZE_X;
+}
+
+unsigned int Interpolation::getGridSizeY()
+{
+ return GRID_SIZE_Y;
+}
+
117 Interpolation.h
@@ -0,0 +1,117 @@
+/*
+*
+COPYRIGHT AND LICENSE
+
+Copyright (c) 2011 The Regents of the University of California.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or
+without modification, are permitted provided that the following
+conditions are met:
+
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above
+copyright notice, this list of conditions and the following
+disclaimer in the documentation and/or other materials provided
+with the distribution.
+
+3. All advertising materials mentioning features or use of this
+software must display the following acknowledgement: This product
+includes software developed by the San Diego Supercomputer Center.
+
+4. Neither the names of the Centers nor the names of the contributors
+may be used to endorse or promote products derived from this
+software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS''
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS
+OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
+USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
+AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
+ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+POSSIBILITY OF SUCH DAMAGE.
+*
+*
+* Based on the notes by Prof. Ramon Arrowsmith(ramon.arrowsmith@asu.edu)
+* Authors: Han S Kim (hskim@cs.ucsd.edu), Sriram Krishnan (sriram@sdsc.edu)
+*
+*/
+
+#ifndef _INTERP_H_
+#define _INTERP_H_
+
+#include <string>
+#include <iostream>
+
+using namespace std;
+
+#include "GridPoint.h"
+#include "CoreInterp.h"
+#include "OutCoreInterp.h"
+#include "InCoreInterp.h"
+
+enum coreMode
+{
+ INCORE = 0,
+ OUTCORE
+};
+
+//class GridPoint;
+
+class Interpolation
+{
+ public:
+ Interpolation(double x_dist, double y_dist, double radius, int _window_size);
+ ~Interpolation();
+
+ int init(char *inputName, int inputFormat);
+ int interpolation(char *inputName, char *outputName, int inputFormat,
+ int outputFormat, unsigned int type);
+ unsigned int getDataCount();
+
+ unsigned int getGridSizeX();
+ unsigned int getGridSizeY();
+
+ // for debug
+ void printArray();
+
+ // depricated
+ void setRadius(double r);
+
+ public:
+ double GRID_DIST_X;
+ double GRID_DIST_Y;
+
+ static const int MAX_POINT_SIZE = 16000000;
+ static const int WEIGHTER = 2;
+
+ // update this to the maximum grid that will fit in memory
+ // as a rule of thumb, memory requirement = MEM_LIMIT*55 bytes
+ static const unsigned int MEM_LIMIT = 200000000;
+
+ private:
+ double min_x;
+ double min_y;
+ double max_x;
+ double max_y;
+
+ unsigned int GRID_SIZE_X;
+ unsigned int GRID_SIZE_Y;
+
+ unsigned int data_count;
+ double radius_sqr;
+ int window_size;
+
+ int core_mode;
+
+ CoreInterp *interp;
+};
+
+#endif
37 LICENSE
@@ -0,0 +1,37 @@
+COPYRIGHT AND LICENSE
+
+Copyright (c) 2011 The Regents of the University of California.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or
+without modification, are permitted provided that the following
+conditions are met:
+
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above
+copyright notice, this list of conditions and the following
+disclaimer in the documentation and/or other materials provided
+with the distribution.
+
+3. All advertising materials mentioning features or use of this
+software must display the following acknowledgement: This product
+includes software developed by the San Diego Supercomputer Center.
+
+4. Neither the names of the Centers nor the names of the contributors
+may be used to endorse or promote products derived from this
+software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS''
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS
+OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
+USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
+AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
+ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+POSSIBILITY OF SUCH DAMAGE.
26 Makefile
@@ -0,0 +1,26 @@
+CC=g++
+
+# flags needed to compile sources
+# ensure that your CPLUS_INCLUDE_PATH is set to pick up the required headers
+
+# on Macs
+# CFLAGS=-O3 -Wall -g -D_GNU_SOURCE -D_FILE_OFFSET_BITS=64 -D_LARGEFILE_SOURCE -D_USE_LARGEFILE64 -Dfopen64=fopen
+
+# on everything else
+CFLAGS=-O3 -m64 -Wall -g -D_GNU_SOURCE -D_FILE_OFFSET_BITS=64 -D_LARGEFILE_SOURCE -D_USE_LARGEFILE64
+
+# ensure that your LIBRARY_PATH is set up to pick up these libraries
+LIB=-llas -lm -lcurl
+
+OBJS=main.o Interpolation.o InCoreInterp.o OutCoreInterp.o GridMap.o GridFile.o
+
+all: points2grid
+
+points2grid: $(OBJS)
+ $(CC) $(CFLAGS) $(LIB) -o $@ $(OBJS)
+
+.cpp.o:
+ $(CC) $(CFLAGS) $(INC) -c $<
+
+clean:
+ rm -f points2grid *.o
1,103 OutCoreInterp.cpp
@@ -0,0 +1,1103 @@
+/*
+*
+COPYRIGHT AND LICENSE
+
+Copyright (c) 2011 The Regents of the University of California.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or
+without modification, are permitted provided that the following
+conditions are met:
+
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above
+copyright notice, this list of conditions and the following
+disclaimer in the documentation and/or other materials provided
+with the distribution.
+
+3. All advertising materials mentioning features or use of this
+software must display the following acknowledgement: This product
+includes software developed by the San Diego Supercomputer Center.
+
+4. Neither the names of the Centers nor the names of the contributors
+may be used to endorse or promote products derived from this
+software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS''
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS
+OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
+USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
+AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
+ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+POSSIBILITY OF SUCH DAMAGE.
+*
+*
+* Based on the notes by Prof. Ramon Arrowsmith(ramon.arrowsmith@asu.edu)
+* Authors: Han S Kim (hskim@cs.ucsd.edu), Sriram Krishnan (sriram@sdsc.edu)
+*
+*/
+
+#include "OutCoreInterp.h"
+
+#include <float.h>
+#include <math.h>
+#include <time.h>
+#include <sys/times.h>
+#include <sys/time.h>
+#include "Interpolation.h"
+#include "Global.h"
+#include <stdio.h>
+#include <string.h>
+
+OutCoreInterp::OutCoreInterp(double dist_x, double dist_y,
+ int size_x, int size_y,
+ double r_sqr,
+ double _min_x, double _max_x,
+ double _min_y, double _max_y,
+ int _window_size)
+{
+ int i;
+
+ GRID_DIST_X = dist_x;
+ GRID_DIST_Y = dist_y;
+
+ GRID_SIZE_X = size_x;
+ GRID_SIZE_Y = size_y;
+
+ radius_sqr = r_sqr;
+
+ min_x = _min_x;
+ max_x = _max_x;
+ min_y = _min_y;
+ max_y = _max_y;
+
+ window_size = _window_size;
+
+ overlapSize = (int)ceil(sqrt(radius_sqr)/GRID_DIST_Y);
+ int window_dist = window_size / 2;
+ if (window_dist > overlapSize) {
+ overlapSize = window_dist;
+ }
+
+ // how many pieces will there be?
+ // numFiles = (GRID_SIZE_X*GRID_SIZE_Y + (2 * overlapSize + 1) * GRID_SIZE_X * numFiles) / MEM_LIMIT;
+ // (2 * overlapSize + 1) means, each component has one more row to handle the points in-between two components.
+ numFiles = (int)ceil((double)GRID_SIZE_X * GRID_SIZE_Y / (Interpolation::MEM_LIMIT - (2 * overlapSize + 1) * GRID_SIZE_X));
+ cout << "numFiles " << numFiles << endl;
+
+ if(numFiles == 0)
+ cout << "The number of files is 0!" << endl;
+
+ // TODO: if one row is too long,
+ // i.e. the shape of decomposition is long, and thin in one direction,
+ // alternative decomposition scheme may be more efficient.
+
+ // row-wise decomposition
+
+ local_grid_size_x = GRID_SIZE_X;
+ local_grid_size_y = (int)ceil((double)GRID_SIZE_Y/numFiles);
+
+
+ // define overlap.. a funtion of (radius_sqr)
+ // overlap size: sqrt(radius_sqr)/GRID_DIST_Y;
+
+ // construct a map indicating which file corresponds which area
+ if((gridMap = new GridMap*[numFiles]) == NULL)
+ cout << "OutCoreInterp::OutCoreInterp() GridMap[] allocation error" << endl;
+
+ for(i = 0; i < numFiles; i++)
+ {
+ int lower_bound = i * local_grid_size_y;
+ int upper_bound = (i+1) * local_grid_size_y - 1;
+
+ if(upper_bound >= GRID_SIZE_Y)
+ upper_bound = GRID_SIZE_Y - 1;
+
+ int overlap_lower_bound = lower_bound - overlapSize;
+ if(overlap_lower_bound < 0)
+ overlap_lower_bound = 0;
+
+ int overlap_upper_bound = upper_bound + overlapSize + 1;
+ if(overlap_upper_bound >= GRID_SIZE_Y)
+ overlap_upper_bound = GRID_SIZE_Y - 1;
+
+ char fname[1024];
+
+ struct timeval tp;
+ gettimeofday(&tp, NULL);
+
+ snprintf(fname, sizeof(fname), "data_%d_%ld", i, tp.tv_usec);
+
+ gridMap[i] = new GridMap(i,
+ GRID_SIZE_X,
+ lower_bound,
+ upper_bound,
+ overlap_lower_bound,
+ overlap_upper_bound,
+ false,
+ fname);
+ if(gridMap[i] == NULL)
+ cout << "OutCoreInterp::OutCoreInterp() GridMap alloc error" << endl;
+
+ cout << "[" << lower_bound << "," << upper_bound << "]" << endl;
+ cout << "[" << overlap_lower_bound << "," << overlap_upper_bound << "]" << endl;
+ }
+
+ // initializing queues
+ qlist = new list<UpdateInfo> [numFiles];
+ if(qlist == NULL)
+ cout << "OutCoreInterp::OutCoreInterp() qlist alloc error" << endl;
+
+ openFile = -1;
+
+ //cout << overlapSize << endl;
+ //cout << "data size: " << (size_x * size_y) << " numFiles: " << numFiles << " local_grid_size_y: " << local_grid_size_y << " overlap: " << overlapSize << endl;
+ //cout << fname << endl;
+ //cout << "GridMap is created: " << i << endl;
+}
+
+OutCoreInterp::~OutCoreInterp()
+{
+ /*
+ if(qlist != NULL)
+ delete qlist;
+ */
+
+ for(int i = 0; i < numFiles; i++)
+ {
+ //cout << "gridMap " << i << " deleted" << endl;
+ if(gridMap[i] != NULL)
+ delete gridMap[i];
+ }
+ if(gridMap != NULL)
+ delete gridMap;
+}
+
+int OutCoreInterp::init()
+{
+ // open up a memory mapped file
+ openFile = 0;
+ gridMap[openFile]->getGridFile()->map();
+
+ return 0;
+}
+
+int OutCoreInterp::update(double data_x, double data_y, double data_z)
+{
+ // update()
+ // push the point info into an appropriate queue
+ // if the length of the queue exceeds a limit,
+ // then
+ // if another file is loaded, write back (lazy update)
+ // if the file we need has not opened,
+ // read a file from disk
+ // update the file (PROC1)
+
+ //
+ // find which file should be updated
+
+ /*
+ if(data_x < 0){
+ cout << "4. !!!! " << data_x << endl;
+ return -1;
+ }
+ */
+
+ int fileNum;
+
+ //fileNum = upper_grid_y / local_grid_size_y;
+ if((fileNum = findFileNum(data_y)) < 0)
+ {
+ cout << "OutCoreInterp::update() findFileNum() error!" << endl;
+ cout << "data_x: " << data_x << " data_y: " << data_y << " grid_y: " << (int)(data_y/GRID_DIST_Y) << " fileNum: " << fileNum << " open file: " << openFile << endl;
+ return -1;
+ }
+
+ if(openFile == fileNum)
+ {
+ // write into memory;
+ updateInterpArray(fileNum, data_x, data_y, data_z);
+
+ }else{
+ UpdateInfo ui(data_x, data_y, data_z);
+ qlist[fileNum].push_back(ui);
+
+ if(qlist[fileNum].size() == QUEUE_LIMIT)
+ {
+ //cout << "erase" << endl;
+ if(openFile != -1)
+ {
+ // if we use mmap and write directly into disk, no need this step.
+ // munmap would be enough for this step..
+
+ // write back to disk
+ gridMap[openFile]->getGridFile()->unmap();
+ openFile = -1;
+ }
+ // upload from disk to memory
+ gridMap[fileNum]->getGridFile()->map();
+ openFile = fileNum;
+
+ // pop every update information
+ list<UpdateInfo>::const_iterator iter;
+
+ for(iter = qlist[openFile].begin(); iter!= qlist[openFile].end(); iter++){
+ updateInterpArray(openFile, (*iter).data_x, (*iter).data_y, (*iter).data_z);
+ }
+ // flush
+ qlist[openFile].erase(qlist[openFile].begin(), qlist[openFile].end());
+ }
+ }
+
+ return 0;
+}
+
+int OutCoreInterp::finish(char *outputName, int outputFormat, unsigned int outputType)
+{
+ int i, j;
+ GridPoint *p;
+ GridFile *gf;
+ int len_y;
+ int offset;
+
+ //struct tms tbuf;
+ clock_t t0, t1;
+
+
+ /*
+ // managing overlap.
+ read adjacent blocks and store ghost cells and update the cells
+
+ merging multiple files
+ */
+ if(openFile != -1){
+ gridMap[openFile]->getGridFile()->unmap();
+ openFile = -1;
+ }
+
+ ////////////////////////////////////////////////////////////
+ // flushing
+ // can be moved inside the communications
+ for(i = 0; i < numFiles; i++)
+ {
+ if(qlist[i].size() != 0)
+ {
+ if((gf = gridMap[i]->getGridFile()) == NULL)
+ {
+ cout << "OutCoreInterp::finish() getGridFile() NULL" << endl;
+ return -1;
+ }
+
+ gf->map();
+ openFile = i;
+
+ // flush UpdateInfo queue
+ // pop every update information
+ list<UpdateInfo>::const_iterator iter;
+
+ for(iter = qlist[i].begin(); iter!= qlist[i].end(); iter++){
+ updateInterpArray(i, (*iter).data_x, (*iter).data_y, (*iter).data_z);
+ }
+ qlist[i].erase(qlist[i].begin(), qlist[i].end());
+
+ gf->unmap();
+ openFile = -1;
+ gf = NULL;
+ }
+ }
+ ////////////////////////////////////////////////////////////
+
+ for(i = 0; i < numFiles - 1; i++)
+ {
+ //////////////////////////////////////////////////////////////
+ // read the upside overlap
+
+ if((gf = gridMap[i]->getGridFile()) == NULL)
+ {
+ cout << "OutCoreInterp::finish() getGridFile() NULL" << endl;
+ return -1;
+ }
+
+ if(gf->map() != -1)
+ {
+ openFile = i;
+ // Sriram's edit to copy over DEM values to overlap also
+ len_y = 2 * (gridMap[i]->getOverlapUpperBound() - gridMap[i]->getUpperBound() - 1);
+
+ if((p = (GridPoint *)malloc(sizeof(GridPoint) * len_y * GRID_SIZE_X)) == NULL)
+ {
+ cout << "OutCoreInterp::finish() malloc error" << endl;
+ return -1;
+ }
+
+ int start = (gridMap[i]->getOverlapUpperBound() - gridMap[i]->getOverlapLowerBound() - len_y) * GRID_SIZE_X;
+ cout << "copy from " << start << " to " << (start + len_y * GRID_SIZE_X) << endl;
+
+ memcpy(p, &(gf->interp[start]), sizeof(GridPoint) * len_y * (GRID_SIZE_X) );
+
+ gf->unmap();
+ gf = NULL;
+ openFile = -1;
+
+ } else {
+ cout << "OutCoreInterp::finish() gf->map() error" << endl;
+ return -1;
+ }
+ //////////////////////////////////////////////////////////////
+
+ //////////////////////////////////////////////////////////////
+ // update (i+1)th component
+ if((gf = gridMap[i+1]->getGridFile()) == NULL)
+ {
+ cout << "OutCoreInterp::finish() getGridFile() NULL" << endl;
+ return -1;
+ }
+
+ if(gf->map() != -1)
+ {
+ offset = 0;
+ openFile = i - 1;
+
+ for(j = 0; j < len_y * GRID_SIZE_X; j++)
+ {
+ if(gf->interp[j + offset].Zmin > p[j].Zmin)
+ gf->interp[j + offset].Zmin = p[j].Zmin;
+
+ if(gf->interp[j + offset].Zmax < p[j].Zmax)
+ gf->interp[j + offset].Zmax = p[j].Zmax;
+
+ gf->interp[j + offset].Zmean += p[j].Zmean;
+ gf->interp[j + offset].count += p[j].count;
+
+ if (p[j].sum != -1) {
+ gf->interp[j + offset].Zidw += p[j].Zidw;
+ gf->interp[j + offset].sum += p[j].sum;
+ } else {
+ gf->interp[j + offset].Zidw = p[j].Zidw;
+ gf->interp[j + offset].sum = p[j].sum;
+ }
+ }
+
+ if(p != NULL){
+ free(p);
+ p = NULL;
+ }
+
+ gf->unmap();
+ gf = NULL;
+ openFile = -1;
+
+ } else {
+ cout << "OutCoreInterp::finish() gf->map() error" << endl;
+ return -1;
+ }
+ //////////////////////////////////////////////////////////////
+ }
+
+ for(i = numFiles -1; i > 0; i--)
+ {
+ //////////////////////////////////////////////////////////////
+ // read the downside overlap
+
+ if((gf = gridMap[i]->getGridFile()) == NULL)
+ {
+ cout << "GridFile is NULL" << endl;
+ return -1;
+ }
+
+ if(gf->map() != -1)
+ {
+ openFile = i;
+ len_y = 2 * (gridMap[i]->getLowerBound() - gridMap[i]->getOverlapLowerBound());
+
+ if((p = (GridPoint *)malloc(sizeof(GridPoint) * len_y * GRID_SIZE_X)) == NULL)
+ {
+ cout << "OutCoreInterp::finish() malloc error" << endl;
+ return -1;
+ }
+
+ memcpy(p, &(gf->interp[0]), len_y * sizeof(GridPoint) * GRID_SIZE_X);
+
+ //finalize();
+
+ gf->unmap();
+ gf = NULL;
+ openFile = -1;