Bonjour à tous,

J'aurais une question d'algèbre matricielle/programmation assez pointue.

Voici le contexte: je travaille actuellement sur un code de déconvolution avec une fonction d'étalement du point (PSF - Point Spread Function) qui varie dans l'espace. Typiquement il s'agit de résoudre un problème du type:



où est une image (floutée) observée (connue), est la PSF (qui est connue), est l'image nette (inconnue) à retrouver, est l'opérateur de convolution et un terme de bruit (que l'on peut ignorer ici).

En général on peut réécrire l'équation ci-dessus sous forme d'une équation linéaire du type:



où est une matrice circulante.


Le fait que soit circulante est très utile: on peut alors facilement et rapidement résoudre le problème ci-dessus en passant dans l'espace de Fourier en profitant de la FFT (Fast Fourier Transform). (voir par exemple: https://fr.wikipedia.org/wiki/Matrice_circulante)

Cependant, si la PSF varie dans l'espace, n'est en générale plus circulante. Donc on ne peut plus profiter directement des avantages de la FFT.
J'ai bien implémenté un code qui résout le problème ci-dessus pour une PSF spatialement variable; mais celui-ci implique une résolution itérative d'un système linéaire qui est un peu trop lente à mon goût.

Je me suis donc un peu renseigné dans la littérature et j'ai notamment trouvé cela:

https://link.springer.com/article/10...041-015-9395-0

Il est apparemment possible de décomposer les matrices carrées* en produit de matrices circulantes et diagonales. Je pense qu'il serait alors à nouveau possible de profiter des avantages de la FFT.

*L'algorithme de déconvolution que j'utilise considère la matrice qui est toujours carrée.

Comme je suis un peu paresseux, je me demandais si l'un d'entre-vous aurait connaissance d'un code (peu importe le langage) ou d'un algorithme (même décrit en anglais) qui effectuerais cette décomposition ou une décomposition similaire.