DELPHI - How to make rotation of shape around point

Suppose this easy class for shape - triangle:

TShape = class( TObject )
  p1 : TPoint;
  p2 : TPoint;
  p3 : TPoint;

  procedure Moving( _iDeltaX, _iDeltaY : integer );
  procedure Rotation( _iAngle : integer );
  procedure RotationAroundXY( _x, _y : integer; _iAngle : integer );
end;

And initialization:

pShape := TShape.Create;
pShape.p1.X := 300;  pShape.p1.Y := 250;
pShape.p2.X := 250;  pShape.p2.Y := 280;
pShape.p3.X := 350;  pShape.p3.Y := 330;

For rotation shape around point use this code:

procedure TShape.RotationAroundXY( _x, _y : integer; _iAngle : integer );
var
  iNewX, iNewY : integer;
  dCos, dSin : double;
  dRadian : double;
  p1x, p1y, p2x, p2y, p3x, p3y : integer;
begin
  { -- angle }

  dRadian := ( 2 * pi ) / ( 360 / _iAngle );
                     
  dCos := cos( dRadian );
  dSin := sin( dRadian );

  { -- point for rotation around (in other case the rotation is around start of coords) }

  p1x := p1.x - _x;
  p1y := p1.y - _y;
  p2x := p2.x - _x;
  p2y := p2.y - _y;
  p3x := p3.x - _x;
  p3y := p3.y - _y;

  { -- calculation of the positions }

  iNewX := round( ( p1x * dCos ) - ( p1y * dSin ) );
  iNewY := round( ( p1y * dCos ) + ( p1x * dSin ) );
  p1.X := iNewX + _x;
  p1.Y := iNewY + _y;

  iNewX := round( ( p2x * dCos ) - ( p2y * dSin ) );
  iNewY := round( ( p2y * dCos ) + ( p2x * dSin ) );
  p2.X := iNewX + _x;
  p2.Y := iNewY + _y;

  iNewX := round( ( p3x * dCos ) - ( p3y * dSin ) );
  iNewY := round( ( p3y * dCos ) + ( p3x * dSin ) );
  p3.X := iNewX + _x;
  p3.Y := iNewY + _y;
end;

Output:

DELPHI - How get greatest common divisor (for two digits)

The function returns greatest common divisor for two digits. When it does not found it, returns 1.

function GCD( _iNumber1, _iNumber2 : integer ) : integer;
var
  iTemp : integer;
begin
  if _iNumber1 < 0 then _iNumber1 := -_iNumber1;
  if _iNumber2 < 0 then _iNumber2 := -_iNumber2;

  repeat

    if _iNumber1 < _iNumber2 then
      begin
        iTemp := _iNumber1;
        _iNumber1 := _iNumber2;
        _iNumber2 := iTemp;
      end;

    _iNumber1 := _iNumber1 mod _iNumber2;

  until ( _iNumber1 = 0 );

  result := _iNumber2;
end;

Calling:

var
  i : integer;
begin
  i := GCD( 12, 16 );

  ShowMessage( IntToStr( i ) );
end;

Output:

DELPHI - How to show chi-square distribution curve

Propability density function (chi-square distribution).

var
  i, iK : integer;
  iX, iY, iGamma, iNumerator, iDenominator : double;
  pSerie : TLineSeries;
  pIntegral : TIntegral;
begin
  Memo1.Clear;

  pSerie := TLineSeries( Chart1.Series[0] );
  pSerie.Clear;

  pIntegral := TIntegral.Create;
  pIntegral.iStep := 0.05;
  pIntegral.SetInterval( 0.001, 20 );

  { -- set param + calc GAMMA function }

  iK := 8;
  iGamma := GetGamma( iK/2 );

  for i := 0 to pIntegral.pX.Count - 1 do
    begin
      iX := pIntegral.pX.GetValue( i );

      try
        iNumerator := power( iX, iK/2-1 ) * exp( ( -1*iX ) /2 );
        iDenominator := power( 2, iK / 2 ) * iGamma;

        if iDenominator <> 0 then
          iY := iNumerator / iDenominator
        else
          iY := 0;
      except
        iY := 0;
      end;
      pIntegral.pY.AddValue( iY );

      { -- add to chart }

      pSerie.AddXY( iX, iY );
    end;
end;

And here is very important Gamma() function with k parameter.

function GammaStirF( X : double ) : double;
var
  y : double;
  w : double;
  v : double;
  stir : double;
begin
  w := 1 / x;
  stir := 7.87311395793093628397E-4;
  stir := -2.29549961613378126380E-4 + w * stir;
  stir := -2.68132617805781232825E-3 + w * stir;
  stir := 3.47222221605458667310E-3 + w * stir;
  stir := 8.33333333333482257126E-2 + w * stir;
  w := 1 + w * stir;
  y := exp(x);
  if x > 143.01608 then
    begin
      v := power( x, 0.5 * x -0.25 );
      y := v *( v / y );
    end
  else
    begin
      y := power( x, x -0.5 ) / y;
    end;
  result := 2.50662827463100050242 * y * w;
end;

function GetGamma( _x : double ) : double;
var
  p : double;
  PP : double;
  q : double;
  QQ : double;
  z : double;
  i : longint;
  SgnGam : double;
begin
  SgnGam := 1;
  q := abs( _x );
  if q > 33.0 then
  begin
    if _x < 0.0 then
    begin
      p := floor(q);
      i := round(p);

      if i mod 2 = 0 then
      begin
        SgnGam := -1;
      end;

      z := q - p;

      if z > 0.5 then
      begin
        p := p+1;
        z := q-p;
      end;

      z := q * sin( pi * z );
      z := Abs(z);
      z := pi / ( z * GammaStirF( q ) );
    end
    else
    begin
      z := GammaStirF( _x );
    end;
    result := SgnGam * z;
    exit;
  end;
  z := 1;
  while _x >= 3 do
  begin
    _x := _x - 1;
    z := z *_x;
  end;
  while _x < 0 do
  begin
    if _x > -0.000000001 then
    begin
      result := z / ( ( 1 + 0.5772156649015329 *_x ) * _x );
      exit;
    end;
    z := z /_x;
    _x := _x + 1;
  end;
  while _x < 2 do
  begin
    if _x < 0.000000001 then
    begin
      result := z / ( ( 1 + 0.5772156649015329 * _x ) *_x );
      exit;
    end;
    z := z /_x;
    _x := _x + 1.0;
  end;
  if _x = 2 then
  begin
    result := z;
    exit;
  end;

  _x := _x - 2.0;
  PP := 1.60119522476751861407E-4;
  PP := 1.19135147006586384913E-3 + _X * PP;
  PP := 1.04213797561761569935E-2 + _X * PP;
  PP := 4.76367800457137231464E-2 + _X * PP;
  PP := 2.07448227648435975150E-1 + _X * PP;
  PP := 4.94214826801497100753E-1 + _X * PP;
  PP := 9.99999999999999996796E-1 + _X * PP;
  QQ := -2.31581873324120129819E-5;
  QQ := 5.39605580493303397842E-4 + _X * QQ;
  QQ := -4.45641913851797240494E-3 + _X * QQ;
  QQ := 1.18139785222060435552E-2 + _X * QQ;
  QQ := 3.58236398605498653373E-2 + _X * QQ;
  QQ := -2.34591795718243348568E-1 + _X * QQ;
  QQ := 7.14304917030273074085E-2 + _X * QQ;
  QQ := 1.00000000000000000320 + _X * QQ;

  result := z * PP / QQ;
  exit;
end;

Output:

k=2

k=8

DELPHI - How to show normal distribution curve (Gaussian)

Typical propability density function (Gaussian distribution).

var
  i : integer;
  iY : double;
  pDN : TDistribution_Normal;
  pSerie : TLineSeries;
begin
  Memo1.Clear;

  /* prepare chart */

  pSerie := TLineSeries( Chart1.Series[0] );
  pSerie.Clear;

  /* prepare distribution data */

  pDN := TDistribution_Normal.Create;
  pDN.Init( 0, 1 );

  for i := 0 to pDN.pIntegral.pY.Count - 1 do
    begin
      { -- get value }

      iY := pDN.pIntegral.pY.GetValue( i );

      pSerie.AddXY( i, iY );
    end;

And here is very important TDistribution_Normal.Init() function with parameters as MEAN and SIGMA, here is generated range MEAN +-3.5 SIGMA.
The sum values under curve is equeal to 1.

function TDistribution_Normal.Init( _iMean, _iSigma : double ) : boolean;
var
  i : integer;
  iFrom, iTo, iValue : double;
begin
  result := false;

  { -- check }

  if _iSigma <= 0 then
    begin
      ShowMessage( 'Distribution_Normal() - standard deviation can`t be 0.' );
      exit;
    end;

  { -- add range -3.5 sigma - Mean - + 3.5 sigma }

  pIntegral.Clear;

  iFrom := _iMean - ( 3.5 * _iSigma );
  iTo := _iMean + ( 3.5 * _iSigma );
  { aproximation - 300 items }
  pIntegral.iStep := ( iTo - iFrom ) / 300;
  pIntegral.SetInterval( iFrom, iTo );

  { -- calc values normal distribution for every X }

  for i := 0 to pIntegral.pX.Count - 1 do
    begin
      iValue := ( 1 / ( _iSigma * sqrt( 2*pi ) ) ) * 
                exp( - ( sqr( pIntegral.pX.GetValue( i ) - _iMean ) ) / 
                ( 2 * sqr( _iSigma ) ) );

      pIntegral.pY.AddValue( iValue );
    end;

  pIntegral.Calc;

  result := true;
end;

Output:

DELPHI - How compute integral - (sin(x))^2 dx from 0 to PI

For integration is used class TIntegral which use for calculation trapezoidal rule.

procedure TForm1.BIntegralClick(Sender: TObject);
var
  i : integer;
  iX, iY : double;
  pIntegral : TIntegral;
  pSerie : TLineSeries;
begin
    pSerie := TLineSeries( Chart1.Series[0] );
    pSerie.Clear;

    { -- prepare integral class - interval 0..PI }

    pIntegral := TIntegral.Create;
    pIntegral.SetInterval( 0, pi );

    for i := 0 to pIntegral.pX.Count - 1 do
      begin
        { calc Y value }
        pIntegral.pY.AddValue( sqr( sin( pIntegral.pX.GetValue( i ) ) ) );

        { -> to graph }
        iX := pIntegral.pX.GetValue( i );
        iY := pIntegral.pY.GetValue( i );
        pSerie.AddXY( iX, iY );
      end;

    ShowMessage( FloatToStr( pIntegral.Calc ) );        // 1.57

    pIntegral.Free;
end;

Here is TIntegral.calc() method:

{ ---------------------------------------------------------------------------
  Method computes integral value (trapezoidal rule).  
  -------------------------------------------------------------------------- }
function TIntegral.Calc( _iCheckSumIntegralValue : double = 999 ) : double;
var
  i, iDeleteXFrom : integer;
  iValue : double;
begin
  result := 0;

  { -- checks }

  if pX.Count < 2 then exit;
  if pX.Count <> pY.Count then exit;

  { -- delete values }

  pIntegral.Clear;
  pCumulation.Clear;

  { -- trapezoidal rule }

  iDeleteXFrom := -1;

  for i := 1 to pX.Count - 1 do
    begin
      { -- compute integral value in interval }

      iValue := (
                  ( pX.GetValueCheckNull( i ) - pX.GetValueCheckNull( i - 1 ) )
                  *
                  ( ( pY.GetValueCheckNull( i ) + pY.GetValueCheckNull( i - 1 ) ) / 2 )
                );

       result := result + iValue;

       { -- check -> is greater than parameter ? }

       if _iCheckSumIntegralValue <> 999 then
       if result > _iCheckSumIntegralValue then
         begin
           { -- subtract }
           result := result - iValue;

           { -- index for removing }
           iDeleteXFrom := i;

           break;
         end;

       { -- record values }

       if bLogIntegral then pIntegral.AddValue( iValue );
       if bLogCumulation then pCumulation.AddValue( result );
     end;

  { -- delete items over test value (if exists) }

  if iDeleteXFrom <> -1 then
    for i := pX.Count - 1 downto iDeleteXFrom do
      begin
        pX.Delete( i );
        pY.Delete( i );
      end;

end;

Output:

And check (same result):