This function computes the Hankel transform defined as integral from 0 to inf f(x, dstruct) J sub order (x r) dx, where J sub order is a Bessel function of the first kind and order is either 0 or 1. The variable dstruct is a structure to pass data to f the (possibly) complex user provided function f. The variables rerr and aerr are the relative and absolute error passed to the Matlab function integral.
The algorithm used is integration between the zero crossings of the Bessel function to obtain partial integrals, followed by their summation using the Padé approximants. This approach is generally more accurate than digital filter algorithms. It also is able to handle well defined but divergent integrals such as f(x) = x.
Alan Chave (2021). Hankel(order,r,f,dstruct,rerr,aerr) (https://www.mathworks.com/matlabcentral/fileexchange/65797-hankel-order-r-f-dstruct-rerr-aerr), MATLAB Central File Exchange. Retrieved .
This code is intended for operation on integrands (outside the Bessel function) that are fairly smooth and without singularities. This occurs in the limit where displacement current is negligible compared to conduction current, which is the situation in electromagnetic geophysics. It works well (and substantially better than digital filter approaches) in that case..
Thanks Prof. Alan. I have tested your function with my Sommerfeld Integrals and I was not satisfied with the accuracy and speed too. But, I will investigate other methods from your book if they will suit my application. Anyway it is a great method and code, and is definitely helpful in other applications.
The method is not based on a filter, but rather direct integration between the zero crossings of the bessel function followed by summation of the ensuing alternating series using Padé approximants to speed up convergence. The original Fortran code is described in Geophysics, 48, 1671-1686, 1983. Other applications of the Padé approximants are described in my book Computational Statistics in the Earth Sciences, Cambridge U Press, 2017.
Is this an FHT filter based? Are there limitations for this method when r>>1 or r<<1?
Currently, I am using the weighted sums along the real axis method for computing the Hankel transforms appearing in Sommerfeld Integrals.
I would like to know, which paper did you use for this implementation?
Find the treasures in MATLAB Central and discover how the community can help you!Start Hunting!