Calculate a skew matrix using DEG angles
@param ax Skew on axis x
@param ay Skew on axis y
@returns {{a: number, b: number, c: number, e: number, d: number, f: number}} Affine matrix
494740
function skewDEG(ax, ay) {
    return skew(ax * Math.PI / 180, ay * Math.PI / 180);
}
808943
public static final DoubleMatrix2D invert2x2Matrix(DoubleMatrix2D A) throws Exception{
		if(2!=A.rows() || A.rows() != A.columns()){
			throw new IllegalArgumentException("matrix is not 2x2");
		}
		//ad - bc
		double s = A.getQuick(0, 0)*A.getQuick(1, 1)-A.getQuick(0, 1)*A.getQuick(1, 0);
		if(Math.abs(s) < 1.e-16){
			throw new Exception("Matrix is singular");
		}
		DoubleMatrix2D ret = new DenseDoubleMatrix2D(2, 2);
		ret.setQuick(0, 0, A.getQuick(1, 1) / s);
		ret.setQuick(1, 1, A.getQuick(0, 0) / s);
		ret.setQuick(0, 1, -A.getQuick(1, 0) / s);
		ret.setQuick(1, 0, A.getQuick(0, 1) / s);
		return ret;
	}
564823
public static DoubleMatrix ger(double alpha, DoubleMatrix x,
			DoubleMatrix y, DoubleMatrix a) {
		NativeBlas.dger(a.rows, a.columns, alpha, x.data, 0, 1, y.data, 0, 1, a.data,
				0, a.rows);
		return a;
	}
979471
public static double elementMinAbs( DMatrix5 a ) {
        double min = Math.abs(a.a1);
        double tmp = Math.abs(a.a1); if( tmp < min ) min = tmp;
        tmp = Math.abs(a.a2); if( tmp < min ) min = tmp;
        tmp = Math.abs(a.a3); if( tmp < min ) min = tmp;
        tmp = Math.abs(a.a4); if( tmp < min ) min = tmp;
        tmp = Math.abs(a.a5); if( tmp < min ) min = tmp;

        return min;
    }
980370
function(matrix) {
      var skew = this._skew;

      if (!skew) {
        skew = this._createElement('skew');
        this._container.appendChild(skew);
        skew.style.behavior = 'url(#default#VML)';
        this._skew = skew;
      }

      // handle skew/translate separately, cause it's broken
      var mt = matrix[0].toFixed(8) + " " + matrix[1].toFixed(8) + " " +
        matrix[2].toFixed(8) + " " + matrix[3].toFixed(8) + " 0 0";
      var offset = Math.floor(matrix[4]).toFixed() + ", " +
        Math.floor(matrix[5]).toFixed() + "";

      var s = this._container.style;
      var l = parseFloat(s.left);
      var t = parseFloat(s.top);
      var w = parseFloat(s.width);
      var h = parseFloat(s.height);

      if (isNaN(l)) l = 0;
      if (isNaN(t)) t = 0;
      if (isNaN(w) || !w) w = 1;
      if (isNaN(h) || !h) h = 1;

      var origin = (-l / w - 0.5).toFixed(8) + " " + (-t / h - 0.5).toFixed(8);

      skew.on = "f";
      skew.matrix = mt;
      skew.origin = origin;
      skew.offset = offset;
      skew.on = true;
    }
484576
public static DMatrixRBlock transpose(DMatrixRBlock A , DMatrixRBlock A_tran )
    {
        if( A_tran != null ) {
            if( A.numRows != A_tran.numCols || A.numCols != A_tran.numRows )
                throw new IllegalArgumentException("Incompatible dimensions.");
            if( A.blockLength != A_tran.blockLength )
                throw new IllegalArgumentException("Incompatible block size.");
        } else {
            A_tran = new DMatrixRBlock(A.numCols,A.numRows,A.blockLength);

        }

        for( int i = 0; i < A.numRows; i += A.blockLength ) {
            int blockHeight = Math.min( A.blockLength , A.numRows - i);

            for( int j = 0; j < A.numCols; j += A.blockLength ) {
                int blockWidth = Math.min( A.blockLength , A.numCols - j);

                int indexA = i*A.numCols + blockHeight*j;
                int indexC = j*A_tran.numCols + blockWidth*i;

                transposeBlock( A , A_tran , indexA , indexC , blockWidth , blockHeight );
            }
        }

        return A_tran;
    }
37126
public static double det( DMatrix4x4 mat ) {

        double  a11 = mat.a22;
        double  a12 = mat.a23;
        double  a13 = mat.a24;
        double  a21 = mat.a32;
        double  a22 = mat.a33;
        double  a23 = mat.a34;
        double  a31 = mat.a42;
        double  a32 = mat.a43;
        double  a33 = mat.a44;

        double ret = 0;
        ret += mat.a11 * ( + a11*(a22*a33 - a23*a32) - a12*(a21*a33 - a23*a31) + a13*(a21*a32 - a22*a31));
        a11 = mat.a21;
        a21 = mat.a31;
        a31 = mat.a41;
        ret -= mat.a12 * ( + a11*(a22*a33 - a23*a32) - a12*(a21*a33 - a23*a31) + a13*(a21*a32 - a22*a31));
        a12 = mat.a22;
        a22 = mat.a32;
        a32 = mat.a42;
        ret += mat.a13 * ( + a11*(a22*a33 - a23*a32) - a12*(a21*a33 - a23*a31) + a13*(a21*a32 - a22*a31));
        a13 = mat.a23;
        a23 = mat.a33;
        a33 = mat.a43;
        ret -= mat.a14 * ( + a11*(a22*a33 - a23*a32) - a12*(a21*a33 - a23*a31) + a13*(a21*a32 - a22*a31));

        return ret;
    }
757325
public static FloatMatrix gemv(float alpha, FloatMatrix a,
			FloatMatrix x, float beta, FloatMatrix y) {
		if (false) {
			NativeBlas.sgemv('N', a.rows, a.columns, alpha, a.data, 0, a.rows, x.data, 0,
					1, beta, y.data, 0, 1);
		} else {
			if (beta != 0.0f) {
				for (int i = 0; i < y.length; i++)
					y.data[i] = beta * y.data[i];
			} else {
				for (int i = 0; i < y.length; i++)
					y.data[i] = 0.0f;
			}

	
			for (int j = 0; j < a.columns; j++) {
				float xj = x.get(j);
				if (xj != 0.0f) {
					for (int i = 0; i < a.rows; i++)
							y.data[i] += alpha * a.get(i, j) * xj;
				}
			}
		}
		return y;
	}
1010264
public Matrix3f rotateXYZ(float angleX, float angleY, float angleZ, Matrix3f dest) {
        float sinX = (float) Math.sin(angleX);
        float cosX = (float) Math.cosFromSin(sinX, angleX);
        float sinY = (float) Math.sin(angleY);
        float cosY = (float) Math.cosFromSin(sinY, angleY);
        float sinZ = (float) Math.sin(angleZ);
        float cosZ = (float) Math.cosFromSin(sinZ, angleZ);
        float m_sinX = -sinX;
        float m_sinY = -sinY;
        float m_sinZ = -sinZ;

        // rotateX
        float nm10 = m10 * cosX + m20 * sinX;
        float nm11 = m11 * cosX + m21 * sinX;
        float nm12 = m12 * cosX + m22 * sinX;
        float nm20 = m10 * m_sinX + m20 * cosX;
        float nm21 = m11 * m_sinX + m21 * cosX;
        float nm22 = m12 * m_sinX + m22 * cosX;
        // rotateY
        float nm00 = m00 * cosY + nm20 * m_sinY;
        float nm01 = m01 * cosY + nm21 * m_sinY;
        float nm02 = m02 * cosY + nm22 * m_sinY;
        dest.m20 = m00 * sinY + nm20 * cosY;
        dest.m21 = m01 * sinY + nm21 * cosY;
        dest.m22 = m02 * sinY + nm22 * cosY;
        // rotateZ
        dest.m00 = nm00 * cosZ + nm10 * sinZ;
        dest.m01 = nm01 * cosZ + nm11 * sinZ;
        dest.m02 = nm02 * cosZ + nm12 * sinZ;
        dest.m10 = nm00 * m_sinZ + nm10 * cosZ;
        dest.m11 = nm01 * m_sinZ + nm11 * cosZ;
        dest.m12 = nm02 * m_sinZ + nm12 * cosZ;
        return dest;
    }
271182
private static double[] balance(DenseMatrix A) {
        double sqrdx = Math.RADIX * Math.RADIX;

        int n = A.nrows();

        double[] scale = new double[n];
        for (int i = 0; i < n; i++) {
            scale[i] = 1.0;
        }

        boolean done = false;
        while (!done) {
            done = true;
            for (int i = 0; i < n; i++) {
                double r = 0.0, c = 0.0;
                for (int j = 0; j < n; j++) {
                    if (j != i) {
                        c += Math.abs(A.get(j, i));
                        r += Math.abs(A.get(i, j));
                    }
                }
                if (c != 0.0 && r != 0.0) {
                    double g = r / Math.RADIX;
                    double f = 1.0;
                    double s = c + r;
                    while (c < g) {
                        f *= Math.RADIX;
                        c *= sqrdx;
                    }
                    g = r * Math.RADIX;
                    while (c > g) {
                        f /= Math.RADIX;
                        c /= sqrdx;
                    }
                    if ((c + r) / f < 0.95 * s) {
                        done = false;
                        g = 1.0 / f;
                        scale[i] *= f;
                        for (int j = 0; j < n; j++) {
                            A.mul(i, j, g);
                        }
                        for (int j = 0; j < n; j++) {
                            A.mul(j, i, f);
                        }
                    }
                }
            }
        }

        return scale;
    }
333305
public static double elementMinAbs( DMatrix3x3 a ) {
        double min = Math.abs(a.a11);
        double tmp = Math.abs(a.a12); if( tmp < min ) min = tmp;
        tmp = Math.abs(a.a13); if( tmp < min ) min = tmp;
        tmp = Math.abs(a.a21); if( tmp < min ) min = tmp;
        tmp = Math.abs(a.a22); if( tmp < min ) min = tmp;
        tmp = Math.abs(a.a23); if( tmp < min ) min = tmp;
        tmp = Math.abs(a.a31); if( tmp < min ) min = tmp;
        tmp = Math.abs(a.a32); if( tmp < min ) min = tmp;
        tmp = Math.abs(a.a33); if( tmp < min ) min = tmp;

        return min;
    }
902024
public static AxisAlignedBB rotate(AxisAlignedBB aabb, int angle)
	{
		return rotate(aabb, angle, Axis.Y);
	}
835649
def _create_affine(x_axis, y_axis, z_axis, image_pos, voxel_sizes):
    

    # Create affine matrix (http://nipy.sourceforge.net/nibabel/dicom/dicom_orientation.html#dicom-slice-affine)

    affine = numpy.array(
        [[x_axis[0] * voxel_sizes[0], y_axis[0] * voxel_sizes[1], z_axis[0] * voxel_sizes[2], image_pos[0]],
         [x_axis[1] * voxel_sizes[0], y_axis[1] * voxel_sizes[1], z_axis[1] * voxel_sizes[2], image_pos[1]],
         [x_axis[2] * voxel_sizes[0], y_axis[2] * voxel_sizes[1], z_axis[2] * voxel_sizes[2], image_pos[2]],
         [0, 0, 0, 1]])
    return affine
509502
public static ComplexFloatMatrix gerc(ComplexFloat alpha, ComplexFloatMatrix x,
			ComplexFloatMatrix y, ComplexFloatMatrix a) {
		NativeBlas.cgerc(a.rows, a.columns, alpha, x.data, 0, 1, y.data, 0, 1, a.data,
				0, a.rows);
		return a;
	}
944245
function skew($angle_x, $angle_y, $x, $y)
    {
        $y = $this->currentPageSize["height"] - $y;

        $tan_x = tan(deg2rad($angle_x));
        $tan_y = tan(deg2rad($angle_y));

        $tm = array(
            1,
            -$tan_y,
            -$tan_x,
            1,
            $tan_x * $y,
            $tan_y * $x,
        );

        $this->transform($tm);
    }
358962
public Matrix4f rotateTowardsXY(float dirX, float dirY, Matrix4f dest) {
        if ((properties & PROPERTY_IDENTITY) != 0)
            return dest.rotationTowardsXY(dirX, dirY);
        float rm00 = dirY;
        float rm01 = dirX;
        float rm10 = -dirX;
        float rm11 = dirY;
        float nm00 = m00 * rm00 + m10 * rm01;
        float nm01 = m01 * rm00 + m11 * rm01;
        float nm02 = m02 * rm00 + m12 * rm01;
        float nm03 = m03 * rm00 + m13 * rm01;
        dest._m10(m00 * rm10 + m10 * rm11);
        dest._m11(m01 * rm10 + m11 * rm11);
        dest._m12(m02 * rm10 + m12 * rm11);
        dest._m13(m03 * rm10 + m13 * rm11);
        dest._m00(nm00);
        dest._m01(nm01);
        dest._m02(nm02);
        dest._m03(nm03);
        dest._m20(m20);
        dest._m21(m21);
        dest._m22(m22);
        dest._m23(m23);
        dest._m30(m30);
        dest._m31(m31);
        dest._m32(m32);
        dest._m33(m33);
        dest._properties(properties & ~(PROPERTY_PERSPECTIVE | PROPERTY_IDENTITY | PROPERTY_TRANSLATION));
        return dest;
    }
157835
function(angle) {
		// TODO this clears any scale, should we care?
		var a = Math.cos(angle);
		var b = Math.sin(angle);
		this.matrix[0] = this.matrix[3] = a;
		this.matrix[1] = (0-b);
		this.matrix[2] = b;
	}