# fresnelSynthNuMaxVarPublic

Compares Institute Fresnel's data set with synthetic data generated by the exact contrast for several approximation orders nuMax. This is useful to find an appropriate approximation order nuMax for the incident field's matching.

Warning: This function clears variables and creates a directory. Be careful when using.

## Contents

- Syntax
- Description
- More About
- References
- See Also
- Code

## Syntax

Start this routine in the folder `code` by

    init;
    fresnelSynthPaper

## Description

This routine is similar to `fresnelSynthPaper`. Therefore we refer to fresnelSynthPaper.html for a description.

In contrast to the routine `fresnelSynthPaper` computations with various values for the approximation order `nuMax` are supported. Therefore the parameter `nuMaxMax` can be set.

In addition, resulting figures are saved in `output/fresnelSynthNuMaxVarPublic`.

Note that a suitable choice of `nuMax` is essentially for a suitable matching of the incident field.

## More About

This file is essentially a copy of the intern routine `fresnelSynthNuMaxVar` with adaptions and documentation. It was used to create Fig. 3.9 in [1]. (The routine `fresnelSynthNuMaxVar` originated from `fresnelSynth`. It is a small part, that was adapted.)

## References

- [1] Florian Bürgel. Effective and Efficient Reconstruction Schemes for the Inverse Medium Problem in Scattering. PhD thesis, Universität Bremen, 2019.

## See Also

- fresnelSynthPaper.html

## Code

Initialization

    close all;
    % clearvars -except inseti htext;
    clear all;

    disp(' ')
    disp('-----')
    disp('fresnelSynth')
    disp(' ')

    % -- initialize
    tic;
    init;

    calcyExact = 0;

    % Make dir to save figures
    dirOutput = 'output/fresnelSynthNuMaxVarPublic/';
    if exist(dirOutput,'dir')~=7
        fprintf('Make dir: %s.\n',dirOutput)
        mkdir(dirOutput);
    end

    % -- set parameters

    % Two possibilities:
    % * Use "-- load paperFresnelSynth" from start to end.
    % or
    % * Use "-- Set parameters from paperFresnelSynth (is a copy)" from start to end.

    % -- load paperFresnelSynth: start --
    % if ~exist('inseti','var')
    %     %inseti = 'setiContrastFresnelSynth'; % synthetic fresnel data
    %     inseti = 'paperFresnelSynth'; % synthetic fresnel data
    % end
    % seti.inseti = inseti;
    % eval(inseti)
    % -- load paperFresnelSynth: end --

    % -- Set parameters from paperFresnelSynth (is a copy): start --

    % -- setGrid
    seti.dim = 2;
    seti.rCD = 0.2;
    seti.nCD = 256; % original

    % -- fresnel data (loadData.m)
    seti.expData = 'fresnel';
    seti.fresnelFreq = 5*1E9; % standard 3 GHz (4 GHz)
    seti.fresnelFile = 'inexpdata/fresnel_opus_1/dielTM_dec8f.exp'; % 1--8 GHz, test: single cylinder, 4 GHz
    %seti.fresnelFile = 'data/fresnel_opus_1/twodielTM_8f.exp'; % 1--8 GHz
    filename = seti.fresnelFile;

    seti.nuMax = 10; seti.ampCalc = 1; % default (B1)

    % -- setKernel
    cLight = 2.99792458E8; % light velocity
    seti.k = 2*pi*seti.fresnelFreq/cLight;
    seti.model = 'helmholtz';

    % -- setMeas (measurements)
    % simulate fresnel_opus_1:
    seti.incType = 'pointSource'; % used in op. 1 and 2 (correct? should be...)
    seti.measType = 'nearField'; % used in op. 1 and 2
    seti.incNb = 36;
    seti.measNb = 72;
    seti.radSrc = 0.72;
    seti.radMeas = 0.76;

    % -- inversion process (regularization)
    %seti.invNo = 3;
    seti.invNo = 6;
    seti.delta = 0; % to test

    %seti = setiAuxRec(seti,usefindpdatau,usefindalpha);
    seti.alpha = 1E8; % old: 2E4 corresponding 0.05 (in old code)
    seti.beta = 0; % because we expect no sharp edges...
    seti.physBounds = [-500,500,0,inf];
    seti.nOut = 30;

    % -- minimization with pda
    seti.pdaN = 50;

    % -- Set parameters from paperFresnelSynth (is a copy): end --

Load Data from Institute Fresnel

    seti = checkConsisExpData(seti);

    for freq = 1:8 % INPUT: frequency in GHz

        seti.fresnelFreq = freq*1E9;

        nuMaxMax = 20; % INPUT: nuMaxMax
        % nuMaxMax = 100; % INPUT: nuMaxMax

        errCmean = zeros(nuMaxMax,1);
        errRe = zeros(nuMaxMax,1);
        errIm = zeros(nuMaxMax,1);

        for nuMax = 1:nuMaxMax

            disp('----------------------------------')
            fprintf('nuMax = %g\n',nuMax);
            disp('----------------------------------')

            seti.nuMax = nuMax;
            disp(['loading fresnel data from "' seti.fresnelFile '" for frequency = ' num2str(seti.fresnelFreq/1E9) ' GHz']) ;
            seti = loadData(seti);
            %seti = evalc('loadData(seti)'); % setGeomSim is called in loadData
            %evalc: computes but no disp output...
            %seti.errC is inside...

            yFresnel = seti.FmeasDelta; % is loaded Fresnel data

            %disp('- yExactInc calculation using incField from Fresnel data')
            %disp('  but incField is NOT exact, so yExactInc might not be Fresnel data')
            FmeasExactInc = mimo(seti, seti.qROIexact, 'simo');
            FmeasExactInc = fresnelMiss(FmeasExactInc,yFresnel);

Compare yFresnel and yExactInc

Compare real-world data from Institute Fresnel, yFresnel, with exact computed data, yExactInc.

            disp('Compare yFresnel and yExactInc')

            errCmean(nuMax) = mean(seti.errC);
            errRe(nuMax) = norm(real(yFresnel-FmeasExactInc))/norm(real(FmeasExactInc));
            errIm(nuMax) = norm(imag(yFresnel-FmeasExactInc))/norm(imag(FmeasExactInc));

            %fprintf('nuMax = %g | ',nuMax);
            %fprintf('rel. err (real part) = %g | ',errRe);
            %fprintf('rel. err (imag part) = %g\n',errIm);

        end

Notation of Values

- `Cmean`: rel. err. of (mean because min = max in our experience)
- `errRe`: rel. err of data y (real part)
- `errIm`: rel. err of data y (imag part)

        disp('----------')
        disp(' nuMax | Cmean in % | errRe in % | errIm in %'); % with *100 below... then %...
        for nuMax = 1:nuMaxMax
            fprintf('%d & %4.1f & %4.1f & %4.1f\n',nuMax,errCmean(nuMax)*100, errRe(nuMax)*100, errIm(nuMax)*100);
        end
        disp('----------')

Figures

        figFontSize = 20; % default in ipscatt publish: 20 (pubFontSize) (does influence the result; tested with 40)
    %    lw = 2; % LineWidth (default for printing; is very thin in PhD thesis)
        lw = 3; % LineWidth (default for talks)

        figure(101);
    %     h = plot(t,step,'b-',t,stepinv,'r:',t,t,'k-.');
    %     h = plot(1:nuMaxMax, errCmean, 1:nuMaxMax, errRe, 1:nuMaxMax, errIm);
        h = plot(1:nuMaxMax,errCmean,'k-.',...
            1:nuMaxMax,errRe,'r-',...
            1:nuMaxMax, errIm,'b:');
        % -. dash-dot line; - solid line, : dotted line

        set(h(1),'LineWidth',lw);
        set(h(2),'LineWidth',lw);
        set(h(3),'LineWidth',lw);

        axis square;
    %     xlim([-5 5]);
        ylim([0 2]);
        set(gca, 'FontSize', figFontSize)
    %     set(gca, 'XTick', [a 0 b]);
    %     set(gca, 'YTick', [a 0 b]);
    %     set(gca, 'XTickLabelMode', 'manual', 'XTickLabel', {'a','0','b'});
    %     set(gca, 'YTickLabelMode', 'manual', 'YTickLabel', {'a','0','b'});

        % legend('mean(errC)','data errRe','data errIm');
        % xlabel('nuMax');
        % title('suitable amplitudes: variation of nuMax')

Save Figures as png, eps and fig

        seti.dirDatetime = datestr(now, 'yyyymmddTHHMMSS');
        % additional information
        [~,name,~] = fileparts(seti.fresnelFile);
        filePng = sprintf('%s%s_%d_GHz_%s.png',dirOutput,name,seti.fresnelFreq/1E9,seti.dirDatetime);
        figNo = 101;
        figNostr = sprintf('-f%d',figNo);
        print(figNostr,'-dpng',filePng);
        fileEpsc = strrep(filePng,'.png','.eps');
        print(figNostr,'-depsc',fileEpsc);
        fileFig = strrep(filePng,'.png','.fig');
        savefig(figNo,fileFig);

    end