clear;
close all;
skellam = Skellam(1024,0,32);
k=1;
xq=linspace(-20,20,1024);
for h=4:0.03:10
    j=1;
    for x=-20:0.7:20;
        suma=0;
        for i=1:length(skellam)
            if ((((skellam(i)/h)-floor(x/h))>=0) && (((skellam(i)/h)-floor(x/h))<1))
               suma=suma+1/h;
            end
        end
        wynik1(k,j)=suma/length(skellam);
        j=j+1;
    end
    j=1;
    for x=-20:0.7:20;
        suma=0;
        for i=1:length(skellam)
            if ((((skellam(i)/h)-(x/h))>=0) && (((skellam(i)/h)-(x/h))<1))
               suma=suma+1/h;
            end
        end
        wynik2(k,j)=suma/length(skellam);
        j=j+1;
    end
    suma1=0;
    for l=1:length(xq)
        suma1 = suma1+(wynik1(k,floor(xq(i)))-f(xq(i)))^2;
    end
    val1(k)=suma1/2*length(xq);
    suma2=0;
    for l=1:length(xq)
        suma2 = suma2+(wynik2(k,floor(xq(i)))-f(xq(i)))^2;
    end
    val2(k)=suma1/2*length(xq);
    k=k+1;
end
[a1 b2] = min(val1)
[a1 b2] = min(val2)