//  Reconstruction d'images par SVD
//  Utilisation de la toolbox IPD
//  Images disponibles en noir et blanc
//  boat.png, cameraman.png, peppers-bw, shepp256.png,
//  Sinogram1.png, taxi3.png

clear();
stacksize('max');

//  Lecture image
    titre = 'cameraman';
    nFic = titre+'.png';
    a = ReadImage(nFic);
    
//  Affichage image
    figure();
    nuPlan = 2 
    if (size(a,3) < 2) then nuPlan = 1;  end; 
    ShowImage(a(:,:,nuPlan),titre);

//  Calcul de la svd
    a = double(a);  //  a contient des valeurs réeelles
    [u,s,v] = svd(a(:,:,nuPlan));
    
//  Plot des valeurs singulières
    lp = min(size(s));    
    figure();
    plot2d([1:lp]',log10(diag(s)));
    xgrid(2);
    xtitle('Graphique des valeurs singulieres. Axe y logarithmique')
    
//  Reconstruction approximative de l'image     
    vs = [40, 60, 80, 100, 120, 140, 160];    //  p = nombre valeurs singulières utilisées
    for i = 1:size(vs,2);
        q = zeros(size(s));
        p = vs(i);
        q = u(:,1:p)*s(1:p,1:p)*v(:,1:p)';
        ps = string(p);
        lbl = titre+' p = '+ps;
        figure();
        ShowImage(uint8(q),lbl);
    end;
    
