Showing posts with label Algorithm. Show all posts
Showing posts with label Algorithm. Show all posts

Mar 26, 2014

Geospatial function: Point in Ploygon Test

There is an excellent post on point in polygon test from Mathematica Stackexchange.

Point in polygon test is one of the most useful functions when processing Geospatial data. Unfortunately, there is no official built-in function from Mathematica yet. Here is an example of showing usefulness of this function.
We have some data on Greenland icesheet thickness in format of {lon,lat, thickness}.


Our goal is to make a map to show the thickness of icesheet.
First try with ListDensityPlot:


You can see it is not working well, we need to limit the plot region inside the boundary of Greenland. Here is the place we can use the point-in-polygon test. In this example, I use inPolyQ2 from the answer by Simon Woods:


Still not right? What's wrong? The trick is to increase MaxPlotPoints  to 100 at least:
ListDensityPlot[data, PlotRange -> All,
 ColorFunction -> (ColorData["Rainbow"][1 - #] &),
 RegionFunction -> (inPolyQ2[greenland_ploygon[[1, 1]], #1, #2] &),
 MaxPlotPoints -> 150, MeshFunctions -> {#3 &}, Mesh -> 10]


It looks like a real map now.

Oct 5, 2008

QHull for Mathematica

mPower.m is a Mathematica package that interfaces with qhull binaries.

Installation guide is not for windows platform.

Here is the quick way to get it running on Windows Vista with Mathematica 6.0.

1. download mPower.m, and put it into the Applications subdirectory of $UserBaseDirectory 

In[1]:= $UserBaseDirectory

Out[1]= C:\Users\...\AppData\Roaming\Mathematica

2. download qhull for windows, put it into the folder C:\qhull

3. open mPower.m, modify the following two lines,

$QHULL=ToFileName[{$UserBaseDirectory, "Applications", "qhull", "src"}]

to: $QHULL="C:\\qhull\\"

qhullFiles={"qconvex",  "qdelaunay", "qhalf","qhull", "qvoronoi"};

to: qhullFiles={"qconvex.exe",  "qdelaunay.exe", "qhalf.exe","qhull.exe", \
"qvoronoi.exe"};

4. then testing, ignore the warning on regtet and pwrvtx, we only need qhull

<< "mPower.m"

points = RandomReal[1, {10, 3}];

ch = convexHull[points, convexHullFormat -> {facetNormals -> True, facets -> True}]

 

Thanks to the mPower developers

May 3, 2008

Self random walk on sphere

Self random on sphere lattice (longitude, latitude).
(red -- starting point, green -- ending point)


Draw on a sphere:Self avoid version:

The longest horizontal lines indicate the jumps between -180 and 180.

Apr 7, 2008

Minimum Bounding Rectangle of a Point Set

Give a 2-D point set, we need to find the minimum bounding rectangle in term of area.

I read it somewhere that the minimum rectangle must always contain at least one edge of the convex hull, however, I can't find a citation now. So the algorithm can be constructed in the following way: first construct the convex hull, then rotate each edges for the convex hull to the position parallel to x-axes, then calculate area of bounding rectangle; find the minimum one of these rotated rectangles and rotate it back to normal.