虽然没有写过类似功能，但是看到好的文章，忍不住要转载，以供自己学习：
原文：isolining package for ActionScript 3
A week or so back I wrote about a package I ported/modified to create the Delaunay triangulation in Flash with a few AS3 classes. As I noted there, such atriangulated irregular network (TIN) allows us to interpolateisolines — lines of constant value (aka isarithms, commonly called contours).
So, given a field of points (weather stations, say)…
…with one or more attributes attached (temperature, say)…
…a TIN can be constructed.
With the above TIN, values can be interpolated along each edge between the points of known values (control points). The interpolation is strictly linear (that is, the value 50 would be interpolated halfway along an edge whose control points were valued 48 and 52).
With a given contouring interval (I’m using 4 degrees F here), we can connect some of these interpolated points, creating our contour lines.
With the previous steps stripped away, this creates a passable isoline map.
The lines are rigid, though, and should be smoothed for presentation. I allow two methods for this. You can use the “simple” method, which just uses the built-in graphics method curveTo between the midpoint of each isoline segment (below with the isoline interval decreased to 3 degrees).
The above looks alright, but the curves are not continuous, closed loops can still have hard corners, and the isolines no longer pass through the interpolated points (we have therefore generalized an already-inaccurate interpolation). My compatriot Andy Woodruff, author of the glorious new Cartogrammar blog, offered to write a nice continuous curve method that ensured isolines would still pass through the interpolated values. You can read about the method inhis post. Here she blows:
Bringing it all together, then, and incorporating the only extra feature I wrote (tinting of isolines), here’s a nice finished isoline map of temperature across the U.S.
My new isolining package for Flash/ActionScript3 accomplishes all of the above, requiring only an array of point data with attribute values attached. The above example, was accomplished with the following lines of code (after drawing the U.S. states from a shapefile).
//first, generate the array of triangles (ITriangle objects) from the point data
var triangles:Array = Delaunay.triangulate(points);
Delaunay.drawDelaunay(triangles, points, triClip, false); //comment this out if you don't want to draw the triangulation
//generate an array of isolines (isoline objects)
var isos:Array = IsoUtils.isoline(triangles, points, triClip, 3, 0);
//create color and class arrays for tinting the isolines
var classesArray:Array = new Array(40, 44, 48, 52, 56, 60, 64, 68, 72, 76);
var colorsArray:Array = new Array(0x051CFD, 0x4602FD, 0x6D0EEB, 0x8400FF, 0xC400FF, 0xEA00FF, 0xFF00E2, 0xFF0095, 0xFF0030, 0xFF0015, 0xFB3507);
//then, actually draw them, using a continuous curve
IsoUtils.drawIsolines(isos, triClip, "continuous", colorsArray, classesArray, .5, .95);
Keep in mind: triangulation is just one interpolation method, and is many ways the least technical (and accurate). More accurate interpolation techniques includeinverse-distance andkriging. ***If you’re having trouble, and your isoline interval is not an integer, check out the comment at line 171 ofisoUtils.as. Please fix that, BTW.
I meant to add other features, but since I started work this past week, I’m posting the package as-is, and invite others to modify. On my wishlist:
- hypsometric tinting, or color between the lines, would allow for more effective terrain or temperature mapping
- support for projections and other coordinate conversions in the drawIsolines method. I have packages for converting lat/long to a number of map projections, but currently the drawIsolines method doesn’t have support for passing a point coordinate conversion method.
- an animated demo. This thing’s lightning-fast, so why not?
- something that would be super wicked would be if someone would implementTanaka’s illuminated contours [pdf] method, that thickens/thins and darkens/lightens lines like so…
…creating beautiful relief maps like the one below
If you add anything to the package, feel free to post a link to your revised version in the comments.
虽然没有细看算法，但是看来该大师做法也是：
1）离散区域呈有限个三角单元；
2）绘制每个三角单元内的isoline；
3）遍历所有三角单元；
4）Done，Enjoy。
但是他的三角单元如何离散得到，且看大师另一篇博文：delaunay triangulation in ActionScript 3，讲述了delaunay算法的前世今生，照抄如下：
update: for a cool usage of Delaunay triangulation, see my isolining package for ActionScript 3
The Delaunay triangulation was invented in 1934 by Boris Delaunay. According to Paul Bourke,
The Delaunay triangulation is closely related geometrically to the Direchlet tesselation also known as the Voronoi or Theissen tesselations. These tesselations split the plane into a number of polygonal regions called tiles. Each tile has one sample point in its interior called a generating point. All other points inside the polygonal tile are closer to the generating point than to any other. The Delauney triangulation is created by connecting all generating points which share a common tile edge. Thus formed, the triangle edges are perpendicular bisectors of the tile edges.
（由于不会在CSDN插入flash，该Demo地在此：http://indiemaps.com/flash/clickTest_v1.swf，抄袭者注）
With the idea of creating client-side isolines and inspired in part by the spirit of Keith Peters’ ongoing MathWorld Problem of the Week, I set out to write an AS3 class to create the above triangulation. After an hour of hacking, I gave up, and decided to just port the algorithm from another language. I choseFlorian Jenett’s Java version, itself a port of Paul Bourke’s original C. The port was easy, and makes it a cinch to create a triangulation of, say, selected weather stations in the U.S.
以下是Java版本的该算法：
/*
* ported from p bourke's triangulate.c
* http://astronomy.swin.edu.au/~pbourke/terrain/triangulate/triangulate.c
*
* fjenett, 20th february 2005, offenbach-germany.
* contact: http://www.florianjenett.de/
*
* run like this:
* javac *.java
* java triangulate
*
* to view the output: http://processing.org/
*
*/
class ITRIANGLE {
int p1, p2, p3;
ITRIANGLE() { ; }
}
class IEDGE {
int p1, p2;
IEDGE() { p1=-1; p2=-1; }
}
class XYZ {
double x, y, z;
XYZ() { ; }
XYZ( double _x, double _y, double _z) {
this.x = _x; this.y = _y; this.z = _z;
}
}
public class triangulate {
public static double EPSILON = 0.000001;
/*
Return TRUE if a point (xp,yp) is inside the circumcircle made up
of the points (x1,y1), (x2,y2), (x3,y3)
The circumcircle centre is returned in (xc,yc) and the radius r
NOTE: A point on the edge is inside the circumcircle
*/
static boolean CircumCircle(
double xp, double yp,
double x1, double y1,
double x2, double y2,
double x3, double y3,
/*double xc, double yc, double r*/
XYZ circle
)
{
double m1,m2,mx1,mx2,my1,my2;
double dx,dy,rsqr,drsqr;
double xc, yc, r;
/* Check for coincident points */
if ( Math.abs(y1-y2) < EPSILON && Math.abs(y2-y3) < EPSILON )
{
System.out.println("CircumCircle: Points are coincident.");
return false;
}
if ( Math.abs(y2-y1) < EPSILON )
{
m2 = - (x3-x2) / (y3-y2);
mx2 = (x2 + x3) / 2.0;
my2 = (y2 + y3) / 2.0;
xc = (x2 + x1) / 2.0;
yc = m2 * (xc - mx2) + my2;
}
else if ( Math.abs(y3-y2) < EPSILON )
{
m1 = - (x2-x1) / (y2-y1);
mx1 = (x1 + x2) / 2.0;
my1 = (y1 + y2) / 2.0;
xc = (x3 + x2) / 2.0;
yc = m1 * (xc - mx1) + my1;
}
else
{
m1 = - (x2-x1) / (y2-y1);
m2 = - (x3-x2) / (y3-y2);
mx1 = (x1 + x2) / 2.0;
mx2 = (x2 + x3) / 2.0;
my1 = (y1 + y2) / 2.0;
my2 = (y2 + y3) / 2.0;
xc = (m1 * mx1 - m2 * mx2 + my2 - my1) / (m1 - m2);
yc = m1 * (xc - mx1) + my1;
}
dx = x2 - xc;
dy = y2 - yc;
rsqr = dx*dx + dy*dy;
r = Math.sqrt(rsqr);
dx = xp - xc;
dy = yp - yc;
drsqr = dx*dx + dy*dy;
circle.x = xc;
circle.y = yc;
circle.z = r;
return ( drsqr <= rsqr ? true : false );
}
/*
Triangulation subroutine
Takes as input NV vertices in array pxyz
Returned is a list of ntri triangular faces in the array v
These triangles are arranged in a consistent clockwise order.
The triangle array 'v' should be malloced to 3 * nv
The vertex array pxyz must be big enough to hold 3 more points
The vertex array must be sorted in increasing x values say
qsort(p,nv,sizeof(XYZ),XYZCompare);
int XYZCompare(void *v1,void *v2)
{
XYZ *p1,*p2;
p1 = v1;
p2 = v2;
if (p1->x < p2->x)
return(-1);
else if (p1->x > p2->x)
return(1);
else
return(0);
}
*/
static int Triangulate ( int nv, XYZ pxyz[], ITRIANGLE v[] )
{
boolean complete[] = null;
IEDGE edges[] = null;
int nedge = 0;
int trimax, emax = 200;
int status = 0;
boolean inside;
//int i, j, k;
double xp, yp, x1, y1, x2, y2, x3, y3, xc, yc, r;
double xmin, xmax, ymin, ymax, xmid, ymid;
double dx, dy, dmax;
int ntri = 0;
/* Allocate memory for the completeness list, flag for each triangle */
trimax = 4*nv;
complete = new boolean[trimax];
for (int ic=0; ic<trimax; ic++) complete[ic] = false;
/* Allocate memory for the edge list */
edges = new IEDGE[emax];
for (int ie=0; ie<emax; ie++) edges[ie] = new IEDGE();
/*
Find the maximum and minimum vertex bounds.
This is to allow calculation of the bounding triangle
*/
xmin = pxyz[0].x;
ymin = pxyz[0].y;
xmax = xmin;
ymax = ymin;
for (int i=1;i<nv;i++)
{
if (pxyz[i].x < xmin) xmin = pxyz[i].x;
if (pxyz[i].x > xmax) xmax = pxyz[i].x;
if (pxyz[i].y < ymin) ymin = pxyz[i].y;
if (pxyz[i].y > ymax) ymax = pxyz[i].y;
}
dx = xmax - xmin;
dy = ymax - ymin;
dmax = (dx > dy) ? dx : dy;
xmid = (xmax + xmin) / 2.0;
ymid = (ymax + ymin) / 2.0;
/*
Set up the supertriangle
This is a triangle which encompasses all the sample points.
The supertriangle coordinates are added to the end of the
vertex list. The supertriangle is the first triangle in
the triangle list.
*/
pxyz[nv+0].x = xmid - 2.0 * dmax;
pxyz[nv+0].y = ymid - dmax;
pxyz[nv+0].z = 0.0;
pxyz[nv+1].x = xmid;
pxyz[nv+1].y = ymid + 2.0 * dmax;
pxyz[nv+1].z = 0.0;
pxyz[nv+2].x = xmid + 2.0 * dmax;
pxyz[nv+2].y = ymid - dmax;
pxyz[nv+2].z = 0.0;
v[0].p1 = nv;
v[0].p2 = nv+1;
v[0].p3 = nv+2;
complete[0] = false;
ntri = 1;
/*
Include each point one at a time into the existing mesh
*/
for (int i=0;i<nv;i++) {
xp = pxyz[i].x;
yp = pxyz[i].y;
nedge = 0;
/*
Set up the edge buffer.
If the point (xp,yp) lies inside the circumcircle then the
three edges of that triangle are added to the edge buffer
and that triangle is removed.
*/
XYZ circle = new XYZ();
for (int j=0;j<ntri;j++)
{
if (complete[j])
continue;
x1 = pxyz[v[j].p1].x;
y1 = pxyz[v[j].p1].y;
x2 = pxyz[v[j].p2].x;
y2 = pxyz[v[j].p2].y;
x3 = pxyz[v[j].p3].x;
y3 = pxyz[v[j].p3].y;
inside = CircumCircle( xp, yp, x1, y1, x2, y2, x3, y3, circle );
xc = circle.x; yc = circle.y; r = circle.z;
if (xc + r < xp) complete[j] = true;
if (inside)
{
/* Check that we haven't exceeded the edge list size */
if (nedge+3 >= emax)
{
emax += 100;
IEDGE[] edges_n = new IEDGE[emax];
for (int ie=0; ie<emax; ie++) edges_n[ie] = new IEDGE();
System.arraycopy(edges, 0, edges_n, 0, edges.length);
edges = edges_n;
}
edges[nedge+0].p1 = v[j].p1;
edges[nedge+0].p2 = v[j].p2;
edges[nedge+1].p1 = v[j].p2;
edges[nedge+1].p2 = v[j].p3;
edges[nedge+2].p1 = v[j].p3;
edges[nedge+2].p2 = v[j].p1;
nedge += 3;
v[j].p1 = v[ntri-1].p1;
v[j].p2 = v[ntri-1].p2;
v[j].p3 = v[ntri-1].p3;
complete[j] = complete[ntri-1];
ntri--;
j--;
}
}
/*
Tag multiple edges
Note: if all triangles are specified anticlockwise then all
interior edges are opposite pointing in direction.
*/
for (int j=0;j<nedge-1;j++)
{
//if ( !(edges[j].p1 < 0 && edges[j].p2 < 0) )
for (int k=j+1;k<nedge;k++)
{
if ((edges[j].p1 == edges[k].p2) && (edges[j].p2 == edges[k].p1))
{
edges[j].p1 = -1;
edges[j].p2 = -1;
edges[k].p1 = -1;
edges[k].p2 = -1;
}
/* Shouldn't need the following, see note above */
if ((edges[j].p1 == edges[k].p1) && (edges[j].p2 == edges[k].p2))
{
edges[j].p1 = -1;
edges[j].p2 = -1;
edges[k].p1 = -1;
edges[k].p2 = -1;
}
}
}
/*
Form new triangles for the current point
Skipping over any tagged edges.
All edges are arranged in clockwise order.
*/
for (int j=0;j<nedge;j++)
{
if (edges[j].p1 == -1 || edges[j].p2 == -1)
continue;
if (ntri >= trimax) return -1;
v[ntri].p1 = edges[j].p1;
v[ntri].p2 = edges[j].p2;
v[ntri].p3 = i;
complete[ntri] = false;
ntri++;
}
}
/*
Remove triangles with supertriangle vertices
These are triangles which have a vertex number greater than nv
*/
for (int i=0;i<ntri;i++)
{
if (v[i].p1 >= nv || v[i].p2 >= nv || v[i].p3 >= nv)
{
v[i] = v[ntri-1];
ntri--;
i--;
}
}
return ntri;
}
public static void main (String[] args)
{
int nv = 20;
if (args.length > 0 && args[0] != null) nv = new Integer(args[0]).intValue();
if (nv <= 0 || nv > 1000) nv = 20;
//System.out.println("Creating " + nv + " random points.");
XYZ[] points = new XYZ[ nv+3 ];
for (int i=0; i<points.length; i++)
points[i] = new XYZ( i*4.0, 400.0 * Math.random(), 0.0 );
ITRIANGLE[] triangles = new ITRIANGLE[ nv*3 ];
for (int i=0; i<triangles.length; i++)
triangles[i] = new ITRIANGLE();
int ntri = Triangulate( nv, points, triangles );
/*
copy-paste the following output into free processing:
http://processing.org/
*/
System.out.println("size(400,400); noFill();");
for (int tt=0; tt<points.length; tt++)
{
System.out.println("rect("+points[tt].x+","+points[tt].y+", 3, 3);");
}
for (int tt=0; tt<ntri; tt++)
{
System.out.println("beginShape(TRIANGLES);");
System.out.println("vertex( "+points[triangles[tt].p1].x+","+points[triangles[tt].p1].y+");");
System.out.println("vertex( "+points[triangles[tt].p2].x+","+points[triangles[tt].p2].y+");");
System.out.println("vertex( "+points[triangles[tt].p3].x+","+points[triangles[tt].p3].y+");");
System.out.println("endShape();");
}
}
}
以下是C语言版本的该算法：
typedef struct {
int p1,p2,p3;
} ITRIANGLE;
typedef struct {
int p1,p2;
} IEDGE;
typedef struct {
double x,y,z;
} XYZ;
/*
Triangulation subroutine
Takes as input NV vertices in array pxyz
Returned is a list of ntri triangular faces in the array v
These triangles are arranged in a consistent clockwise order.
The triangle array 'v' should be malloced to 3 * nv
The vertex array pxyz must be big enough to hold 3 more points
The vertex array must be sorted in increasing x values say
qsort(p,nv,sizeof(XYZ),XYZCompare);
:
int XYZCompare(void *v1,void *v2)
{
XYZ *p1,*p2;
p1 = v1;
p2 = v2;
if (p1->x < p2->x)
return(-1);
else if (p1->x > p2->x)
return(1);
else
return(0);
}
*/
int Triangulate(int nv,XYZ *pxyz,ITRIANGLE *v,int *ntri)
{
int *complete = NULL;
IEDGE *edges = NULL;
int nedge = 0;
int trimax,emax = 200;
int status = 0;
int inside;
int i,j,k;
double xp,yp,x1,y1,x2,y2,x3,y3,xc,yc,r;
double xmin,xmax,ymin,ymax,xmid,ymid;
double dx,dy,dmax;
/* Allocate memory for the completeness list, flag for each triangle */
trimax = 4 * nv;
if ((complete = malloc(trimax*sizeof(int))) == NULL) {
status = 1;
goto skip;
}
/* Allocate memory for the edge list */
if ((edges = malloc(emax*(long)sizeof(EDGE))) == NULL) {
status = 2;
goto skip;
}
/*
Find the maximum and minimum vertex bounds.
This is to allow calculation of the bounding triangle
*/
xmin = pxyz[0].x;
ymin = pxyz[0].y;
xmax = xmin;
ymax = ymin;
for (i=1;i<nv;i++) {
if (pxyz[i].x < xmin) xmin = pxyz[i].x;
if (pxyz[i].x > xmax) xmax = pxyz[i].x;
if (pxyz[i].y < ymin) ymin = pxyz[i].y;
if (pxyz[i].y > ymax) ymax = pxyz[i].y;
}
dx = xmax - xmin;
dy = ymax - ymin;
dmax = (dx > dy) ? dx : dy;
xmid = (xmax + xmin) / 2.0;
ymid = (ymax + ymin) / 2.0;
/*
Set up the supertriangle
This is a triangle which encompasses all the sample points.
The supertriangle coordinates are added to the end of the
vertex list. The supertriangle is the first triangle in
the triangle list.
*/
pxyz[nv+0].x = xmid - 20 * dmax;
pxyz[nv+0].y = ymid - dmax;
pxyz[nv+0].z = 0.0;
pxyz[nv+1].x = xmid;
pxyz[nv+1].y = ymid + 20 * dmax;
pxyz[nv+1].z = 0.0;
pxyz[nv+2].x = xmid + 20 * dmax;
pxyz[nv+2].y = ymid - dmax;
pxyz[nv+2].z = 0.0;
v[0].p1 = nv;
v[0].p2 = nv+1;
v[0].p3 = nv+2;
complete[0] = FALSE;
*ntri = 1;
/*
Include each point one at a time into the existing mesh
*/
for (i=0;i<nv;i++) {
xp = pxyz[i].x;
yp = pxyz[i].y;
nedge = 0;
/*
Set up the edge buffer.
If the point (xp,yp) lies inside the circumcircle then the
three edges of that triangle are added to the edge buffer
and that triangle is removed.
*/
for (j=0;j<(*ntri);j++) {
if (complete[j])
continue;
x1 = pxyz[v[j].p1].x;
y1 = pxyz[v[j].p1].y;
x2 = pxyz[v[j].p2].x;
y2 = pxyz[v[j].p2].y;
x3 = pxyz[v[j].p3].x;
y3 = pxyz[v[j].p3].y;
inside = CircumCircle(xp,yp,x1,y1,x2,y2,x3,y3,&xc,&yc,&r);
if (xc < xp && ((xp-xc)*(xp-xc)) > r)
complete[j] = TRUE;
if (inside) {
/* Check that we haven't exceeded the edge list size */
if (nedge+3 >= emax) {
emax += 100;
if ((edges = realloc(edges,emax*(long)sizeof(EDGE))) == NULL) {
status = 3;
goto skip;
}
}
edges[nedge+0].p1 = v[j].p1;
edges[nedge+0].p2 = v[j].p2;
edges[nedge+1].p1 = v[j].p2;
edges[nedge+1].p2 = v[j].p3;
edges[nedge+2].p1 = v[j].p3;
edges[nedge+2].p2 = v[j].p1;
nedge += 3;
v[j] = v[(*ntri)-1];
complete[j] = complete[(*ntri)-1];
(*ntri)--;
j--;
}
}
/*
Tag multiple edges
Note: if all triangles are specified anticlockwise then all
interior edges are opposite pointing in direction.
*/
for (j=0;j<nedge-1;j++) {
for (k=j+1;k<nedge;k++) {
if ((edges[j].p1 == edges[k].p2) && (edges[j].p2 == edges[k].p1)) {
edges[j].p1 = -1;
edges[j].p2 = -1;
edges[k].p1 = -1;
edges[k].p2 = -1;
}
/* Shouldn't need the following, see note above */
if ((edges[j].p1 == edges[k].p1) && (edges[j].p2 == edges[k].p2)) {
edges[j].p1 = -1;
edges[j].p2 = -1;
edges[k].p1 = -1;
edges[k].p2 = -1;
}
}
}
/*
Form new triangles for the current point
Skipping over any tagged edges.
All edges are arranged in clockwise order.
*/
for (j=0;j<nedge;j++) {
if (edges[j].p1 < 0 || edges[j].p2 < 0)
continue;
if ((*ntri) >= trimax) {
status = 4;
goto skip;
}
v[*ntri].p1 = edges[j].p1;
v[*ntri].p2 = edges[j].p2;
v[*ntri].p3 = i;
complete[*ntri] = FALSE;
(*ntri)++;
}
}
/*
Remove triangles with supertriangle vertices
These are triangles which have a vertex number greater than nv
*/
for (i=0;i<(*ntri);i++) {
if (v[i].p1 >= nv || v[i].p2 >= nv || v[i].p3 >= nv) {
v[i] = v[(*ntri)-1];
(*ntri)--;
i--;
}
}
skip:
free(edges);
free(complete);
return(status);
}
/*
Return TRUE if a point (xp,yp) is inside the circumcircle made up
of the points (x1,y1), (x2,y2), (x3,y3)
The circumcircle centre is returned in (xc,yc) and the radius r
NOTE: A point on the edge is inside the circumcircle
*/
int CircumCircle(double xp,double yp,
double x1,double y1,double x2,double y2,double x3,double y3,
double *xc,double *yc,double *rsqr)
{
double m1,m2,mx1,mx2,my1,my2;
double dx,dy,drsqr;
double fabsy1y2 = fabs(y1-y2);
double fabsy2y3 = fabs(y2-y3);
/* Check for coincident points */
if (fabsy1y2 < EPSILON && fabsy2y3 < EPSILON)
return(FALSE);
if (fabsy1y2 < EPSILON) {
m2 = - (x3-x2) / (y3-y2);
mx2 = (x2 + x3) / 2.0;
my2 = (y2 + y3) / 2.0;
*xc = (x2 + x1) / 2.0;
*yc = m2 * (*xc - mx2) + my2;
} else if (fabsy2y3 < EPSILON) {
m1 = - (x2-x1) / (y2-y1);
mx1 = (x1 + x2) / 2.0;
my1 = (y1 + y2) / 2.0;
*xc = (x3 + x2) / 2.0;
*yc = m1 * (*xc - mx1) + my1;
} else {
m1 = - (x2-x1) / (y2-y1);
m2 = - (x3-x2) / (y3-y2);
mx1 = (x1 + x2) / 2.0;
mx2 = (x2 + x3) / 2.0;
my1 = (y1 + y2) / 2.0;
my2 = (y2 + y3) / 2.0;
*xc = (m1 * mx1 - m2 * mx2 + my2 - my1) / (m1 - m2);
if (fabsy1y2 > fabsy2y3) {
*yc = m1 * (*xc - mx1) + my1;
} else {
*yc = m2 * (*xc - mx2) + my2;
}
}
dx = x2 - *xc;
dy = y2 - *yc;
*rsqr = dx*dx + dy*dy;
dx = xp - *xc;
dy = yp - *yc;
drsqr = dx*dx + dy*dy;
// Original
//return((drsqr <= *rsqr) ? TRUE : FALSE);
// Proposed by Chuck Morris
return((drsqr - *rsqr) <= EPSILON ? TRUE : FALSE);
}
作为一个喜欢画蛇添足的人，最后给出自己一段类似的算法，程序是根据有限元中的节点MPoint（FEModel.model.vo.primitives.MPoint，不好意思，当时为了时髦，用了pureMVC架构，名称太长，事实证明，不需要架构时用架构可能入不敷出）和单元（FEModel.model.vo.primitives.MElement）拓扑关系，亦即网格，如下图所示，根据网格绘制Contourr图形，我们的网格是根据图形而剖分，非前述由离散点根据Delaunay算法得到网格。由于年代久远（大约2010年9月才完成，耗时1.5年多才完成)，后处理的代码也许是个半成品，代码暂时不做解释了，况且笔者也不能很快读懂，所以养成一个好的写代码做注释的习惯是很重要的。好吧，开始给貂续狗狗尾吧了，以自己的程序结束本文：
package FEModel.util.visual.contour
{
import FEModel.model.vo.primitives.MElement;
import FEModel.model.vo.primitives.MPoint;
import FEModel.util.math.algebra.NumbericUtil;
import FEModel.util.math.geometry.PointUtil;
import flash.display.Graphics;
public class Contour2DUtil
{
public static function draw2DTContour(g:Graphics,
mpV:Vector.<MPoint>,
tsFun:Function,
colorKey:Vector.<uint>,
valueKey:Vector.<Number>,
mode:uint=0):void
{
var pLstA:Vector.<MPoint>=PointUtil.lineContainValues(mpV[0],mpV[2],valueKey);
var pLstB:Vector.<MPoint>=PointUtil.lineContainValues(mpV[0],mpV[1],valueKey);
for each(var tmpMPoint:MPoint in pLstA)
{
if(tmpMPoint.T==mpV[1].T)
{
pLstB.push(mpV[1]);
break;
}
}
var pLstC:Vector.<MPoint>=PointUtil.lineContainValues(mpV[1],mpV[2],valueKey);
if(pLstB!=null)
{
if((pLstC!=null)&&(pLstC.length))
{
for each(var mPoint:MPoint in pLstC)
pLstB.push(mPoint);
}
}
else
{
pLstB=pLstC;
}
pLstC=null;
var triMode:uint=triangleMode(mpV,pLstA);
if(triMode==0)
{
pLstA.push(mpV[0]);
pLstB.push(mpV[0]);
pLstA.push(mpV[2]);
pLstB.push(mpV[1]);
}
else if(triMode==1)
{
pLstA.unshift(mpV[0]);
pLstB.unshift(mpV[0]);
pLstA.push(mpV[2]);
pLstB.push(mpV[1]);
}
else if(triMode==2)
{
pLstA.unshift(mpV[0]);
pLstB.unshift(mpV[1]);
pLstA.push(mpV[2]);
pLstB.push(mpV[2]);
}
else if(triMode==3)
{
pLstA.unshift(mpV[0]);
pLstB.unshift(mpV[0]);
pLstA.push(mpV[2]);
pLstB.push(mpV[2]);
}
if(pLstA.length!=pLstB.length)
return;
var pos:MPoint;
var path:Vector.<Number>;
var pathCmd:Vector.<int>=new Vector.<int>;
pathCmd.push(1,2,2,2,2);
for(var i:uint=0;i<pLstA.length-1;i++)
{
path=new Vector.<Number>;
pos=tsFun(pLstA[i]);
path.push(pos.x,pos.y);
pos=tsFun(pLstB[i]);
path.push(pos.x,pos.y);
if(triMode==3)
{
if(NumbericUtil.NumberContain(mpV[1].T,pLstB[i].T,pLstB[i+1].T,false))
{
pos=tsFun(mpV[1]);
path.push(pos.x,pos.y);
pathCmd.push(2);
}
}
pos=tsFun(pLstB[i+1]);
path.push(pos.x,pos.y);
pos=tsFun(pLstA[i+1]);
path.push(pos.x,pos.y);
pos=tsFun(pLstA[i]);
path.push(pos.x,pos.y);
var color:uint=GetCColor((pLstA[i].T+pLstA[i+1].T)/2,colorKey,valueKey);
g.beginFill(color);
g.drawPath(pathCmd,path);
g.endFill();
}
}
public static function draw2DContour(g:Graphics,
mpV:Vector.<MPoint>,
meV:Vector.<MElement>,
tsFun:Function,
colorKey:Vector.<uint>,
valueKey:Vector.<Number>,
mode:uint=0):void
{
g.clear();
for each(var me:MElement in meV)
{
var triangle:Vector.<MPoint>=new Vector.<MPoint>(3,true);
for(var i:uint=0;i<3;i++)
{
triangle[i]=mpV[me.p[i]];
}
triangle.sort(MPoint.sortMPoint);
draw2DTContour(g,triangle,tsFun,colorKey,valueKey,mode);
triangle=null;
}
}
public static function triangleMode(mpV:Vector.<MPoint>,contourLine:Vector.<MPoint>):uint
{
if(contourLine.length==0)
return 0;
var hiCnt:uint=0,loCnt:uint=0;
var lo:Number=contourLine[0].T;
var hi:Number=0;
if(contourLine.length==1)
{
for each(var mp0:MPoint in mpV)
{
if(mp0.T<lo)
loCnt++;
}
return loCnt;
}
hi=contourLine[contourLine.length-1].T;
for each(var mp:MPoint in mpV)
{
if(mp.T>hi) hiCnt++;
if(mp.T<lo) loCnt++;
}
if(loCnt==hiCnt)
return 3;
else
return loCnt;
}
public static function GetCColor(value:Number,colorKey:Vector.<uint>,valueKey:Vector.<Number>):uint
{
if(value<=valueKey[0])
return colorKey[0];
else if(value>=valueKey[valueKey.length-1])
return colorKey[colorKey.length-1];
for(var i:uint=0;i<valueKey.length-1;i++)
{
if(NumbericUtil.NumberContain(value,valueKey[i],valueKey[i+1]))
return colorKey[i+1];
}
return 0;
}
}
}另外，本文所转载的isoline AS3 API的作者Zachary Forest Johnson可能是个艺术青年，把直线做了样条曲线的艺术处理，可能会造成误差，普通青年注意了，2青年随意。
2736
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
