function seti = setWavelet(seti)
% set wavelet
% uses nInv (and not nROI)
% if gscale is not set or 0, then nInv = nROI and nCDinv = nCD
% is set in checkConsistency

if seti.useWavelet == 1

    reshapeVec = seti.nROI*ones(1,seti.dim);
    G  = @(x) reshape(x,reshapeVec); %  G: n x 1  -> matrix
    iG = @(x) x(:);                  % iG: matrix -> n x 1

    seti = Wfuncs(seti,G,iG);

    if seti.wavWeight == 1
        w = whiteNoiseCalc(seti,seti.W,G,iG,seti.nCD,seti.nROI); %weight VECTOR for ROI
        w = normsApproxIsom(w,seti.W);
    else
        w = 1; % no weight
    end
    seti.wavWght = w;

    % wavelet transforms with weight (then seti.wavWght is not needed)
    seti.wW = @(x) w.*seti.W(x);
    seti.wiW = @(x) seti.iW(1./w.*x);
    seti.wWstar = @(x) seti.Wstar(w.*x);
    seti.wiWstar = @(x) 1./w.*seti.iWstar(x);

    % norms approximateley isometric?
    % norm(randq,1)
    % norm(w.*W(randq),1)

    if isfield(seti,'wavIsom')
        if strcmp(seti.wavIsom, 'W1p')
            % construct weights for W^{1,p}-norm
            if seti.gscale == 1
                %seti.omegaW1pInv = normsApproxW1Isom(seti.nCDinv, seti);
                seti.omegaW1p = normsApproxW1Isom(seti.nInv, seti);
                seti.omegaW1p = seti.omegaW1p(:);
            else
                seti.omegaW1p = normsApproxW1Isom(seti.nROI, seti);
                seti.omegaW1p = seti.omegaW1p(:);
            end
            seti.omegaW1p = seti.omegaW1p/abs(min(seti.omegaW1p));
        else
            error('Error in setWavelet.m: seti.wavIsom not implemented.')
        end
    end

    %figure(102);
    %imagesc(G(iW(w))); colormap(seti.colormap);
    %axis xy;
    %error('stop')
end

if seti.useWavelet == 1 % && seti.gscale == 1
% similar to if before... for nInv...
% definitions e.g. Gi instead of G, iGi instead of iG etc.

    reshapeVeci = seti.nInv*ones(1,seti.dim);

    Gi  = @(x) reshape(x,reshapeVeci); %  Gi: n x 1  -> matrix
    iGi = @(x) x(:);                  % iGi: matrix -> n x 1

    seti = Wifuncs(seti,Gi,iGi);

    if seti.wavWeight == 1
        wi = whiteNoiseCalc(seti,seti.Wi,Gi,iGi,seti.nCDinv,seti.nInv); %weight VECTOR for down scaled ROI
        wi = normsApproxIsom(wi,seti.Wi);
    else
        wi = 1; % no weight
    end
    seti.wavWghti = wi;

    % wavelet transforms with weight (then seti.wavWghti is not needed)
    seti.wWi = @(x) wi.*seti.Wi(x);
    seti.wiWi = @(x) seti.iWi(1./wi.*x);
    seti.wWstari = @(x) seti.Wstari(wi.*x);
    seti.wiWstari = @(x) 1./wi.*seti.iWstari(x);

end

end

function seti = Wfuncs(seti,G,iG)
seti.W      = @(x) iG(ldwt(G(x),'forward',seti.wavelet,seti.extension,seti.smin));
seti.iW     = @(x) iG(ldwt(G(x),'inverse',seti.wavelet,seti.extension,seti.smin));
seti.Wstar  = @(x) iG(ldwt(G(x),'adjoint',seti.wavelet,seti.extension,seti.smin));
seti.iWstar = @(x) iG(ldwt(G(x),'adjointinverse',seti.wavelet,seti.extension,seti.smin));
end

function seti = Wifuncs(seti,Gi,iGi)
seti.Wi      = @(x) iGi(ldwt(Gi(x),'forward',seti.wavelet,seti.extension,seti.smin));
seti.iWi     = @(x) iGi(ldwt(Gi(x),'inverse',seti.wavelet,seti.extension,seti.smin));
seti.Wstari  = @(x) iGi(ldwt(Gi(x),'adjoint',seti.wavelet,seti.extension,seti.smin));
seti.iWstari = @(x) iGi(ldwt(Gi(x),'adjointinverse',seti.wavelet,seti.extension,seti.smin));
end

function w = normsApproxIsom(w,W)
% -- make norms approximateley isometric
randq = rand(size(w))+1i*rand(size(w));
isom = norm(randq,1)/norm(w.*W(randq),1);
w = isom.*w;
end

function weights = normsApproxW1Isom(n, seti)
% compute weights for approximate W^{1,p*} norm with p* = 3/2
% works only in 2D !!!

wavMat = zeros(n,n);
weights = wavMat;
pAst = 3/2;

for jj=1:n
    for ll=1:n
        wavMat(jj,ll)=1;

        oriMat = ldwt(wavMat, 'inverse',seti.wavelet,seti.extension,seti.smin);
        [oriMatX,oriMatY] = gradient(oriMat,seti.rCD/sqrt(2)/n,seti.rCD/sqrt(2)/n);

        normOri = ( normroi(oriMat(:),seti,pAst)^pAst ...
                  + normroi(oriMatX(:),seti,pAst)^pAst ...
                  + normroi(oriMatY(:),seti,pAst)^pAst )^(1/pAst);

        weights(jj,ll) = normOri/normroi(wavMat,seti,1);

        wavMat(jj,ll)=0;
    end
end

% figure(101)
% surf(w); colorbar;

end