separability - Fast/efficient way to decompose separable integer 2D filter coefficients

I would like to be able to quickly determine whether a given 2D kernel of integer coefficients is separable into two 1D kernels with integer coefficients. E.g.

 2   3   2
 4   6   4
 2   3   2

is separable into

 2   3   2

and

 1
 2
 1

The actual test for separability seems to be fairly straightforward using integer arithmetic, but the decomposition into 1D filters with integer coefficients is proving to be a more difficult problem. The difficulty seems to lie in the fact that ratios between rows or columns may be non-integer (rational fractions), e.g. in the above example we have ratios of 2, 1/2, 3/2 and 2/3.

I don't really want to use a heavy duty approach like SVD because (a) it's relatively computationally expensive for my needs and (b) it still doesn't necessarily help to determine integer coefficients.

Any ideas ?

---

FURTHER INFORMATION

Coefficients may be positive, negative or zero, and there may be pathological cases where the sum of either or both 1D vectors is zero, e.g.

-1   2  -1
 0   0   0
 1  -2   1

is separable into

 1  -2   1

and

-1
 0
 1

Answer

I have taken @Phonon's answer and modified it somewhat so that it uses the GCD approach on just the top row and left column, rather than on row/column sums. This seems to handle pathological cases a little better. It can still fail if the top row or left column are all zeroes, but these cases can be checked for prior to applying this method.

function [X, Y, valid] = separate(M)    % separate 2D kernel M into X and Y vectors 
  X = M(1, :);                          % init X = top row of M
  Y = M(:, 1);                          % init Y = left column of M

  nx = numel(X);                        % nx = no of columns in M
  ny = numel(Y);                        % ny = no of rows in M
  gx = X(1);                            % gx = GCD of top row
  for i = 2:nx
    gx = gcd(gx, X(i));
  end
  gy = Y(1);                            % gy = GCD of left column
  for i = 2:ny
    gy = gcd(gy, Y(i));
  end

  X = X / gx;                           % scale X by GCD of X
  Y = Y / gy;                           % scale Y by GCD of Y
  scale = M(1, 1) / (X(1) * Y(1));      % calculate scale factor
  X = X * scale;                        % apply scale factor to X
  valid = all(all((M == Y * X)));       % result valid if we get back our original M
end

Many thanks to @Phonon and @Jason R for the original ideas for this.