# factMeth

Computes an image of the shape of the scatterer from far field pattern using the factorization method.

## Contents

- Syntax
- Description
- Input Arguments
- Output Arguments
- More About
- References
- See Also
- Code: function: factMeth
- Code: function: factMeth: Construct F_# and Compute the SVD
- Code: function: factMeth: Numerical Implementation of Factorization Method

## Syntax

    z = factMeth(seti);

## Description

`z = factMeth(seti)` computes a real matrix `z`.

## Input Arguments

- `seti` : structure array

- `seti.incNb` : number of transmitters (number of incident fields)
- `seti.delta` : Relative (artificial) noise level of data (default: 0.01)
- `seti.FmeasDelta` : scattered field at receivers positions (i.e. 'the data') (uScaRX = uTotRX - uIncRX, i.e. total field minus incident field) (complex matrix of size seti.measNb x seti.incNb)

Note that the factorization method is restricted to plane waves and far field data.

## Output Arguments

- `z` : Result of the numercial implementation of the factorization method. (real matrix of size seti.nROI x seti.nROI)

## More About

Note that the factorization method's result separates the obstacle from the background but does not distinguish whether the contrast is different from zero in the real or imaginary part. This is important for the interpretation of the result.

## References

- [1] Armin Lechleiter. Factorization Methods for Photonics and Rough Surfaces. PhD thesis, Universität Karlsruhe, 2008. URL: https://doi.org/10.5445/KSP/1000008735.
- [2] Willy Dörfler, Armin Lechleiter, Michael Plum, Guido Schneider, and Christian Wieners. Photonic Crystals: Mathematical Analysis and Numerical Approximation, volume 42 of Oberwolfach Seminars. Birkhäuser, Basel, 2011. URL: https://doi.org/10.1007/978-3-0348-0113-3.
- [3] Andreas Kirsch and Natalia Grinberg. The Factorization Method for Inverse Problems, volume 36 of Oxford Lecture Series in Mathematics and its Applications. Oxford University Press, 2008. URL: https://doi.org/10.1093/acprof:oso/9780199213535.001.0001.

## See Also

- factorizationGuide.html

## Code: function: factMeth

    function z = factMeth(seti)

Load data F (farField)

Notation: F: data (far field)

    farField = seti.FmeasDelta; % requires far field data (as factorization method is restricted to far field)

## Code: function: factMeth: Construct F_# and Compute the SVD

Construct F_# = `Re(F)` + Im(F) and compute its SVD.

Construct F_#: 1st step: Re(F) and Im(F)

Note that the real and imaginary parts of a bounded linear operator on a Hilbert space over are the self-adjoint operators, see [3],

These definitions are in accordance with the corresponding definition for complex numbers, see [1]. That does mean that

This does not mean that in the numerical implementation the real part only contains real numbers and the imaginary part only imaginary ones for a matrix .

    realPart = (farField + farField')/2;      % Re(F)
    imagPart = (farField - farField')/(2*1i); % Im(F)

Construct F_#: 2nd step: `Re(F)`

    [V,D] = eig(realPart);
    absRealPart = V*abs(D)*V'; % |Re(F)|

Construct: F_#: 3rd step: F_#

    farFieldSharp = absRealPart + imagPart; % F_# = |Re(F)| + Im(F)

Compute SVD of F_#

    [V,D,~] = svd(farFieldSharp); % V will be used for \psi

## Code: function: factMeth: Numerical Implementation of Factorization Method

The factorization method decides whether a point belongs to the image or not. In theory, we have (4.32) in Th. 4.2.6 in [2]. As we have to deal with a finite number of receivers and perturbed data we follow (4.33) with a regularization due to noise, i.e.

Preparations

    k = seti.k; % wave number

    s = zeros(size(seti.G(seti.qROIexact))); % ROI as matrix

    directions = seti.incPnts;
    % The directions are the coordinates of transmitters as well as receivers.
    % Actually it does not matter because they have the same positions.
    % Anyway, the same number of transmitters and receivers is needed (incNb = measNb).

    % Use grid points from ROI
    x1ROI = seti.gridROI(2,1:seti.nROI);
    gridPnts = x1ROI;

Numercial Implementation of Series of Factorization Method

    reg = seti.facRegFac*seti.delta; % regularization parameter
    Neig = floor(seti.measNb*seti.facNeigFac); % numbers of eigenvalues used in series before truncated
    if Neig > seti.measNb
        error('Neig is greater than seti.measNb. Set seti.facNeigFac not bigger than 1.');
    end

    for l = 1:length(gridPnts)
        for m = 1:length(gridPnts)
            testfunc = exp(-1i*k* (directions(1,:)*gridPnts(l)+directions(2,:)*gridPnts(m))); % \Phi_\infty
            % requires: incNb = measNb (e.g. 64)
            % testfunc: (e.g. size(testfunc) is 1 x 64)
            % V: eigenvectors in cols (e.g. size(V) is 64 x 64)
            % Use less eigenvalues and eigenvectors:
            %   Reduce number of eigenvalues and eigenvectors.
            D = D(1:Neig,1:Neig);
            V = V(:,1:Neig);
            % Series of factorization method:
            s(l,m) = 1./sum(diag(D)./(diag(D)+reg).^2.*transpose(abs(testfunc*V)).^2); % (4.33) in [2]
        end
    end
    z = s;

Rotate the result Rotate the result such that it fits (probably because of the grid)

    z = transpose(z); %
    z = rot90(rot90((z))); % rotate 180 degrees

    end