function [ker_refined] = single_level(dir,fn, im, gamma, segind, blursize)
    
im = im.^gamma;

dsp_ratio = 1/sqrt(2);
dsp_level = ceil( log(max(blursize)/3) / log(1/dsp_ratio) ) + 1;
dsp_pry = dsp_ratio.^(dsp_level-1:-1:0);    

dsp_level = numel(dsp_pry);
mu = 5e-2;
mu_mulfac = (5e-3/mu)^(1/(dsp_level-1));

for ids = 1:dsp_level

    maxiter = 10;

    fprintf('\nlevel:%2d\n', ids);

    imB = imresize(im,dsp_pry(ids));
    [Bx1,By1] = gradient(imB,1);
    [h,w] = size(imB);

    ksize = round(blursize*dsp_pry(ids));
    ksize = ksize + (1-mod(ksize,2));

    %% main
    if ids == 1
        % init kernel image
        ker = zeros(ksize,ksize);
        [Bx0,By0] = gradient(im,1);
        if sum(sum(abs(Bx0))) < sum(sum(abs(By0)))
            ker(2,:) = [1,1,1]/3;
        elseif sum(sum(abs(Bx0))) > sum(sum(abs(By0)))
            ker(:,2) = [1;1;1]/3;
        else
            ker(2,2) = 1;
        end
        % init latent image
        imL_est = imB;
        imL_est(imL_est>1) = 1;
        imL_est(imL_est<0) = 0;
    else
        % upsample kernel image
        ker = imresize(ker, [ksize,ksize]);
        % ker(ker>1) = 1; % uncomment on 22-Oct-2014
        ker(ker<0) = 0;
        ker = ker./sum(sum(ker));
        % upsample latent image
        imL_est = imresize(imL_est, size(imB));
        imL_est(imL_est>1) = 1;
        imL_est(imL_est<0) = 0;
    end
    tol = 1e-4;
    thres_r = 0.8;
    for iter = 1:maxiter

        %% latent image update
        imL_est_old = imL_est;
        imL_est = admmfft(imB,ker,mu,1);
        imL_est(imL_est>1) = 1; imL_est(imL_est<0) = 0;

        %% kernel im/params update
        [Px,Py] = gradient(imL_est,1);

        br = 0.05;
        Px = ignore_border(Px,br);
        Py = ignore_border(Py,br);
        Bx = ignore_border(Bx1,br);
        By = ignore_border(By1,br);

        % rmap Oct-20-2014
        [Px,Py,rmap] = filter_rmap(Bx,By,Px,Py,ksize(1),thres_r);
        thres_r = thres_r/1.1;
        % imwrite(rmap,sprintf('%s/%d_rmap_scale_%d_iter_%d.png',filename,segind,ids,iter));


        ker = ifft2( ( conj(fft2(Px)).*fft2(Bx) + conj(fft2(Py)).*fft2(By) ) ./ ( conj(fft2(Px)).*fft2(Px) + conj(fft2(Py)).*fft2(Py) + 10 ) ) ;
        ker = fftshift(ker);
        center = ceil([(h+0.1)/2 (w+0.1)/2]);
        ker = ker(center(1)-floor(ksize/2):center(1)+floor(ksize/2),...
            center(2)-floor(ksize/2):center(2)+floor(ksize/2));
        % ker(ker>1) = 1; % uncomment on 22-Oct-2014
        ker(ker<0) = 0;
        ker = ker./sum(sum(ker));


        %% termination
        err = sqrt(sum(sum((imL_est-imL_est_old).^2)))/size(imL_est,1)/size(imL_est,2);
        if err<tol
            break;
        end

    end

    % at the end of each level, update weight of tv prior
    mu = mu*mu_mulfac;

    kerout = uint8( (ker-min(min(ker))) / (max(max(ker))-min(min(ker))) *256);
    imwrite(imresize(ind2rgb(kerout,gray(256)),[blursize,blursize],'nearest'), sprintf('%s/%s_seg_%d_psf_scale_%d.png',dir,fn,segind,ids));
    imwrite(imL_est.^(1/gamma), sprintf('%s/%s_seg_%d_imL_est_scale_%d.png',dir,fn,segind,ids));
end

% refine using ISD
ker_refined = ker;%kernel_refinement_isd(imL_est, imB, ker);
    

