procedure RotateModel(ModelInstance: TModelInstance; const Scale: Double);
var
  OrigVector: TModelVector;
  RotatedVector: TModelVector;
  Vector: TVector;
  ScaleMatrix: TMatrix;
  V: TVector;
begin
  ModelInstance.RotationMatrixX := CreateRotationMatrixX(Sin(DegToRad(ModelInstance.RotX)), Cos(DegToRad(ModelInstance.RotX)));
  ModelInstance.RotationMatrixY := CreateRotationMatrixY(Sin(DegToRad(ModelInstance.RotY - 90)), Cos(DegToRad(ModelInstance.RotY - 90)));
  ModelInstance.RotationMatrixZ := CreateRotationMatrixZ(Sin(DegToRad(-ModelInstance.RotZ)), Cos(DegToRad(-ModelInstance.RotZ)));
  V[0] := Scale;
  V[1] := Scale;
  V[2] := Scale;
  ScaleMatrix := CreateScaleMatrix(V);
  ModelInstance.RotationMatrix := MatrixMultiply(MatrixMultiply(ModelInstance.RotationMatrixX, ModelInstance.RotationMatrixZ), ModelInstance.RotationMatrixY);
  ModelInstance.RotationMatrix := MatrixMultiply(ScaleMatrix, ModelInstance.RotationMatrix);

  ModelInstance.OrigModelVertices.First;
  while not ModelInstance.OrigModelVertices.EOF do
  begin
    OrigVector := TModelVector(ModelInstance.OrigModelVertices.CurrentNode);
    Vector[0] := -OrigVector.X;
    Vector[1] := -OrigVector.Y;
    Vector[2] := OrigVector.Z;
    Vector[3] := 0;

    Vector := VectorTransform(Vector, ModelInstance.RotationMatrix);

    RotatedVector := TModelVector.Create(nil);
    RotatedVector.Index := OrigVector.Index;
    RotatedVector.X := Vector[0];
    RotatedVector.Y := Vector[1];
    RotatedVector.Z := Vector[2];

    ModelInstance.RotatedModelVertices.Add(RotatedVector);

    ModelInstance.OrigModelVertices.Next;
  end;
end;