| <?php | <?php |
| /*======================================================================= | /*======================================================================= |
| // File: JPGRAPH_REGSTAT.PHP | // File: JPGRAPH_REGSTAT.PHP |
| // Description: Regression and statistical analysis helper classes | // Description: Regression and statistical analysis helper classes |
| // Created: 2002-12-01 | // Created: 2002-12-01 |
| // Ver: $Id: jpgraph_regstat.php 1131 2009-03-11 20:08:24Z ljp $ | // Ver: $Id: jpgraph_regstat.php 1131 2009-03-11 20:08:24Z ljp $ |
| // | // |
| // Copyright (c) Aditus Consulting. All rights reserved. | // Copyright (c) Aditus Consulting. All rights reserved. |
| //======================================================================== | //======================================================================== |
| */ | */ |
| //------------------------------------------------------------------------ | //------------------------------------------------------------------------ |
| // CLASS Spline | // CLASS Spline |
| // Create a new data array from an existing data array but with more points. | // Create a new data array from an existing data array but with more points. |
| // The new points are interpolated using a cubic spline algorithm | // The new points are interpolated using a cubic spline algorithm |
| //------------------------------------------------------------------------ | //------------------------------------------------------------------------ |
| class Spline { | class Spline { |
| // 3:rd degree polynom approximation | // 3:rd degree polynom approximation |
| private $xdata,$ydata; // Data vectors | private $xdata,$ydata; // Data vectors |
| private $y2; // 2:nd derivate of ydata | private $y2; // 2:nd derivate of ydata |
| private $n=0; | private $n=0; |
| function __construct($xdata,$ydata) { | function __construct($xdata,$ydata) { |
| $this->y2 = array(); | $this->y2 = array(); |
| $this->xdata = $xdata; | $this->xdata = $xdata; |
| $this->ydata = $ydata; | $this->ydata = $ydata; |
| $n = count($ydata); | $n = count($ydata); |
| $this->n = $n; | $this->n = $n; |
| if( $this->n !== count($xdata) ) { | if( $this->n !== count($xdata) ) { |
| JpGraphError::RaiseL(19001); | JpGraphError::RaiseL(19001); |
| //('Spline: Number of X and Y coordinates must be the same'); | //('Spline: Number of X and Y coordinates must be the same'); |
| } | } |
| // Natural spline 2:derivate == 0 at endpoints | // Natural spline 2:derivate == 0 at endpoints |
| $this->y2[0] = 0.0; | $this->y2[0] = 0.0; |
| $this->y2[$n-1] = 0.0; | $this->y2[$n-1] = 0.0; |
| $delta[0] = 0.0; | $delta[0] = 0.0; |
| // Calculate 2:nd derivate | // Calculate 2:nd derivate |
| for($i=1; $i < $n-1; ++$i) { | for($i=1; $i < $n-1; ++$i) { |
| $d = ($xdata[$i+1]-$xdata[$i-1]); | $d = ($xdata[$i+1]-$xdata[$i-1]); |
| if( $d == 0 ) { | if( $d == 0 ) { |
| JpGraphError::RaiseL(19002); | JpGraphError::RaiseL(19002); |
| //('Invalid input data for spline. Two or more consecutive input X-values are equal. Each input X-value must differ since from a mathematical point of view it must be a one-to-one mapping, i.e. each X-value must correspond to exactly one Y-value.'); | //('Invalid input data for spline. Two or more consecutive input X-values are equal. Each input X-value must differ since from a mathematical point of view it must be a one-to-one mapping, i.e. each X-value must correspond to exactly one Y-value.'); |
| } | } |
| $s = ($xdata[$i]-$xdata[$i-1])/$d; | $s = ($xdata[$i]-$xdata[$i-1])/$d; |
| $p = $s*$this->y2[$i-1]+2.0; | $p = $s*$this->y2[$i-1]+2.0; |
| $this->y2[$i] = ($s-1.0)/$p; | $this->y2[$i] = ($s-1.0)/$p; |
| $delta[$i] = ($ydata[$i+1]-$ydata[$i])/($xdata[$i+1]-$xdata[$i]) - | $delta[$i] = ($ydata[$i+1]-$ydata[$i])/($xdata[$i+1]-$xdata[$i]) - |
| ($ydata[$i]-$ydata[$i-1])/($xdata[$i]-$xdata[$i-1]); | ($ydata[$i]-$ydata[$i-1])/($xdata[$i]-$xdata[$i-1]); |
| $delta[$i] = (6.0*$delta[$i]/($xdata[$i+1]-$xdata[$i-1])-$s*$delta[$i-1])/$p; | $delta[$i] = (6.0*$delta[$i]/($xdata[$i+1]-$xdata[$i-1])-$s*$delta[$i-1])/$p; |
| } | } |
| // Backward substitution | // Backward substitution |
| for( $j=$n-2; $j >= 0; --$j ) { | for( $j=$n-2; $j >= 0; --$j ) { |
| $this->y2[$j] = $this->y2[$j]*$this->y2[$j+1] + $delta[$j]; | $this->y2[$j] = $this->y2[$j]*$this->y2[$j+1] + $delta[$j]; |
| } | } |
| } | } |
| // Return the two new data vectors | // Return the two new data vectors |
| function Get($num=50) { | function Get($num=50) { |
| $n = $this->n ; | $n = $this->n ; |
| $step = ($this->xdata[$n-1]-$this->xdata[0]) / ($num-1); | $step = ($this->xdata[$n-1]-$this->xdata[0]) / ($num-1); |
| $xnew=array(); | $xnew=array(); |
| $ynew=array(); | $ynew=array(); |
| $xnew[0] = $this->xdata[0]; | $xnew[0] = $this->xdata[0]; |
| $ynew[0] = $this->ydata[0]; | $ynew[0] = $this->ydata[0]; |
| for( $j=1; $j < $num; ++$j ) { | for( $j=1; $j < $num; ++$j ) { |
| $xnew[$j] = $xnew[0]+$j*$step; | $xnew[$j] = $xnew[0]+$j*$step; |
| $ynew[$j] = $this->Interpolate($xnew[$j]); | $ynew[$j] = $this->Interpolate($xnew[$j]); |
| } | } |
| return array($xnew,$ynew); | return array($xnew,$ynew); |
| } | } |
| // Return a single interpolated Y-value from an x value | // Return a single interpolated Y-value from an x value |
| function Interpolate($xpoint) { | function Interpolate($xpoint) { |
| $max = $this->n-1; | $max = $this->n-1; |
| $min = 0; | $min = 0; |
| // Binary search to find interval | // Binary search to find interval |
| while( $max-$min > 1 ) { | while( $max-$min > 1 ) { |
| $k = ($max+$min) / 2; | $k = ($max+$min) / 2; |
| if( $this->xdata[$k] > $xpoint ) | if( $this->xdata[$k] > $xpoint ) |
| $max=$k; | $max=$k; |
| else | else |
| $min=$k; | $min=$k; |
| } | } |
| // Each interval is interpolated by a 3:degree polynom function | // Each interval is interpolated by a 3:degree polynom function |
| $h = $this->xdata[$max]-$this->xdata[$min]; | $h = $this->xdata[$max]-$this->xdata[$min]; |
| if( $h == 0 ) { | if( $h == 0 ) { |
| JpGraphError::RaiseL(19002); | JpGraphError::RaiseL(19002); |
| //('Invalid input data for spline. Two or more consecutive input X-values are equal. Each input X-value must differ since from a mathematical point of view it must be a one-to-one mapping, i.e. each X-value must correspond to exactly one Y-value.'); | //('Invalid input data for spline. Two or more consecutive input X-values are equal. Each input X-value must differ since from a mathematical point of view it must be a one-to-one mapping, i.e. each X-value must correspond to exactly one Y-value.'); |
| } | } |
| $a = ($this->xdata[$max]-$xpoint)/$h; | $a = ($this->xdata[$max]-$xpoint)/$h; |
| $b = ($xpoint-$this->xdata[$min])/$h; | $b = ($xpoint-$this->xdata[$min])/$h; |
| return $a*$this->ydata[$min]+$b*$this->ydata[$max]+ | return $a*$this->ydata[$min]+$b*$this->ydata[$max]+ |
| (($a*$a*$a-$a)*$this->y2[$min]+($b*$b*$b-$b)*$this->y2[$max])*($h*$h)/6.0; | (($a*$a*$a-$a)*$this->y2[$min]+($b*$b*$b-$b)*$this->y2[$max])*($h*$h)/6.0; |
| } | } |
| } | } |
| //------------------------------------------------------------------------ | //------------------------------------------------------------------------ |
| // CLASS Bezier | // CLASS Bezier |
| // Create a new data array from a number of control points | // Create a new data array from a number of control points |
| //------------------------------------------------------------------------ | //------------------------------------------------------------------------ |
| class Bezier { | class Bezier { |
| /** | /** |
| * @author Thomas Despoix, openXtrem company | * @author Thomas Despoix, openXtrem company |
| * @license released under QPL | * @license released under QPL |
| * @abstract Bezier interoplated point generation, | * @abstract Bezier interoplated point generation, |
| * computed from control points data sets, based on Paul Bourke algorithm : | * computed from control points data sets, based on Paul Bourke algorithm : |
| * http://local.wasp.uwa.edu.au/~pbourke/geometry/bezier/index2.html | * http://local.wasp.uwa.edu.au/~pbourke/geometry/bezier/index2.html |
| */ | */ |
| private $datax = array(); | private $datax = array(); |
| private $datay = array(); | private $datay = array(); |
| private $n=0; | private $n=0; |
| function __construct($datax, $datay, $attraction_factor = 1) { | function __construct($datax, $datay, $attraction_factor = 1) { |
| // Adding control point multiple time will raise their attraction power over the curve | // Adding control point multiple time will raise their attraction power over the curve |
| $this->n = count($datax); | $this->n = count($datax); |
| if( $this->n !== count($datay) ) { | if( $this->n !== count($datay) ) { |
| JpGraphError::RaiseL(19003); | JpGraphError::RaiseL(19003); |
| //('Bezier: Number of X and Y coordinates must be the same'); | //('Bezier: Number of X and Y coordinates must be the same'); |
| } | } |
| $idx=0; | $idx=0; |
| foreach($datax as $datumx) { | foreach($datax as $datumx) { |
| for ($i = 0; $i < $attraction_factor; $i++) { | for ($i = 0; $i < $attraction_factor; $i++) { |
| $this->datax[$idx++] = $datumx; | $this->datax[$idx++] = $datumx; |
| } | } |
| } | } |
| $idx=0; | $idx=0; |
| foreach($datay as $datumy) { | foreach($datay as $datumy) { |
| for ($i = 0; $i < $attraction_factor; $i++) { | for ($i = 0; $i < $attraction_factor; $i++) { |
| $this->datay[$idx++] = $datumy; | $this->datay[$idx++] = $datumy; |
| } | } |
| } | } |
| $this->n *= $attraction_factor; | $this->n *= $attraction_factor; |
| } | } |
| /** | /** |
| * Return a set of data points that specifies the bezier curve with $steps points | * Return a set of data points that specifies the bezier curve with $steps points |
| * @param $steps Number of new points to return | * @param $steps Number of new points to return |
| * @return array($datax, $datay) | * @return array($datax, $datay) |
| */ | */ |
| function Get($steps) { | function Get($steps) { |
| $datax = array(); | $datax = array(); |
| $datay = array(); | $datay = array(); |
| for ($i = 0; $i < $steps; $i++) { | for ($i = 0; $i < $steps; $i++) { |
| list($datumx, $datumy) = $this->GetPoint((double) $i / (double) $steps); | list($datumx, $datumy) = $this->GetPoint((double) $i / (double) $steps); |
| $datax[$i] = $datumx; | $datax[$i] = $datumx; |
| $datay[$i] = $datumy; | $datay[$i] = $datumy; |
| } | } |
| $datax[] = end($this->datax); | $datax[] = end($this->datax); |
| $datay[] = end($this->datay); | $datay[] = end($this->datay); |
| return array($datax, $datay); | return array($datax, $datay); |
| } | } |
| /** | /** |
| * Return one point on the bezier curve. $mu is the position on the curve where $mu is in the | * Return one point on the bezier curve. $mu is the position on the curve where $mu is in the |
| * range 0 $mu < 1 where 0 is tha start point and 1 is the end point. Note that every newly computed | * range 0 $mu < 1 where 0 is tha start point and 1 is the end point. Note that every newly computed |
| * point depends on all the existing points | * point depends on all the existing points |
| * | * |
| * @param $mu Position on the bezier curve | * @param $mu Position on the bezier curve |
| * @return array($x, $y) | * @return array($x, $y) |
| */ | */ |
| function GetPoint($mu) { | function GetPoint($mu) { |
| $n = $this->n - 1; | $n = $this->n - 1; |
| $k = 0; | $k = 0; |
| $kn = 0; | $kn = 0; |
| $nn = 0; | $nn = 0; |
| $nkn = 0; | $nkn = 0; |
| $blend = 0.0; | $blend = 0.0; |
| $newx = 0.0; | $newx = 0.0; |
| $newy = 0.0; | $newy = 0.0; |
| $muk = 1.0; | $muk = 1.0; |
| $munk = (double) pow(1-$mu,(double) $n); | $munk = (double) pow(1-$mu,(double) $n); |
| for ($k = 0; $k <= $n; $k++) { | for ($k = 0; $k <= $n; $k++) { |
| $nn = $n; | $nn = $n; |
| $kn = $k; | $kn = $k; |
| $nkn = $n - $k; | $nkn = $n - $k; |
| $blend = $muk * $munk; | $blend = $muk * $munk; |
| $muk *= $mu; | $muk *= $mu; |
| $munk /= (1-$mu); | $munk /= (1-$mu); |
| while ($nn >= 1) { | while ($nn >= 1) { |
| $blend *= $nn; | $blend *= $nn; |
| $nn--; | $nn--; |
| if ($kn > 1) { | if ($kn > 1) { |
| $blend /= (double) $kn; | $blend /= (double) $kn; |
| $kn--; | $kn--; |
| } | } |
| if ($nkn > 1) { | if ($nkn > 1) { |
| $blend /= (double) $nkn; | $blend /= (double) $nkn; |
| $nkn--; | $nkn--; |
| } | } |
| } | } |
| $newx += $this->datax[$k] * $blend; | $newx += $this->datax[$k] * $blend; |
| $newy += $this->datay[$k] * $blend; | $newy += $this->datay[$k] * $blend; |
| } | } |
| return array($newx, $newy); | return array($newx, $newy); |
| } | } |
| } | } |
| // EOF | // EOF |
| ?> | ?> |