0
|
9570 Views
|
2 Replies
|
1 Total Likes
View groups...
Share
GROUPS:

# fitting the Airy diffraction pattern with FindFit

Posted 9 years ago
 Hello! I have 2D data (images) which are results of Gaussian type beam diffraction. According to the theory of diffraction the pattern can be discribed by function ~ J1(r)/r. I tried to fit my data with this function using the FIndFit routine, but had an error : Blockquote Power::infy: "Infinite expression 1/0. encountered" Blockquote data01 = Import["366.0000cm_660.0000nm.png"]; data02 = ImageData[data01]; min02 = Min[data02]; max02 = Max[data02]; (*Dimensionalities of images*) dim01a = Import["366.0000cm_660.0000nm.png", "ImageSize"] []; dim01b = Import["366.0000cm_660.0000nm.png", "ImageSize"] []; dim01 = Min[dim01a, dim01b]; ReducPix = 5; data03 = Table[ (data02[[i]][[j]][] - min02)/max02, {i, 1, dim01b, ReducPix}, {j, 1, dim01a, ReducPix}]; (*Dimensionalities*) dim03a = Dimensions[data03][]; dim03b = Dimensions[data03][]; dim03 = Min[dim03a, dim03b] ListPlot3D[data03, PlotRange -> All] s[x_, y_] := Sqrt[a2 (x - a1)^2 + b2 (y - b1)^2] Bessel3D = Piecewise[ {{1/4, x == 0 && y == 0}, {(a0 BesselJ[1, s[x, 0]]/s[x, 0]), y == 0}, {(a0 BesselJ[1, s[0, y]]/s[0, y]), x == 0}}, {(a0 BesselJ[1, s[x, y]]/s[x, y])}] data3D = Flatten[MapIndexed[{#2[], #2[], #1} &, data03, {2}], 1]; fit = NonlinearModelFit[ data3D, (Bessel3D)^2, {{a0, 1}, {a1, 50}, {a2, 10}, {b1, 50}, {b2, 10}}, {x, y}];  I tried to use different method in the function (NMinimize, Gradient, so on ) but it did not help. It is really puzzling for me. Will appreciate any suggestions!!
2 Replies
Sort By:
Posted 9 years ago
 Hi Alex,This is an interesting problem. The airy function has a singularity at 0 which can be removed by defining the value there as 1, which is the limit. I do that in the code below, and also make use of the radial symmetry to fit for only the central amplitude and the center. (Which is not all that is needed.)I have found that the code below works intermittently. It seems still dependent on choices of values made during the fit function, even though I think I have removed the singularity. I Hopefully someone with more knowledge will speak up. In:= (2 BesselJ[1, r]/r)^2 /. r -> 0 During evaluation of In:= Power::infy: Infinite expression 1/0^2 encountered. >> During evaluation of In:= Infinity::indet: Indeterminate expression 0 ComplexInfinity encountered. >> Out= Indeterminate In:= Limit[(2 BesselJ[1, r]/r)^2, r -> 0] Out= 1 In:= airy[r_] := (2 BesselJ[1, r]/r)^2 /; r != 0. In:= airy[0.] = 1 Out= 1 In:= airyData = Flatten[Table[{x, y, 5 airy[Sqrt[(x - 1)^2 + (y - 2)^2]]}, {x, -10, 10, .1}, {y, -10, 10, .1}], 1]; In:= ListPlot3D[airyData]; In:= NonlinearModelFit[airyData, a airy[Sqrt[(x - h)^2 + (y - k)^2]], {a, h, k}, {x, y}] 
Posted 9 years ago
 Hi David! thanks for your interest to this problem! You are right, there is the singularity at zero point. I tried to eliminate it using a Picewise definiton of the model function in my code but for some reason it did no help and I still had the error. Looks like your definition of the functionairy[r_] := (2 BesselJ[1, r]/r)^2 /; r != 0.airy[0.] = 1;is able to overcome that issue. Great! An addition specification of the range for parameters in NonlinearModelFit helps as well. However when I'm trying to fit my data the result is not accurate. A new error pops up: NonlinearModelFit::nrlnum: The function value {2.3360210^119-3.4376510^103 I,3.2420410^118+0. I,4.6252810^117+6.1297110^101 I,6.7870110^116-1.18210^101 I,<<44>>,1.1806110^96+0. I,1.15086*10^96+0. I,<<12238>>} is not a list of real numbers with dimensions {12288} at {a0,a1,a2,b1,b2} = {0.206555,69.9985,-2.9634,50.0133,-2.99061}. >> not sure what does it mean.. The 3D plots with parameters from the fitting are not very similiar to the original data at all and as you mentioned very sensitive to the initial values (guess) in the fitting. Thank you anyway!Thanks for your suggestions!