function y = quadp(hankel,a,b,rho)

y = 2*pi*quad(@nested_func,a,b);
    function z = nested_func(r)
        z = hankel(r,rho);
    end
end