Moyen rapide / efficace de décomposer des coefficients de filtre 2D entiers séparables

Je voudrais pouvoir déterminer rapidement si un noyau 2D donné de coefficients entiers est séparable en deux noyaux 1D avec des coefficients entiers. Par exemple

 2   3   2
 4   6   4
 2   3   2

est séparable en

 2   3   2

et

 1
 2
 1

Le test réel de séparabilité semble être assez simple en utilisant l'arithmétique entière, mais la décomposition en filtres 1D avec des coefficients entiers s'avère être un problème plus difficile. La difficulté semble résider dans le fait que les ratios entre les lignes ou les colonnes peuvent être non entiers (fractions rationnelles), par exemple dans l'exemple ci-dessus, nous avons des ratios de 2, 1/2, 3/2 et 2/3.

Je ne veux pas vraiment utiliser une approche lourde comme SVD parce que (a) c'est relativement coûteux en calcul pour mes besoins et (b) ça n'aide pas toujours nécessairement à déterminer des coefficients entiers .

Des idées ?

PLUS D'INFORMATIONS

Les coefficients peuvent être positifs, négatifs ou nuls, et il peut y avoir des cas pathologiques où la somme de l'un ou des deux vecteurs 1D est nulle, par exemple

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

est séparable en

 1  -2   1

et

-1
 0
 1

J'ai pris @Phononla réponse de et l ' ai quelque peu modifiée pour qu'elle utilise l' approche GCD uniquement sur la ligne du haut et la colonne de gauche, plutôt que sur les sommes de ligne / colonne. Cela semble mieux gérer les cas pathologiques. Il peut toujours échouer si la ligne supérieure ou la colonne de gauche sont toutes des zéros, mais ces cas peuvent être vérifiés avant d'appliquer cette méthode.

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

Un grand merci à @Phononet @Jason Rpour les idées originales pour cela.

Je l'ai! Publication du code MATLAB, affichera une explication ce soir ou demain

% Two original arrays
N = 3;
range = 800;
a = round( range*(rand(N,1)-0.5) )
b = round( range*(rand(1,N)-0.5) )

% Create a matrix;
M = a*b;
N = size(M,1);

% Sanity check
disp([num2str(rank(M)) ' <- this should be 1!']);

% Sum across rows and columns
Sa = M * ones(N,1);
Sb = ones(1,N) * M;

% Get rid of zeros
SSa = Sa( Sa~=0 );
SSb = Sb( Sb~=0 );

if isempty(SSa) | isempty(SSb)
    break;
end

% Sizes of array without zeros
Na = numel(SSa);
Nb = numel(SSb);

% Find Greatest Common Divisor of Sa and Sb.
Ga = SSa(1);
Gb = SSb(1);

for l=2:Na
    Ga = gcd(Ga,SSa(l));
end

for l=2:Nb
    Gb = gcd(Gb,SSb(l));
end

%Divide by the greatest common divisor
Sa = Sa / Ga;
Sb = Sb / Gb;

%Scale one of the vectors
MM = Sa * Sb;
Sa = Sa * (MM(1) / M(1));

disp('Two arrays found:')
Sa
Sb
disp('Sa * Sb = ');
Sa*Sb
disp('Original = ');
M

Peut-être que je banalise le problème, mais il semble que vous pourriez:

• $N$$M$$\mathbf{A}$$\mathbf{a_i}$$i = 0, 1, \ldots , N-1$

• $j > 0$

  • $\mathbf{a_j}$$\mathbf{a_0}$$j$$\mathbf{r_j}$
  • $\mathbf{r_j}$
  • $\mathbf{r_j}$$\mathbf{a_j}$$\mathbf{a_0}$$j$$0$$\mathbf{x}$
  • $\mathbf{a_j}$$\mathbf{a_0}$

• $\mathbf{x}$

• $\mathbf{x_{norm}}$ $$\mathbf{x_{scaled}} = K\mathbf{x_{norm}}, K = 1, 2, \ldots, M$$ $K$$M$

Ce n'est pas la méthode la plus élégante, et il y a probablement une meilleure façon, mais cela devrait fonctionner, est assez simple à implémenter et devrait être relativement rapide pour les matrices de taille modeste.

$x \otimes y \otimes z$$A$$|A - x \otimes y \otimes z|$
$x\ y\ z$
$y$$z$$x$$x\ y\ z\ x\ y\ z\ ...$ à son tour.

(De approximative-a-convolution-as-a-sum-of-separable-convolutions sur math.stackexchange.)