DECLARE FUNCTION j1! (p!) pi = 4 * ATN(1!) 20 REM function j1 calculates j1(p) 30 REM the bessel function of 1st kind 1st order 40 REM series taken from crc handbook of math p445 50 REM g e copeland 60 REM - updated 2/2/2001 REM updated again improved bessel function evaluations REM 3/26/2001 REM GE Copeland REM Dept of Physics REM Old Dominion University REM Norfolk Va 23504 REM 200 DIM I(300), t(300) 201 PRINT CHR$(12) 210 PRINT "Circular Aperature---Fraunhofer diffraction pattern" 220 PRINT 230 PRINT " A circular opening in a solid plane be the aperature" 240 PRINT "through which parallel light is incident. The radius of the " 250 PRINT "aperature is r . Light of wavelength,l,is brought" 260 PRINT "to a focus on a screen located a distance f away from the" 270 PRINT "aperature." 280 PRINT 290 PRINT "We wish to find the intensity pattern on that screen. Let " 300 PRINT " p=k*r*sin(0)," 310 PRINT " where k is the wavevector = 2 *pi /wavelength." 322 PRINT " and 0 is the angle of observation on the screen as" 324 PRINT "measured from the geometrical center of the aperature" 326 PRINT "as projected along the optical axis." 330 PRINT "Then the intensity is given by:" 340 PRINT 350 PRINT , "I (p) = i0 *[ 2*J1(p)/p ]^2" 360 PRINT "where J1 is the Bessel function of 1st kind 1st order." 370 PRINT 380 PRINT "Now let's calculate the pattern for a case of concern." 390 PRINT "Let's keep the wavelengths between; 100 and 1000 nm.; "; "" 400 PRINT "Input the wavelength in nanometers"; 410 INPUT l 420 IF l < 100 THEN GOTO 390 430 IF l > 1000 THEN GOTO 390 440 REM 450 REM use wave lengths that are reasonaable 460 REM use optical wavelengths 470 l = l * 1E-09 480 490 k = 2 * pi / l 500 PRINT "This corresponds to a wavevector of k ="; k; " m-1" 510 PRINT 520 PRINT "Input the radius of the aperature "; 530 INPUT rmm 540 r = rmm * 1000! 550 PRINT 560 n1 = 0 570 REM counter of the values 580 t1 = 1.22 * l / (2 * r) 590 PRINT "The angular radius of the first ring ="; PRINT USING "##.###############"; t1 PRINT " or "; t1 PRINT " radians" 600 PRINT " or "; PRINT USING "##.#############"; t1 * 180! / pi; PRINT " or"; t1 * 180 / pi; PRINT " degrees." 602 REM 604 REM n=1 order of bessel function 606 n = 1 610 FOR a = 0! TO 2.9 * t1 STEP t1 / 75! 620 n1 = n1 + 1 630 t(n1) = a 640 p1 = k * r * SIN(a) 650 IF p1 = 0 THEN GOTO 680 652 z = p1 654 p = j1(z) 656 I(n1) = (2! * p / p1) ^ 2 670 GOTO 700 680 REM case p1=0 690 I(n1) = 1! 700 REM 710 REM print n1,t(n1),i(n1) 720 NEXT a 721 PRINT "Do you want to see some results of the calculation?"; 722 INPUT s$ 723 IF s$ = "YES" OR s$ = "Yes" OR s$ = "yes" THEN 727 724 IF s$ = "NO" OR s$ = "No" OR s$ = "no" THEN 740 725 PRINT "Please Answer Yes or No." 726 GOTO 721 727 PRINT 728 PRINT "Angle", "Intensity", "Log Intensity" PRINT " Step #" 729 FOR I = 1 TO n1 STEP INT(n1 / 10) 730 PRINT USING " ####.##### "; t(I) / t1; I(I); LOG(I(I)) 731 NEXT I 732 PRINT 740 PRINT "Do you want to plot of the data?"; 750 INPUT a$ 760 IF a$ = "yes" OR a$ = "YES" OR a$ = "Yes" THEN 810 770 IF a$ = "no" OR a$ = "NO" OR a$ = "No" THEN 800 780 PRINT "Please input yes or no." 790 GOTO 740 800 STOP 810 REM print "Do you want plot of the diffraction pattern"; 950 REM graphics plot 960 OPEN "plttek.dat" FOR OUTPUT AS #2 970 PRINT #2, 1 980 PRINT #2, 1 990 PRINT #2, n1 1000 PRINT #2, " Relative intensity" 1010 PRINT #2, "Angular departure from Optical Axis." 1020 PRINT #2, "Circular diffraction radius ="; rmm; " mm" 1030 PRINT #2, "Wavelength ="; l * 1000000!; " nm" 1040 FOR I = 1 TO n1 1050 PRINT #2, t(I) / t1, I(I) 1060 NEXT I CLOSE #2 1070 PRINT "Now return and do the plot. We have chained to convert the" PRINT "plotting file." CHAIN "tek2wplt.bas" 1080 END FUNCTION j1 (p) REM use case where 0 3 in Bessel" GOTO 111 ELSE END IF x2 = x * x x4 = x2 * x2 x6 = x2 * x4 x8 = x4 * x4 x10 = x2 * x8 X12 = x6 * x6 REM a = .5# - .56249985# * x2 + .21093573# * x4 a = a - .03954289# * x6 + .00443319# * x8 a = a - .00031761# * x10 + .00001109# * X12 j1 = p * a EXIT FUNCTION 111 REM do the other series y = 3! / p y2 = y * y y3 = y * y2 y4 = y2 * y2 y5 = y * y4 y6 = y4 * y2 theta1 = p - 2.35619449# + .12499612# * y theta1 = theta1 + .0000565# * y2 - .00637879# * y3 theta1 = theta1 + .00074348# * y4 + .00079824# * y5 theta1 = theta1 - .00029166# * y6 f1 = .79788456# + .00000156# * y + .01659667# * y2 f1 = f1 + .00017105# * y3 - .00249511# * y4 + .00113653# * y5 - .00020033# * y6 j1 = f1 * COS(theta1) / SQR(p) EXIT FUNCTION 80 x1 = p * .5 90 x2 = x1 * x1 100 x3 = x2 * x1 110 x5 = x2 * x3 120 x7 = x2 * x5 130 x9 = x7 * x2 140 x0 = x1 - x3 / (2!) + x5 / (2 * 3 * 2) - x7 / (3 * 2 * 4 * 3 * 2) + x9 / (4 * 3 * 2 * 5 * 4 * 3 * 2) 150 x0 = x0 - (x9 * x2) / (6 * (5 * 4 * 3 * 2) ^ 2) 160 x0 = x0 + x9 * x2 * x2 / (7! * (6 * 5 * 4 * 3 * 2) ^ 2) 170 x0 = x0 - x9 * x5 * x1 / (8! * (7 * 6 * 5 * 4 * 3 * 2) ^ 2) j1 = x0 END FUNCTION