Original CONREC code by Paul Bourke.
This PHP implementation by Ashley Norris.
X Axis: -9.425 to 9.425, Y Axis: -9.425 to 9.425
X Axis: -3.142 to 9.425, Y Axis: -6.283 to 6.283
X Axis: -5 to 10, Y Axis: -6.283 to 6.283
X Axis: 0 to 10, Y Axis: 0 to 10
End of CONREC demonstrations.
The library code is shown below.
<?php
// ------------------------------------------------------------------
// ------------------------------------------------------------------
// CONREC is a straightforward method of contouring some surface represented
// as a regular triangular mesh.
//
// Ported from C / Fortran code by Paul Borke.
// http://paulbourke.net/papers/conrec/
//
// PHP implementation by Ashley Norris.
// http://norris.org.au/
//
// Contouring aids in visualizing three dimensional surfaces on a two dimensional
// medium (on paper or in this case a computer graphics screen). Two most common
// applications are displaying topological features of an area on a map or the air
// pressure on a weather map. In all cases some parameter is plotted as a function
// of two variables, the longitude and latitude or x and y axis. One problem with
// computer contouring is the process is usually CPU intensive and the algorithms
// often use advanced mathematical techniques making them susceptible to error.
//
// ------------------------------------------------------------------
// ------------------------------------------------------------------
//
// This PHP implementation accepts a 2D scalar field (a matrix) and returns
// a set of line segments that would typically then be rendered in some way by the user.
//
// The form of the returned data has been made as simple as possible, and could
// be easily edited (to minimise memory use) by anyone if they so wish.
//
// The returned data is an array of this structure:
// [top] (array of contours) // An array containing all of the contours created by this algorithm
// -> [contour] (mixed array)
// -> 'value' (double) // The scalar value (height) of the contour
// -> ['segments'] (array of segments) // An array containing all the segments (not necessarily continuous) that make up the contour
// -> 'x1' (double)
// -> 'y1' (double)
// -> 'x2' (double)
// -> 'y2' (double)
//
// Best thing to do is to use the function and then do a "print_r();" on its returned value to explore the structure returned.
//
// ------------------------------------------------------------------
// ------------------------------------------------------------------
//
// This library file produces no output, and contains a definition for one useful function: contour();
//
// There are some other functions defined at the end of the file for internal use, which would be useless for external calls.
//
// ------------------------------------------------------------------
// ------------------------------------------------------------------
// Function contour($d, $x, $y, $z, $SKIP_QC = false);
//
// === Input Parameters ===
//
// $d => 2D matrix of values, typically of type double
// $x => 1D array of the x-coordinate for each column
// $y => 1D array of the y-coordinate for each row
// $z => 1D array of the Z value to use for each contour **OR** an integer number of contours to automatically create
// $SKIP_QC => boolean, set to true to skip the Quality Control stage (speed benefit for large datasets)
//
// Some notes:
// - $d is a 2D jagged array, with X on the first dim and Y on the second.
// I say jagged as PHP does not support a true native multidimensional array
// type in the same way as other languages do. Thus, (C, Fortran, etc...) d[i,j] === (PHP) $d[i][j]
// - $x must be the same length as the width (first dim) of $d
// - $y must be the same length as the height (second dim) of $d
// - The $x and $y arrays are not specifically "needed" to create contours, but are
// used to return contour line segments in that coordinate system, rather than a unit system of col and row number.
// - $z takes advantage of PHP's flexible typing (yes, it really is good for some things!)
// Either, pass in an array of values and the number of contours returned will equal the length of $z
// Or, pass in an integer, and that many contours will be returned.
// ------------------------------------------------------------------
// ------------------------------------------------------------------
function CONREC_contour($d, $x, $y, $z, $SKIP_QC = false)
{
// ==================================================
// QC the inputs (turn off for large datasets)
// If wanting to improve performance, but maintain QC, then feel free to weave the QC checks into the actual algorithm code below
// ==================================================
if (!$SKIP_QC) {
// Basic checks first
if (!is_array($d)) { die('Function contour(); failed QC. Parameter $d is not an array.'); }
if (!is_array($x)) { die('Function contour(); failed QC. Parameter $x is not an array.'); }
if (!is_array($y)) { die('Function contour(); failed QC. Parameter $y is not an array.'); }
// NOTE: the count() function in PHP is NOT recursive by default
if (count($d) <= 1) { die('Function contour(); failed QC. Dimension 1 of $d has length smaller than or equal to 1.'); }
if (count($x) <= 1) { die('Function contour(); failed QC. Array $x has length smaller than or equal to 1.'); }
if (count($y) <= 1) { die('Function contour(); failed QC. Array $y has length smaller than or equal to 1.'); }
// Check first dim of $d against $x
if (count($d) != count($x)) { die('Function contour(); failed QC. Dimension 1 of $d ['.count($d).'] does not match length of $x ['.count($x).'].'); }
// Loop first dim of $d and check second dim for all rows (as a side-effect this ensures $d is regular)
foreach($d as $index=>$row) {
if (!is_int($index)) { die('Function contour(); failed QC. Row "'.$index.'" of $d must have an integer index.'); }
if (count($row) != count($y)) { die('Function contour(); failed QC. Dimension 2 of row "'.$index.'" of $d does not match length of $y.'); }
// Check indexes of row
foreach(array_keys($row) as $col_index) {
if (!is_int($col_index)) { die('Function contour(); failed QC. On row "'.$index.'", col "'.$col_index.'" of $d must have an integer index.'); }
}
}
// $z is either an array, or an integer
if (!is_array($z) AND !is_int($z)) { die('Function contour(); failed QC. Parameter $z must be either an array or an integer.'); }
}
// ==================================================
// ==================================================
// While PHP does not require us to declare our variables, we do so here anyway, as per the original C code
$x1 = 0.0;
$x2 = 0.0;
$y1 = 0.0;
$y2 = 0.0;
// ==================================================
// For these arrays we can not specify a fixed length, we just init them with zeros and start using them
// ==================================================
$h = array(0,0,0,0,0); // Length 5
$sh = array(0,0,0,0,0); // Length 5
$xh = array(0,0,0,0,0); // Length 5
$yh = array(0,0,0,0,0); // Length 5
// ==================================================
// When getting the range of the X and Y values we can use the first row of the matrix thanks to the QC step above
// ==================================================
// NOTE: the count() function in PHP is NOT recursive by default
$ilb = 0; // Can assume this thanks to QC step
$iub = count($d) - 1; // Counts first dim
$jlb = 0; // Can assume this thanks to QC step
$jub = count($d[0]) - 1; // Counts second dim of first row, which we can do thanks to QC step
if (is_array($z)) {
$nc = count($z);
}
else if (is_int($z)) {
$nc = $z;
}
else {
// Sensible hard-coded default value, in case QC was skipped
$nc = 10;
}
// ==================================================
// Automatic contours (in case only the number has been defined)
// ==================================================
if (!is_array($z)) {
$z = CONREC_auto_z($nc,$d);
}
// ==================================================
// Sort $z
sort($z);
// ==================================================
// The indexing of im and jm should be noted as it has to start from zero
// unlike the fortran counter part
$im = array(0, 1, 1, 0 );
$jm = array( 0, 0, 1, 1 );
// ==================================================
// Note that castab is arranged differently from the FORTRAN code because
// Fortran and C/C++ arrays are transposed of each other, in this case
// it is more tricky as castab is in 3 dimensions
$castab = array(
array( array( 0, 0, 8 ), array( 0, 2, 5 ), array( 7, 6, 9 ) ),
array( array( 0, 3, 4 ), array( 1, 3, 1 ), array( 4, 3, 0 ) ),
array( array( 9, 6, 7 ), array( 5, 2, 0 ), array( 8, 0, 0 ) )
);
// ==================================================
// We are now ready to perform the contouring, init the return array
$return_array = array();
foreach($z as $index=>$value) {
$return_array[$index] = array('value'=>$value,'segments'=>array());
}
// ==================================================
// Primary loop: i (row, or $d first dimension) loop begins here
for ($i = $ilb; $i <= $iub - 1; $i++)
{
// ==================================================
// Primary loop: j (col, or $d second dimension) loop begins here
for ($j = $jub - 1; $j >= $jlb; $j--)
{
// ==================================================
// Find the max and min value of the corners
$corners = array($d[$i][$j], $d[$i][$j + 1], $d[$i + 1][$j], $d[$i + 1][$j + 1]);
$dmin = min($corners);
$dmax = max($corners);
// Perform elimination of trivial cases
if ($dmax >= $z[0] && $dmin <= $z[$nc - 1])
{
// Start k (contour level) loop
for ($k = 0; $k < $nc; $k++)
{
// Test against this contour level
if ($z[$k] >= $dmin && $z[$k] <= $dmax)
{
// Loop over $m (corners)
for ($m = 4; $m >= 0; $m--)
{
if ($m > 0)
{
// The indexing of im and jm should be noted as it has to
// start from zero
$h[$m] = $d[$i + $im[$m - 1]][ $j + $jm[$m - 1]] - $z[$k];
$xh[$m] = $x[$i + $im[$m - 1]];
$yh[$m] = $y[$j + $jm[$m - 1]];
}
else
{
// Compute average
$h[0] = 0.25 * ($h[1] + $h[2] + $h[3] + $h[4]);
$xh[0] = 0.5 * ($x[$i] + $x[$i + 1]);
$yh[0] = 0.5 * ($y[$j] + $y[$j + 1]);
}
if ($h[$m] > 0.0)
{
$sh[$m] = 1;
}
else if ($h[$m] < 0.0)
{
$sh[$m] = -1;
}
else
{
$sh[$m] = 0;
}
} // end $m loop
// Note: at this stage the relative heights of the corners and the
// centre are in the h array, and the corresponding coordinates are
// in the xh and yh arrays. The centre of the box is indexed by 0
// and the 4 corners by 1 to 4 as shown below.
// Each triangle is then indexed by the parameter m, and the 3
// vertices of each triangle are indexed by parameters m1,m2,and
// m3.
// It is assumed that the centre of the box is always vertex 2
// though this is important only when all 3 vertices lie exactly on
// the same contour level, in which case only the side of the box
// is drawn.
//
// vertex 4
// +-------------------+ vertex 3
// | \ / |
// | \ m-3 / |
// | \ / |
// | \ / |
// | m=2 X m=2 | the centre is vertex 0
// | / \ |
// | / \ |
// | / m=1 \ |
// | / \ |
// +-------------------+ vertex 2
// vertex 1
//
// Scan each triangle in the box
for ($m = 1; $m <= 4; $m++)
{
$m1 = $m;
$m2 = 0;
$m3 = 0;
if ($m != 4)
{
$m3 = $m + 1;
}
else
{
$m3 = 1;
}
$caseValue = $castab[$sh[$m1] + 1][$sh[$m2] + 1][$sh[$m3] + 1];
if ($caseValue != 0)
{
switch ($caseValue)
{
case 1: // Line between vertices 1 and 2
$x1 = $xh[$m1];
$y1 = $yh[$m1];
$x2 = $xh[$m2];
$y2 = $yh[$m2];
break;
case 2: // Line between vertices 2 and 3
$x1 = $xh[$m2];
$y1 = $yh[$m2];
$x2 = $xh[$m3];
$y2 = $yh[$m3];
break;
case 3: // Line between vertices 3 and 1
$x1 = $xh[$m3];
$y1 = $yh[$m3];
$x2 = $xh[$m1];
$y2 = $yh[$m1];
break;
case 4: // Line between vertex 1 and side 2-3
$x1 = $xh[$m1];
$y1 = $yh[$m1];
$x2 = CONREC_sect($m2, $m3, $h, $xh);
$y2 = CONREC_sect($m2, $m3);
break;
case 5: // Line between vertex 2 and side 3-1
$x1 = $xh[$m2];
$y1 = $yh[$m2];
$x2 = CONREC_sect($m3, $m1, $h, $xh);
$y2 = CONREC_sect($m3, $m1, $h, $yh);
break;
case 6: // Line between vertex 3 and side 1-2
$x1 = $xh[$m3];
$y1 = $yh[$m3];
$x2 = CONREC_sect($m1, $m2, $h, $xh);
$y2 = CONREC_sect($m1, $m2, $h, $yh);
break;
case 7: // Line between sides 1-2 and 2-3
$x1 = CONREC_sect($m1, $m2, $h, $xh);
$y1 = CONREC_sect($m1, $m2, $h, $yh);
$x2 = CONREC_sect($m2, $m3, $h, $xh);
$y2 = CONREC_sect($m2, $m3, $h, $yh);
break;
case 8: // Line between sides 2-3 and 3-1
$x1 = CONREC_sect($m2, $m3, $h, $xh);
$y1 = CONREC_sect($m2, $m3, $h, $yh);
$x2 = CONREC_sect($m3, $m1, $h, $xh);
$y2 = CONREC_sect($m3, $m1, $h, $yh);
break;
case 9: // Line between sides 3-1 and 1-2
$x1 = CONREC_sect($m3, $m1, $h, $xh);
$y1 = CONREC_sect($m3, $m1, $h, $yh);
$x2 = CONREC_sect($m1, $m2, $h, $xh);
$y2 = CONREC_sect($m1, $m2, $h, $yh);
break;
default:
break;
}
// ==================================================
// OUTPUT to return array structure
// ==================================================
$return_array[$k]['segments'][] = array('x1'=>$x1,'y1'=>$y1,'x2'=>$x2,'y2'=>$y2);
// NOTE: If wishing to improve performance, or reduce memory use, just put your
// rendering routine right here.
// ==================================================
} // end if: caseValue != 0
} // end for: m loop
}
} // end k loop
} // end if: trivial case elimination (outside contouring range)
} // end for: primary j loop
} // end for: primary i loop
// Return the computed contours
return $return_array;
} // end function: conrec();
// ----------------------------------------------------------------
// Internal functions
// ----------------------------------------------------------------
// As we don't have access to the $xh and $yh parameters we pass them into this single function to perform the segment computation
function CONREC_sect($p1_index, $p2_index, $h, $dimh) {
return ($h[$p2_index] * $dimh[$p1_index] - $h[$p1_index] * $dimh[$p2_index]) / ($h[$p2_index] - $h[$p1_index]);
}
function CONREC_auto_z($nc, $array) {
// ==================================================
// NOTE: Simply dividing the range by the number of desired contours doesn't
// actually work very well, as the max and min values do not produce
// any contours at all. To correct for this we add 2 to the desired number of
// contours, divide the range by that number and then use that step to return the
// actual contour values to use.
// ==================================================
$z = array();
$d_min = CONREC_arraymin($array);
$d_max = CONREC_arraymax($array);
$d_range = $d_max - $d_min;
$contour_step = $d_range / ($nc + 2 - 1); // See above for an explanation of this algebra :)
for($tmp_i = 1;$tmp_i <= $nc;$tmp_i++) {
$z[] = $d_min + ($contour_step * $tmp_i);
}
// $z should now contain $nc elements
return $z;
}
function CONREC_arraymax($array){
if (is_array($array)){
foreach($array as $key => $value) {
$array[$key] = CONREC_arraymax($value);
}
return max($array);
} else{
return $array;
}
}
function CONREC_arraymin($array){
if (is_array($array)) {
foreach($array as $key => $value) {
$array[$key] = CONREC_arraymin($value);
}
return min($array);
} else{
return $array;
}
}
// ----------------------------------------------------------------
// End of file.
// ----------------------------------------------------------------
?>
The test generation code (this file) is shown below.
<!DOCTYPE HTML>
<html>
<head><title>CONREC PHP Demonstration Page - Ashley Norris, 2013</title>
<body>
<h1>CONREC PHP Demonstration Page</h1>
<p>Original <a href="http://paulbourke.net/papers/conrec/">CONREC</a> code by <a href="http://paulbourke.net/">Paul Bourke</a>.</p>
<p>This <a href="http://www.php.net/">PHP</a> implementation by <a href="http://norris.org.au/">Ashley Norris</a>.</p>
<?php
// ------------------------------------------------
// ------------------------------------------------
// Include the CONREC library - this assumes the two scripts are side by side on the server.
// ------------------------------------------------
include('conrec.inc.php');
// ------------------------------------------------
// ------------------------------------------------
// Show a couple of examples
// ------------------------------------------------
// Use of the 'eval' function here is for convenience - general use can not be recommended.
// ------------------------------------------------
$examples = array();
// ------------------------------------------------
$examples[] = array('nc' => 10,
'x_min' => -3 * pi(),
'x_range' => 6 * pi(),
'y_min' => -3 * pi(),
'y_range' => 6 * pi(),
'string' => 'f = X.(sin(X) + sin(Y))',
'eval' => '$val = $x_val * (sin($x_val) + sin($y_val));');
// ------------------------------------------------
$examples[] = array('nc' => 10,
'x_min' => -1 * pi(),
'x_range' => 4 * pi(),
'y_min' => -2 * pi(),
'y_range' => 4 * pi(),
'string' => 'f = sin((X^2 + Y^2)^0.5)',
'eval' => '$val = sin(pow($x_val*$x_val + $y_val*$y_val,0.5));');
// ------------------------------------------------
$examples[] = array('nc' => 15,
'x_min' => -5,
'x_range' => 15,
'y_min' => -2 * pi(),
'y_range' => 4 * pi(),
'string' => 'f = ContouringPaper[Figure18.6]',
'eval' => '$val = sin(pow($x_val*$x_val + $y_val*$y_val,0.5)) + (1 / (pow(pow($x_val - 5,2) + $y_val*$y_val,0.5)));');
// ------------------------------------------------
$examples[] = array('nc' => 6,
'x_min' => 0,
'x_range' => 10,
'y_min' => 0,
'y_range' => 10,
'string' => 'f = random(0,100)',
'eval' => '$val = rand(0,100);');
// ------------------------------------------------
// ------------------------------------------------
// ------------------------------------------------
$N = 100;
$img_size = 900; // The pixel size of the outputted image: be aware of reduced contour quality due to oversampling, keep $N <= ($img_size / 3)
// ------------------------------------------------
// ------------------------------------------------
foreach($examples as $ex) {
$d = array();
$x = array();
$y = array();
$x_min = $ex['x_min'];
$y_min = $ex['y_min'];
$x_range = $ex['x_range'];
$y_range = $ex['y_range'];
// ------------------------------------------------
// ------------------------------------------------
// Create the sampled data based on the function chosen
for($i = 0;$i < $N;$i++) {
for($j = 0;$j < $N ;$j++) {
// x
$x_val = $x_min + ($x_range * $i / ($N - 1));
// y
$y_val = $y_min + ($y_range * $j / ($N - 1));
// val (from eval)
$val = 0;
eval($ex['eval']);
$error = error_get_last();
print_r($error);
// Null or NaN check (crude!)
if (!$val) { $val = 0; }
$x[$i] = $x_val;
$y[$j] = $y_val;
$d[$i][$j] = $val;
}
}
// ------------------------------------------------
// ------------------------------------------------
// ------------------------------------------------
// INVOKE CONREC
// ------------------------------------------------
$contours = CONREC_contour($d,$x,$y,$ex['nc']);
// ------------------------------------------------
// ------------------------------------------------
// ------------------------------------------------
// ------------------------------------------------
// Create output image using GD
// ------------------------------------------------
if (!function_exists('imagecreatetruecolor')) { die('Fatal error: no GD module support. Please enable this in php.ini to use this script.'); }
$img = imagecreatetruecolor($img_size,$img_size);
$white = imagecolorallocate($img,255,255,255);
$black = imagecolorallocate($img,0,0,0);
imagefilledrectangle($img,0,0,$img_size,$img_size,$white);
// ------------------------------------------------
$num_contours = count($contours);
$R_base = rand(0,255);
$G_base = rand(0,255);
$B_base = rand(0,255);
foreach($contours as $index=>$c) {
// make a colour
$R = $R_base + ((255 - $R_base) * ($index / $num_contours));
$G = $G_base * ($index / $num_contours);
$B = $B_base - ($B_base * ($index / $num_contours));
$this_contour = imagecolorallocate($img,$R,$G,$B);
// ---
foreach($c['segments'] as $sgm) {
// Convert X and Y to pixels
$x1 = (($sgm['x1'] - $x_min) / $x_range) * $img_size;
$x2 = (($sgm['x2'] - $x_min) / $x_range) * $img_size;
$y1 = $img_size - (($sgm['y1'] - $y_min) / $y_range) * $img_size;
$y2 = $img_size - (($sgm['y2'] - $y_min) / $y_range) * $img_size;
imageline($img,$x1,$y1,$x2,$y2,$this_contour);
}
// ---
imagefilledrectangle($img,0,($index + 1) * 20,20,($index+2)*20,$this_contour);
imagerectangle($img,0,($index + 1)*20,20,($index + 2)*20,$black);
imagestring($img,2,25,(($index + 1) * 20)+3,round($c['value'],3),$black);
}
// ------------------------------------------------
// Label the function
imagefilledrectangle($img,0,0,300,19,$white);
imagestring($img,2,3,3,$ex['string'],$black);
// ------------------------------------------------
// ------------------------------------------------
// Output PNG file to browser
// ------------------------------------------------
ob_start();
imagepng($img);
$image_data = ob_get_contents();
ob_end_clean();
imagedestroy($img);
// Base64 encode (for embedded URI scheme)
$image_data = base64_encode($image_data);
print('<hr />');
print('<h1>'.$ex['string'].'</h1>');
print('<img style="border:1px solid #555;text-align:padding:1em;" src="data:image/png;base64,'.$image_data.'" />');
print('<p>X Axis: '.round($x_min,3).' to '.round($x_min + $x_range,3).', Y Axis: '.round($y_min,3).' to '.round($y_min + $y_range,3).'</p>');
// ------------------------------------------------
} // end of major example loop
?>
<hr />
<p><i>End of CONREC demonstrations.</i></P>
<hr />
<hr />
<h1>CONREC PHP File</h1>
<p>The library code is shown below.</p>
<hr />
<pre>
<?php
$file = str_replace(array('<','>'),array('<','>'),file_get_contents('conrec.inc.php'));
print($file);
?>
</pre>
<hr />
<h1>Test PHP File (This File)</h1>
<p>The test generation code (this file) is shown below.</p>
<hr />
<pre>
<?php
$file = str_replace(array('<','>'),array('<','>'),file_get_contents('index.php'));
print($file);
?>
</pre>
<hr />
<p><i>End of document.</i></p>
</html>
End of document.