Title :PINT Keywords :EELS - planar interface calculations Computer :VAX 780 Operating System :VAX/VMS version 4.7 Programming Language :FORTRAN 77 Hardware Requirements :Tektronics 4010 series graphic terminal or emulator desirable but not essential. Author(s) :C.A. Walsh - based on a BASIC program written by A. Howie. Correspondence Address :Microstructural Physics Group, Cavendish Laboratory, Madingley Road, Cambridge, CB3 0HE, England. Abstract: The program uses the result obtained by Garcia-Molina, Gras-Marti, Howie and Ritchie (J. Phys. C18 (1985) pp5335) to calculate the energy loss of electrons travelling parallel to a planar interface. They obtained an expression for the differential probability dP/(dzdhf), for an electron travelling in the z direction and losing energy hf. The probability of losing energy between hf-dhf/2 and hf+dhf/2 after travelling a distance dz parallel to the interface is given by dP/(dzdhf) x dzdhf. The program also uses this result to calculate the energy loss of electrons incident at an angle to an interface (see Howie, Milne and Walls,(1985) Inst.Phys. Conf. Ser. No. 78: Chapter 4, EMAG'85). The program requires dielectric data for the media on either side of the interface. A limited amount of data is stored inside the program for some materials. Other data can be supplied to the program by data files. Title :PINT.DOC Keywords :EELS - planar interface calculations Computer :VAX 780 Operating System :VAX/VMS version 4.7 Programming Language :FORTRAN 77 Hardware Requirements :Tektronics 4010 series graphic terminal or emulator desirable but not essential. Author(s) :C.A. Walsh - based on a BASIC program written by A. Howie. Correspondence Address :Microstructural Physics Group, Cavendish Laboratory, Madingley Road, Cambridge, CB3 0HE, England. Documentation: The program uses the result obtained by Garcia-Molina, Gras-Marti, Howie and Ritchie (J. Phys. C18 (1985) pp5335) to calculate the energy loss of electrons travelling parallel to a planar interface. They obtained an expression for the differential probability dP/(dzdhf), for an electron travelling in the z direction and losing energy hf. The probability of losing energy between hf-dhf/2 and hf+dhf/2 after travelling a distance dz parallel to the interface is given by dP/(dzdhf) x dzdhf. In practise, for a specimen of thickness T, the calculations are split into x number of steps of length dz=T/x << plasmon mean free path. This allows multiple scattering to be taken into account more effectively. The same result can also be used to calculate the energy loss of electrons incident at an angle to an interface (see Howie, Milne and Walls,(1985) Inst. Phys. Conf. Ser. No. 78: Chapter 4, EMAG'85). The calculations are made by considering the electron path to consist of a series of steps parallel to the interface. Such calculations are generally made in stages: 1. Start with the electron a long way from the interface and bring it in closer to the interface. Further from the interface the probability of loss is smaller, hence fewer steps and longer step lengths can be used. 2. With the electron initially closer to the interface move the electron to the interface. In this stage a short step length should be used so that the distance between successive steps (in the direction perpendicular to the interface) is small (1 Angstrom? try for yourself!). 3. From the interface into the medium 1. If you wish to compare results with experimental data to find the best penetration depth, then you can keep stopping at different distances inside the medium. Printed and graphical output are produced after each stage. N.B. A factor of 2 is already included in the step length to account for losses incurred during the return path of the electron out of the medium again. (This is NOT true when the electrons travel parallel to the interface (THETA=0).). An example of this is given in the listings at the end of this file showing the questions asked by the program and answers given by the user. The dielectric data for aluminium that was used is also given and the ouput file showing the calculated spectra are also shown. No graphical output was obtained. IMPLEMENTATION OF THE PROGRAM. The program is written in FORTRAN 77. It calls graphics routines from the library of routines NGRAPH supplied by Nestor Zaluzec. N.B. _________________________ N.B. ------------------------- This program is written using output unit 6 to write to the terminal/screen. The graphics routines written by Nestor Zaluzec write to the terminal on unit 7. Hence it will be necessary to change all WRITE statements to the terminal in either this program and/or in the NGRAPH routines to whatever unit is used on your own computer for output to the terminal. ------------------------- END OF N.B. ------------------------- The following is the compilation and linking procedure used on a VAX system for this program: $ FORTRAN PINT $ FORTRAN/EXTEND NGRAPH $ LINK PINT,NGRAPH $ RUN PINT The program then assumes the existance of the file TERMIN.DAT in your filespace. This file is accessed by Nestor Zaluzec's graphics routines and gives details of the graphics terminal type being used and the baud rate. A copy of this file is given here: GRAPHICS TERMINAL TYPE:10 BAUD RATE FOR GRAPHICS: 9600 -------------------------------------------------------- TERMINAL DEFINITIONS: --------------------- 0 = TEKTRONICS 4010-1,4014-1 (WITH HARDWARE CURSOR) 1 = TEKTRONICS 4006 ,4010 (NO HARDWARE CURSOR) 2 = LEAR SIEGLER ADM-3,ADM-5 (WITH RETROGRAPHICS RG-512) 3 = PERITEK VCG-512 BIT MAP COLOR GRAPHICS -------------------------------------------------------- BAUD RATE RANGE FOR GRAPHICS (110 - 19200) -------------------------------------------------------- N.B. The first two lines are read in using a formatted read statement. Hence, do not change the number of characters on these two lines if you need to change the terminal type or baud rate. DIELECTRIC DATA The program requires dielectric data for the media on either side of the interface. A limited amount of data is stored inside the program for some materials. Calculations made using this data are then made at 1eV energy intervals which will correspond, experimentally, to using 1eV wide energy channels for the collectin of the spectra. Users can supply their own dielectric data for use by the program. The data must be supplied with a constant energy interval. The data must be in a file with the following format: First line: Exactly 24 characters as a title. Second line: No. of data items energy interval between succesive items (eV) Third line: Energy (eV) real part of epsilon imaginary part of epsilon Fourth line: etc. MISCELLANEOUS INFORMATION 1. The calculations are made using 10E6 incident electrons. These are placed either in the zero loss channel, or as a gaussian distribution about the zero loss channel with the specified full width half height = 2xstandard deviation. Spectra written to a file in the format suitable for being read by Nestor Zaluzec's NELS program include 100 channels corresponding to negative energy loss. 2. The angle of incidence requested by the program is the angle between the electron beam direction and the plane of the interface. If this angle is set to zero the electon beam travels parallel to the interface. 3. THE INNER POTENTIAL In the glancing angle configuration, electrons are refracted when they pass from one medium to another. For an incident angle i the refracted angle r is given by r x r = i x i + V/E where V is positive and is equal to the inner potential of medium 2 when the electrons pass from a vacuum (medium 1) into medium 2. If neither media are a vacuum the refraction of the electrons is ignored in this program. V generally lies between 5 and 30 eV. A default value of 18 eV is specified in this program. The inner potential is used to recalculate the angle the electron beam makes to the interface when it passes from medium 2 into medium 1. However, for this to be done one stage of the calculations must finish with the electron at the interface. 4. When specifying the distance of the electrons from the interface it should be noted that positive values lie in medium 2 and negative values lie in medium 1. 5. The units of the scattering vectors Ky: A new variable = Kyv/w is used in the program (where v is the electron velocity and h-cross x w is the energy loss.) The integration in the energy loss calculation is made using this variable. The required units in whuch RKYMAX and RKYMIN, the maximum and minimum collected scattering vectors, must be specified are those for which Kyv/w' where h-cross x w' = e, the electronic charge. In terms of angular scattering this becomes semi-angle = Ky/K = (Kyv/w')/(Kv/w') = (Kyv/w')/2E where E, v and K are the energy, velocity and wavevector of the incident electrons. Hence, RKYMAX = Semi-angle (rads) x 2 x E (eV). For 100keV electrons 200 units = 1mrad. 6. It should be noted that the integration occurs over Ky only. The integration over Kx (where x is the direction perpendicular to the interface) has been implicitly been taken to infinity. Consequently the values of RKYMAX and RKYMIN do not quite correspond to the upper and lower collection angles in an electron microscope. In particular, it should be noted that for the case of an electron travelling through a bulk material the probability of energy loss is proportional to Im(-1/epsilon)x ln(qv/w) (Landau and Lifschitz, 'Electrodynamics of Continuous Media'), where q is the upper cut-off wavenumber, whereas the expression used in this program yields the result that the energy loss probability is proportional to Im(-1/epsilon)ln(2qv/w). 7. The output from this program consists of the calculated spectrum as number of electrons in each energy channel. The zero loss intensity specified is the number if electrons left in the zero loss channel. The probability given is that calculated for the last step. The spectrum gives the final number of electrons in each energy loss channel. A SAMPLE SET OF DIELECTRIC DATA ..Aluminium............. 50 1.0 1.000 -0.1539D+03 0.3021D+02 2.000 -0.5424D+02 0.1950D+02 3.000 -0.2497D+02 0.5258D+01 4.000 -0.1390D+02 0.2203D+01 5.000 -0.8617D+01 0.1120D+01 6.000 -0.5700D+01 0.6208D+00 7.000 -0.3923D+01 0.3753D+00 8.000 -0.2762D+01 0.2381D+00 9.000 -0.1962D+01 0.1560D+00 10.000 -0.1386D+01 0.1042D+00 11.000 -0.9562D+00 0.6971D-01 12.000 -0.6241D+00 0.5190D-01 13.000 -0.3701D+00 0.4588D-01 14.000 -0.1712D+00 0.4011D-01 15.000 -0.7825D-02 0.3834D-01 16.000 0.1195D+00 0.4174D-01 17.000 0.2232D+00 0.4013D-01 18.000 0.3103D+00 0.3885D-01 19.000 0.3835D+00 0.3745D-01 20.000 0.4460D+00 0.3582D-01 21.000 0.4875D+00 0.3438D-01 22.000 0.5290D+00 0.3294D-01 23.000 0.5706D+00 0.3151D-01 24.000 0.6121D+00 0.3007D-01 25.000 0.6536D+00 0.2863D-01 26.000 0.6762D+00 0.2730D-01 27.000 0.6988D+00 0.2597D-01 28.000 0.7213D+00 0.2463D-01 29.000 0.7439D+00 0.2330D-01 30.000 0.7665D+00 0.2197D-01 31.000 0.7806D+00 0.2117D-01 32.000 0.7946D+00 0.2037D-01 33.000 0.8087D+00 0.1957D-01 34.000 0.8227D+00 0.1877D-01 35.000 0.8368D+00 0.1797D-01 36.000 0.8460D+00 0.1744D-01 37.000 0.8552D+00 0.1691D-01 38.000 0.8643D+00 0.1639D-01 39.000 0.8735D+00 0.1586D-01 40.000 0.8827D+00 0.1533D-01 41.000 0.8892D+00 0.1487D-01 42.000 0.8957D+00 0.1441D-01 43.000 0.9023D+00 0.1396D-01 44.000 0.9088D+00 0.1350D-01 45.000 0.9153D+00 0.1304D-01 46.000 0.9200D+00 0.1271D-01 47.000 0.9247D+00 0.1238D-01 48.000 0.9295D+00 0.1204D-01 49.000 0.9342D+00 0.1171D-01 50.000 0.9389D+00 0.1138D-01 A SAMPLE RUN FOR THE REFLECTION CASE: GLANCING ANGLE PROGRAM FOR VACUUM 0 (MEDIUM 2 ONLY) MgO 1 SILICON 2 SILICA 3 CARBON 4 DIAMOND 5 ALUMINIUM 6 ALUMINA 7 GaAs 8 GaP 9 GaP" 10 DIAMOND" 11 SILICA" 12 USER DEFINED DATA -1 INPUT MEDIUM 1 AND MEDIUM 2: -1 0 FILENAME OF DATA FOR MEDIUM 1: ALS.DAT MEDIUM 2 IS A VACUUM DO YOU WANT TO SEE THE DIELECTRIC DATA (Y/N)? N DO YOU WANT TO SEE THE PROBABILITY INFO? (Y/N) (Reply no unless you know what you are doing.) N PARAMETER CHECK. PRESS Y IF O.K., OR N IF YOU WANT TO ENTER A NEW VALUE. ENERGY OF INCIDENT ELECTRONS (KeV) = E = 100.000 Y ENERGY OF INCIDENT ELECTRONS (KeV) = E = 100.000 KYMIN AND KYMAX ARE THE MINIMUM AND MAXIMUM COLLECTION SEMIANGLES WITH A SCALE OF 1 unit = 0.5/E(eV) mrad. For 100keV electrons 200 = 1mrad. KYMIN = 0.000 Y KYMIN = 0.000 KYMAX = 1600.000 Y KYMAX = 1600.000 DKY0 is the first step length used in the integration over Ky. ALPHA is used to increase the step length with Ky during the integration. FIRST STEP DKY0 = 0.100 Y FIRST STEP DKY0 = 0.100 ALPHA = STEP/Ky = 0.500 Y ALPHA = STEP/Ky = 0.500 If you wish to make the calculation non-relativistically, set BETA=0. BETA = 0.548 Y BETA = 0.548 INITIAL DISTANCE OF THE ELECTRON FROM THE INTERFACE (Angstroms) = 100.000 Y INITIAL DISTANCE OF THE ELECTRON FROM THE INTERFACE (Angstroms) = 100.000 INPUT THETA = ANGLE OF INCIDENCE OF THE ELECTRONS TO THE INTERFACE (radians) (CAN BE 0): 0.003 THETA = 0.003 INNER POTENTIAL OF MEDIUM 1 = 18.00 eV. Y INNER POTENTIAL OF MEDIUM 1 = 18.00 DO YOU WANT A GRAPH (Y/N)? N ARE THE PARAMETERS O.K.? (Y/N): Y DO YOU WISH TO SPECIFY THE ENERGY WIDTH OF THE INCIDENT BEAM (Y/N)? Y INPUT FULL WIDTH AT HALF HEIGHT OF THE INCIDENT BEAM (eV): 1 INPUT OUTPUT FILENAME FOR CALCULATED SPECTRA: ALOUT.DAT THE ELECTRON IS NOW AT A DISTANCE OF 100.000 ANGSTROMS FROM THE INTERFACE. DO YOU WISH TO CONTINUE (Y/N)? Y FINAL DISTANCE OF THE ELECTRONS FROM THE INTERFACE (ANGSTROMS)= 20 KYMAX= 1600.000 O.K.? (Y/N) Y KYMAX= 1600.000 ANGLE OF ELECTRON PATH TO THE INTERFACE= 0.003 radians. INPUT NUMBER OF STEPS INTO WHICH THE ELECTRON PATH IS TO BE DIVIDED: 10 please wait.......calculation in progress. THE ELECTRON IS NOW AT A DISTANCE OF 20.000 ANGSTROMS FROM THE INTERFACE. DO YOU WISH TO CONTINUE (Y/N)? Y FINAL DISTANCE OF THE ELECTRONS FROM THE INTERFACE (ANGSTROMS)= 0 KYMAX= 1600.000 O.K.? (Y/N) Y KYMAX= 1600.000 ANGLE OF ELECTRON PATH TO THE INTERFACE= 0.003 radians. INPUT NUMBER OF STEPS INTO WHICH THE ELECTRON PATH IS TO BE DIVIDED: 10 please wait.......calculation in progress. THE ELECTRON IS NOW AT A DISTANCE OF 0.000 ANGSTROMS FROM THE INTERFACE. DO YOU WISH TO CONTINUE (Y/N)? Y FINAL DISTANCE OF THE ELECTRONS FROM THE INTERFACE (ANGSTROMS)= -5 KYMAX= 1600.000 O.K.? (Y/N) Y KYMAX= 1600.000 ANGLE OF ELECTRON PATH TO THE INTERFACE= 0.014 radians. INPUT NUMBER OF STEPS INTO WHICH THE ELECTRON PATH IS TO BE DIVIDED: 10 please wait.......calculation in progress. THE ELECTRON IS NOW AT A DISTANCE OF -5.000 ANGSTROMS FROM THE INTERFACE. DO YOU WISH TO CONTINUE (Y/N)? N THE OUTPUT FILE PRODUCED BY THE ABOVE PROGRAM RUN: ZERO LOSS PEAK INTENSITY= 10109.134 ENERGY PROBABILITY FINAL LOSS OF LAST STEP SPECTRUM 1.000 0.1292D-02 1500.0079 2.000 0.2862D-02 255.9052 3.000 0.2516D-02 231.6028 4.000 0.2825D-02 242.4906 5.000 0.3512D-02 285.5603 6.000 0.4632D-02 364.5952 7.000 0.6912D-02 534.0717 8.000 0.1243D-01 1021.0975 9.000 0.3464D-01 5561.1537 10.000 0.4600D+00 27391.1123 11.000 0.2045D+00 10718.6308 12.000 0.1514D-01 2349.7541 13.000 0.5692D-02 1175.3002 14.000 0.2989D-02 1044.3196 15.000 0.1986D-02 1123.4728 16.000 0.1655D-02 1365.2915 17.000 0.1286D-02 1885.3434 18.000 0.1041D-02 3292.3165 19.000 0.8610D-03 11145.1720 20.000 0.7197D-03 42450.2752 21.000 0.6207D-03 24062.5074 22.000 0.5357D-03 8148.4303 23.000 0.4626D-03 3287.8088 24.000 0.3994D-03 2409.4517 25.000 0.3448D-03 2375.5717 26.000 0.3043D-03 2735.6897 27.000 0.2682D-03 3590.6723 28.000 0.2360D-03 5819.8024 29.000 0.2073D-03 15463.6870 30.000 0.1816D-03 49103.3315 31.000 0.1641D-03 34155.4438 32.000 0.1481D-03 14986.4219 33.000 0.1336D-03 6421.7827 34.000 0.1203D-03 3998.8188 35.000 0.1083D-03 3589.8275 36.000 0.9927D-04 3914.2073 37.000 0.9095D-04 4904.1614 38.000 0.8335D-04 7465.8988 39.000 0.7629D-04 16975.8512 40.000 0.6979D-04 47019.1014 41.000 0.6425D-04 37180.1564 42.000 0.5911D-04 19307.5265 43.000 0.5439D-04 9139.1923 44.000 0.4999D-04 5345.3286 45.000 0.4592D-04 4360.3841 46.000 0.4263D-04 4493.9081 47.000 0.3956D-04 5390.5201 48.000 0.3668D-04 7776.5507 49.000 0.3403D-04 15809.5247 50.000 0.3155D-04 39351.1456 Parameters used: Media 1 and 2:..Aluminium.............VACUUM Reflection case: THETA = 0.003000 radians Inner potential = 18.000000 eV Initial distance from interface = 100.00 Angstroms Final distance from interface = 20.00 Angstroms Number of steps used = 10 Minimum scattering angle = 0.0000 mrad Maximum scattering angle = 8.0000 mrad Beta = 0.5482 Incident electron energy = 100.0000 keV Energy FWHH of zero loss beam = 1.0000 eV ZERO LOSS PEAK INTENSITY= 108.317 ENERGY PROBABILITY FINAL LOSS OF LAST STEP SPECTRUM 1.000 0.4482D-03 16.6065 2.000 0.1205D-02 3.9167 3.000 0.1257D-02 3.7336 4.000 0.1593D-02 4.1015 5.000 0.2159D-02 5.0301 6.000 0.3024D-02 6.6432 7.000 0.4673D-02 9.9798 8.000 0.8391D-02 19.2006 9.000 0.2135D-01 98.8373 10.000 0.2076D+00 496.0777 11.000 0.5285D+00 350.5868 12.000 0.2308D-01 78.8898 13.000 0.7863D-02 37.2522 14.000 0.4068D-02 34.5248 15.000 0.2732D-02 38.4986 16.000 0.2332D-02 48.1337 17.000 0.1869D-02 67.3689 18.000 0.1565D-02 115.4255 19.000 0.1342D-02 346.3256 20.000 0.1164D-02 1285.4867 21.000 0.1044D-02 1405.6383 22.000 0.9372D-03 701.6077 23.000 0.8412D-03 239.8859 24.000 0.7546D-03 159.3391 25.000 0.6767D-03 157.5016 26.000 0.6215D-03 183.2397 27.000 0.5699D-03 240.5385 28.000 0.5213D-03 374.0283 29.000 0.4759D-03 847.5983 30.000 0.4333D-03 2475.1503 31.000 0.4070D-03 3314.3079 32.000 0.3819D-03 2424.3597 33.000 0.3579D-03 1161.7719 34.000 0.3350D-03 568.8990 35.000 0.3130D-03 463.3257 36.000 0.2981D-03 492.5145 37.000 0.2837D-03 604.7317 38.000 0.2700D-03 864.9973 39.000 0.2565D-03 1632.7324 40.000 0.2436D-03 3932.5686 41.000 0.2327D-03 5797.5903 42.000 0.2221D-03 5296.9913 43.000 0.2120D-03 3334.7937 44.000 0.2020D-03 1755.1708 45.000 0.1923D-03 1132.6663 46.000 0.1851D-03 1059.2090 47.000 0.1780D-03 1206.3395 48.000 0.1710D-03 1598.7762 49.000 0.1643D-03 2640.4734 50.000 0.1577D-03 5438.2855 Parameters used: Media 1 and 2:..Aluminium.............VACUUM Reflection case: THETA = 0.003000 radians Inner potential = 18.000000 eV Initial distance from interface = 20.00 Angstroms Final distance from interface = 0.00 Angstroms Number of steps used = 10 Minimum scattering angle = 0.0000 mrad Maximum scattering angle = 8.0000 mrad Beta = 0.5482 Incident electron energy = 100.0000 keV Energy FWHH of zero loss beam = 1.0000 eV ZERO LOSS PEAK INTENSITY= 48.132 ENERGY PROBABILITY FINAL LOSS OF LAST STEP SPECTRUM 1.000 0.2518D-04 7.3960 2.000 0.6706D-04 1.7806 3.000 0.6857D-04 1.7058 4.000 0.8483D-04 1.8819 5.000 0.1114D-03 2.3151 6.000 0.1500D-03 3.0634 7.000 0.2206D-03 4.6048 8.000 0.3745D-03 8.8448 9.000 0.9145D-03 45.1784 10.000 0.9505D-02 227.0718 11.000 0.1337D-01 167.4173 12.000 0.9107D-03 37.7720 13.000 0.8354D-03 18.0945 14.000 0.2478D-02 19.5403 15.000 0.4551D-01 34.8126 16.000 0.4793D-02 27.3620 17.000 0.1470D-02 33.8005 18.000 0.7662D-03 56.3234 19.000 0.4963D-03 164.9234 20.000 0.3578D-03 608.0609 21.000 0.2905D-03 694.0019 22.000 0.2387D-03 357.5462 23.000 0.1981D-03 125.8178 24.000 0.1656D-03 99.7879 25.000 0.1394D-03 158.9804 26.000 0.1246D-03 156.4820 27.000 0.1112D-03 141.7431 28.000 0.9921D-04 197.1828 29.000 0.8841D-04 426.0785 30.000 0.7865D-04 1219.7952 31.000 0.7304D-04 1691.8847 32.000 0.6779D-04 1285.0553 33.000 0.6285D-04 643.5295 34.000 0.5821D-04 370.3551 35.000 0.5383D-04 464.2263 36.000 0.5102D-04 517.2145 37.000 0.4834D-04 467.4153 38.000 0.4580D-04 514.8103 39.000 0.4333D-04 881.2480 40.000 0.4095D-04 2040.8602 41.000 0.3907D-04 3078.7154 42.000 0.3723D-04 2914.6082 43.000 0.3547D-04 1923.7395 44.000 0.3375D-04 1150.0147 45.000 0.3207D-04 1076.7891 46.000 0.3086D-04 1209.7550 47.000 0.2969D-04 1165.1959 48.000 0.2851D-04 1148.7031 49.000 0.2739D-04 1572.4680 50.000 0.2629D-04 3012.1919 Parameters used: Media 1 and 2:..Aluminium.............VACUUM Reflection case: THETA = 0.014000 radians Inner potential = 18.000000 eV Initial distance from interface = 0.00 Angstroms Final distance from interface = -5.00 Angstroms Number of steps used = 10 Minimum scattering angle = 0.0000 mrad Maximum scattering angle = 8.0000 mrad Beta = 0.5482 Incident electron energy = 100.0000 keV Energy FWHH of zero loss beam = 1.0000 eV C Title :PINT C Keywords :EELS - planar interface calculations C Computer :VAX 780 C Operating System :VAX/VMS version 4.7 C Programming Language :FORTRAN 77 C Hardware Requirements :Tektronics 4010 series graphic terminal or emulator C desirable but not essential. C Author(s) :C.A. Walsh - based on a BASIC program written by C A. Howie. C Correspondence Address :Microstructural Physics Group, Cavendish Laboratory, C Madingley Road, Cambridge, CB3 0HE, England. C C Source Code: C C The program uses the result obtained by Garcia-Molina, Gras-Marti, Howie C and Ritchie (J. Phys. C18 (1985) pp5335) to calculate the energy loss of C electrons travelling parallel to a planar interface. They obtained an C expression for the differential probability dP/(dzdhf), for an electron C travelling in the z direction and losing energy hf. The probability of C losing energy between hf-dhf/2 and hf+dhf/2 after travelling a distance dz C parallel to the interface is given by dP/(dzdhf) x dzdhf. C The program also uses this result to calculate the energy loss of electrons C incident at an angle to an interface (see Howie, Milne and Walls,(1985) Inst. C Phys. Conf. Ser. No. 78: Chapter 4, EMAG'85). C The program requires dielectric data for the media on either side of C the interface. A limited amount of data is stored inside the program for some C materials. Other data can be supplied to the program in data files with the C following format: C C First line: Exactly 24 characters as a title. C Second line: No. of data items energy interval between succesive items (eV) C Third line: Energy (eV) real part of epsilon imaginary part of epsilon C Fourth line: etc. C C This program also calls upon routines in the library NGRAPH produced by C Nestor Zaluzec. C C C PROGRAM FOR THE CALCULATION OF ELECTRON ENERGY LOSS SPECTRA FOR C GLANCING ANGLE OR TRANSMISSION CONFIGURATION, ALLOWING FOR MULTIPLE C SCATTERING. C C FIRST WRITTEN IN BASIC BY A. HOWIE, EXTENDED AND REWRITTEN IN FORTRAN C BY C.A. WALSH, CAVENDISH LABORATORY, CAMBRIDGE, ENGLAND. 1988. C C DEFINITION OF SOME VARIABLES: C C CHARACTER VARIABLES: C C Y or N - answer to question concerning writing the diagnostic C probability information. C G Y or N - is Y if want graphical output. C ANS Y or N - contains the answer to numerous questions put to the user C TITLE array containing the titles of the two media used in the calculation C FN1 Contains the filenames containing user defined dielectric data for C FN2 media 1 and 2. C FN3 Contains the filename of a file opened for output of spectra in the C NELS program data format. C C INTEGER VARIABLES C C M An array containing integers identifying media 1 and 2. C M(i) = -1,0,1,......or 12, for i=1 or 2. The numbers are defined C subroutine INDAT. C NMAT The number of different sets of internally held dielectric data. C This also dictates the maximum allowable value of M(i). C NPOINT The number of dielectric data points supplied by one of the C subroutines when called by subroutine INDAT. The internally C held dielectric data is then interpolated and extrapolated to C provide data at the required energy interval. C NRI Number of data points in the final intensity spectrum. C NZ The number of dielectric data points. C NS Number of steps used in the calculations. C IOPEN 0/1 - is 1 if a file is open on unit 15 for output of spectra. C C FLOATING POINT VARIABLES C C E Energy of the incident electrons in keV. C C In the following, positive values lie in medium 2, negative in medium 1: C X0 Initial distance of the electron from the interface (Angstroms). C In the glancing angle case the calculations are made in stages, each C stage using different path lengths and different numbers of steps (NS). C For each stage: C XL1 Initial distance of the electron from interface (Angstroms). C Xl2 Final distance of the electron from interface (Angstroms). C In the transmission case: C X0 Impact parameter (Angstroms). C XL1 = 0. C XL2 Thickness of the medium through which the electrons travel (Angs.) C C T Step length parallel to the interface (Angs.). ( T=XL2/NS in the C transmission case.) C X Distance from the interface at which the probability of energy C loss is to be calculated in subroutine PCINT. C BETA = electron velocity / speed of light if calculations are to be C made relativistically. BETA is set to 0 if the calculations are C to be made non-relativistically. C BETAP = electron velocity / speed of light. C XSCALE = electron speed * Planck's costant/2/pi/electronic charge * C 1E10 (SI units). It is used to scale distances in the C integration to correspond to the units of RKYMIN and RKYMAX. C ZI Energy interval at which calculations are made (eV). C RL Gives the required constant factor in front of the integral in C the calculation of the probability. C C THE ARRAYS C C ER, EI Arrays used in INDAT to hold the real and imaginary parts of the C internally stored dielectric data before interpolating and C extrapolating it to the desired energy interval and range. C IDE Array used when reading the internally stored dielectric data. C It contains the energy of the corresponding set of dielectric C data in ER and EI. C E1R, E1I Hold the real and imaginary parts of the dielectric data for C medium 1. C E2R, E2I Hold the real and imaginary parts of the dielectric data for C medium 2. C EP Contains the calculated probability of energy loss in an energy C interval ZI, after travelling a distance T/NS, as a function of C energy. C RI Contains the spectrum assuming 10E6 incident electrons. C IMPLICIT REAL*8 (A-H,O-Z) CHARACTER C, G, ANS CHARACTER*4 TITLE(18) CHARACTER*20 FN1,FN2,FN3 DIMENSION M(2) COMMON FN1,FN2,FN3 COMMON /RAYS/ ER(0:20),EI(0:20),E1R(501),E1I(501),E2R(501), + E2I(501),EP(501),RI(-100:801),IDE(0:20) COMMON /SINGS/ RMAX,E,X0,ALPHA,RKYMAX,DKY0,RHE,THETA,XL2,XL1, + THETAP,T,X,BETA,BETAP,RKYMIN,XSCALE,RL,ZI COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI COMMON /CRACT/ C,G,ANS COMMON /GRPH/ TITLE COMMON /WRTFIL/ IOPEN DO 98 L1=13,18 TITLE(L1)=' ' 98 CONTINUE FN2 = ' ' FN3 = ' ' NPOINT = 0 NMAT = 12 NREP = 100 IOPEN = 0 RMAX = 0.0 E = 100.0 X0 = 100.0 ALPHA= 0.5 RKYMAX = 1600.0 DKY0 = 0.1 DO 90 L1=1,12 TITLE(L1)=' ' 90 CONTINUE C C ZI=ENERGY INTERVAL BETWEEN DATA POINTS. C NZ IS THE NO. OF DIELECTRIC DATA POINTS. C NRI IS THE NUMBER OF DATA POINTS IN THE FINAL INTENSITY SPECTRUM C ZI = 1.0 NZ = 40 NRI = 80 C C GET DIELECTRIC DATA C CALL INDAT( M(1) ) C C ASK IF WANT TO SEE DATA C WRITE( 6,2 ) READ( 5,12 ) C IF ( C.EQ.'N' ) GO TO 250 DO 100 L1=1,NZ ERGY = L1*ZI WRITE( 6,3 ) ERGY, E1R(L1), E1I(L1), E2R(L1), E2I(L1) 100 CONTINUE DO 200 L1=1,NPOINT WRITE( 6,4 ) IDE(L1), ER(L1), EI(L1) 200 CONTINUE C C ASK IF WANT TO SEE PROBABILITY DATA C 250 WRITE( 6,14 ) READ( 5,12 ) C IF ( C.NE.'Y'.AND.C.NE.'N' ) GO TO 250 300 CALL PRAMS( M(1),VINPOT ) DO 400 L1=-100,801 RI(L1) = 0.0 400 CONTINUE RI(0) = 1E6 420 WRITE ( 6,15 ) READ ( 5,12 ) ANS IF ( ANS.EQ.'Y' ) THEN CALL DISTR( RI(-100),SIGMA ) ELSE IF ( ANS.NE.'N' ) THEN GO TO 420 ENDIF C C OPEN OUTPUT FILE FOR SPECTRA C CALL OFILE C C INITIALISE GRAPHICS PACKAGE C XL2 = X0 IF ( ABS(THETA).LT.1D-10 ) XL2 = 0 DO 1100 L1=1,NREP XL1 = XL2 IF ( ABS(THETA).GE.1D-10 ) WRITE( 6,5 ) XL1 C C ASK IF WANT TO CONTINUE C 500 WRITE( 6,6 ) READ( 5,12 ) ANS IF ( ANS.EQ.'N') GO TO 1200 IF ( ANS.NE.'Y') GO TO 500 510 IF ( ABS(THETA).LT.1E-10 ) THEN WRITE( 6,7 ) X0 ELSE WRITE( 6,13 ) ENDIF READ( 5,* ) XL2 C C CHECK VALUE OF KYMAX C 700 WRITE( 6,8 ) RKYMAX CALL RNEW( RKYMAX,ANS ) WRITE( 6,1 ) RKYMAX 800 THETAP = THETA IF ( XL2.LT.0.0 ) THETAP = SQRT(THETA*THETA+VINPOT) 850 WRITE( 6,9 ) THETAP C C INPUT NUMBER OF STEPS C READ( 5,* ) NS IF ( NS.LT.1 ) GO TO 850 IF ( ABS(THETA).LT.1E-10 ) THEN T = ( XL2-XL1 )/NS ELSE T = 2.0*ABS( XL2-XL1 )/NS/TAN( THETAP ) ENDIF WRITE(6,16) DO 1000 L2=1,NS X = XL1 + ( XL1-XL2 )*( 0.5-L2 )/NS IF ( ABS(THETA).LT.1D-10 ) X=X0 IF ( ABS(THETA).LT.1D-10.AND.L2.GT.1 ) GO TO 900 CALL PCINT 900 CALL FOLD( RI(-100),EP(1) ) 1000 CONTINUE WRITE( 15,10 ) RI(0) NMAX = NZ IF ( NRI.GT.NMAX ) NMAX = NRI DO 1020 L2=1,NMAX ENERGY = L2*ZI ET = 0.0D0 IF ( L2.LE.NZ ) ET = EP(L2) WRITE( 15,11 ) ENERGY, ET, RI(L2) 1020 CONTINUE SA = RKYMAX/E/2.0 WRITE(15,17) (TITLE(J),J=1,18) IF ( ABS(THETA).LT.1D-10 ) THEN WRITE(15,18) X0,XL2 ELSE WRITE(15,19) VINPOT,THETAP,XL1,XL2 ENDIF SAM = RKYMIN/E/2.0 TSIGMA = SIGMA*2.0 WRITE(15,20) NS,SAM,SA,BETA,E,TSIGMA IF ( G.EQ.'Y' ) CALL GRAPH( RI,SIGMA ) IF ( ABS(THETA).LT.1D-10 ) GO TO 1200 1100 CONTINUE 1200 CONTINUE STOP 1 FORMAT (' KYMAX=', F10.3 ) 2 FORMAT (' DO YOU WANT TO SEE THE DIELECTRIC DATA (Y/N)?' ) 3 FORMAT ( 1X, F8.3, 4( 2X, F10.5 ) ) 4 FORMAT ( 1X, I3, 2( 2X, F10.5 ) ) 5 FORMAT ( ' THE ELECTRON IS NOW AT A DISTANCE OF ',F10.3, + ' ANGSTROMS FROM THE INTERFACE.' ) 6 FORMAT ( ' DO YOU WISH TO CONTINUE (Y/N)?' ) 7 FORMAT ( ' THETA=0 - ELECTRONS TRAVELLING PARALLEL TO ', + 'INTERFACE.'/' IMPACT PARAMETER =',F10.3, + /' INPUT THE TOTAL DISTANCE TRAVELLED BY THE ', + 'ELECTRONS IN THE MEDIUM:' ) 8 FORMAT ( ' KYMAX=', F10.3/' O.K.? (Y/N)' ) 9 FORMAT ( ' ANGLE OF ELECTRON PATH TO THE INTERFACE=', F10.3, + ' radians.'/' INPUT NUMBER OF STEPS INTO WHICH THE', + ' ELECTRON PATH IS TO BE DIVIDED:' ) 10 FORMAT ( ' ZERO LOSS PEAK INTENSITY=', F10.3/ 5X, + ' ENERGY PROBABILITY FINAL '/ + 5X,' LOSS OF LAST STEP SPECTRUM' ) 11 FORMAT ( 3X, F8.3, D12.4, F12.4 ) 12 FORMAT ( 1A ) 13 FORMAT( ' FINAL DISTANCE OF THE ELECTRONS FROM THE INTERFACE', + ' (ANGSTROMS)=' ) 14 FORMAT ( ' DO YOU WANT TO SEE THE PROBABILITY INFO? (Y/N)'/ + ' (Reply no unless you know what you are doing.)' ) 15 FORMAT ( ' DO YOU WISH TO SPECIFY THE ENERGY WIDTH OF THE '/ + ' INCIDENT BEAM (Y/N)?' ) 16 FORMAT ( ' please wait.......calculation in progress.' ) 17 FORMAT ( ' Parameters used:'/' Media 1 and 2:',18A4 ) 18 FORMAT ( ' Transmission case (THETA=0):'/ + ' Impact parameter = ',F10.2,' Angstroms'/ + ' Distance travelled past interface =', F10.2, + ' Angstroms'/ ) 19 FORMAT ( ' Reflection case: THETA = ',F10.6,' radians'/ + ' Inner potential = ',F10.6,' eV'/ + ' Initial distance from interface = ',F10.2, + ' Angstroms'/ + ' Final distance from interface = ',F10.2, + ' Angstroms' ) 20 FORMAT ( ' Number of steps used = ',I6/ + ' Minimum scattering angle = ',F10.4,' mrad'/ + ' Maximum scattering angle = ',F10.4,' mrad'/ + ' Beta = ',F10.4/ + ' Incident electron energy = ',F10.4,' keV'/ + ' Energy FWHH of zero loss beam = ',F10.4,' eV' ) END C From C.A. Walsh, Cavendish Laboratory, Cambridge, England. C C The routine INDAT in the program for calculating the electron C energy loss of electrons passing a planar interface contained C a bug. When extrapolating the internally held dielectric data C to LOW energies, this bug caused the imaginary part C to be set equal to zero, for energies approx equal to zero. In C fact, the imaginary part should be set to some large number C (100 in this program). This bug will only cause errors when a C metallic medium is used from the list of dielectric data offered C by the program (ie aluminium). However, the error this causes in C the calculated spectra is SERIOUS, producing extraneous peaks in C the low energy region. C C Below is a corrected version of the offending routine. C C I hope this doesn't cause you too much inconvenience. C C I would recommend using your own dielectric data, rather than the C data inside the program, when trying to do any serious work. C C CCCCCCCCCCCCCCCC C C C C C INDAT C C C C C CCCCCCCCCCCCCCCC SUBROUTINE INDAT( M ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C This routine puts the dielectric data into the appropiate C C arrays. In the case of internally stored data it also C C interpolates and extrapolates it to the required energy C C interval and range. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION M(2),IU(2) CHARACTER C, G, ANS CHARACTER*4 TITLE(18) COMMON /GRPH/ TITLE COMMON /RAYS/ ER(0:20),EI(0:20),E1R(501),E1I(501),E2R(501), + E2I(501),EP(501),RI(-100:801),IDE(0:20) COMMON /SINGS/ RMAX,E,X0,ALPHA,RKYMAX,DKY0,RHE,THETA,XL2,XL1, + THETAP,T,X,BETA,BETAP,RKYMIN,XSCALE,RL,ZI COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI COMMON /CRACT/ C,G,ANS EI(0) = 0.0D0 100 WRITE( 6,1 ) READ( 5,* ) M(1), M(2) C C CHECK LIMITS ON M C IF ( M(1).LT.1.AND.M(1).NE.-1 ) GO TO 100 IF ( M(1).GT.NMAT ) GO TO 100 IF ( M(2).LT.-1.OR.M(2).GT.NMAT ) GO TO 100 C C OPEN USER DEFINED DATA FILES AND FIND NZ AND ZI C C PUT UNIT NUMBER OF FILE FOR MEDIUM L1 IN IU(L1) C IU(1) = 0 IU(2) = 0 DO 150 L1=1,2 IF ( M(L1).EQ.-1 ) CALL OPUSFL( IU(1),L1 ) 150 CONTINUE C C PUT DIELECTRIC DATA INTO ARRAYS ER AND EI C DO 600 L1=1,2 IF ( M(L1).EQ.-1 ) THEN IF ( L1.EQ.2.AND.RMAX.GT.4.0 ) GO TO 600 CALL USDE( IU(L1),L1,E2R(1),E2I(1) ) GO TO 450 ELSE IF ( M(L1).EQ.0 ) THEN CALL VACUUM( L1,E2R(1),E2I(1),TITLE(1) ) GO TO 600 ELSE IF ( M(L1).EQ.1 ) THEN CALL MGO( L1,IDE(0),ER(0),EI(0),TITLE ) ELSE IF ( M(L1).EQ.2 ) THEN CALL SILICN( L1,IDE(0),ER(0),EI(0),TITLE ) ELSE IF ( M(L1).EQ.3 ) THEN CALL SILICA( L1,IDE(0),ER(0),EI(0),TITLE ) ELSE IF ( M(L1).EQ.4 ) THEN CALL CARBON( L1,IDE(0),ER(0),EI(0),TITLE ) ELSE IF ( M(L1).EQ.5 ) THEN CALL DIAMD( L1,IDE(0),ER(0),EI(0),TITLE ) ELSE IF ( M(L1).EQ.6 ) THEN CALL ALMNM( L1,IDE(0),ER(0),EI(0),TITLE ) ELSE IF ( M(L1).EQ.7 ) THEN CALL ALMINA( L1,IDE(0),ER(0),EI(0),TITLE ) ELSE IF ( M(L1).EQ.8 ) THEN CALL GAAS( L1,IDE(0),ER(0),EI(0),TITLE ) ELSE IF ( M(L1).EQ.9 ) THEN CALL GAP( L1,IDE(0),ER(0),EI(0),TITLE ) ELSE IF ( M(L1).EQ.10 ) THEN CALL GAPP( L1,IDE(0),ER(0),EI(0),TITLE ) ELSE IF ( M(L1).EQ.11 ) THEN CALL DIAMDP( L1,IDE(0),ER(0),EI(0),TITLE ) ELSE IF ( M(L1).EQ.12 ) THEN CALL SILCAP( L1,IDE(0),ER(0),EI(0),TITLE ) ENDIF C C INTERPOLATE DIELECTRIC DATA C IDE(0)= 0 ER(0) = ER(1) IF ( ER(0).LT.0.0D0 ) ER(0) = -100.0D0 EI(0) = 0.00001D0 IF ( ER(0).LT.0.0D0 ) EI(0) = 100.0D0 J = 1 ZZ = ZI DO 300 L2=1,NPOINT RJI = IDE(L2) - IDE(L2-1) 200 IF ( (ZZ-IDE(L2) ).GT.1D-4 ) GO TO 300 E2R(J) = ER(L2-1)+(ER(L2)-ER(L2-1))*(ZZ-IDE(L2-1))/RJI E2I(J) = EI(L2-1)+(EI(L2)-EI(L2-1))*(ZZ-IDE(L2-1))/RJI IF ( E2I(J).LT.0.001D0 ) E2I(J) = 0.001D0 ZZ = ZZ + ZI J = J + 1 GO TO 200 300 CONTINUE C C EXTRAPOLATE THE DIELECTRIC DATA C A = (E2R(J-1)-1)*(J-1)*(J-1) B = E2I(J-1)*(J-1)*(J-1) DO 400 L2=J,NZ E2R(L2) = A/L2/L2 + 1.0 E2I(L2) = B/L2/L2 400 CONTINUE 450 IF ( L1.EQ.1 ) THEN DO 500 L2=1,NZ E1R(L2) = E2R(L2) E1I(L2) = E2I(L2) 500 CONTINUE ENDIF 600 CONTINUE RETURN 1 FORMAT (' GLANCING ANGLE PROGRAM FOR'/ + ' VACUUM 0 (MEDIUM 2 ONLY)'/' MgO 1'/ + ' SILICON 2'/' SILICA 3'/' CARBON 4'/ + ' DIAMOND 5'/' ALUMINIUM 6'/' ALUMINA 7'/ + ' GaAs 8'/' GaP 9'/' GaP" 10'/ + ' DIAMOND" 11'/' SILICA" 12'/ + ' USER DEFINED DATA -1'/ + ' INPUT MEDIUM 1 AND MEDIUM 2:' ) END C CCCCCCCCCCCCCCCC C C C C C PRAMS C C C C C CCCCCCCCCCCCCCCC SUBROUTINE PRAMS( M,VINPOT ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C This routine asks the user for various parameters. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) CHARACTER C, G, ANS DIMENSION M(2) COMMON /RAYS/ ER(0:20),EI(0:20),E1R(501),E1I(501),E2R(501), + E2I(501),EP(501),RI(-100:801),IDE(0:20) COMMON /SINGS/ RMAX,E,X0,ALPHA,RKYMAX,DKY0,RHE,THETA,XL2,XL1, + THETAP,T,X,BETA,BETAP,RKYMIN,XSCALE,RL,ZI COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI COMMON /CRACT/ C,G,ANS 100 GS = ( 1.0+(E/511.0) )*( 1.0+(E/511.0) ) BETAP = SQRT( 1.0-(1.0/GS) ) BETA = BETAP RKYMIN = 0.0 XSCALE = 1973.0*BETA RL = 511.0*0.5292*3.14159265*BETA*BETA/20.0 200 WRITE( 6,1 ) C C CHECK VALUE OF ELECTRON ENERGY, E C WRITE( 6,2 ) E CALL RNEW( E,ANS ) IF ( ANS.EQ.'N' ) GO TO 100 WRITE( 6,2 ) E C C CHECK ON KYMIN, RKYMIN C WRITE ( 6,3 ) RKYMIN CALL RNEW( RKYMIN,ANS ) WRITE ( 6,22 ) RKYMIN C C CHECK ON KYMAX, RKYMAX C WRITE ( 6,4 ) RKYMAX CALL RNEW( RKYMAX,ANS ) WRITE ( 6,4 ) RKYMAX C C CHECK ON DKY0 C WRITE( 6,5 ) DKY0 CALL RNEW( DKY0,ANS ) WRITE( 6,18 ) DKY0 C C CHECK ON ALPHA C WRITE( 6,6 ) ALPHA CALL RNEW( ALPHA,ANS ) WRITE( 6,6 ) ALPHA C C CHECK ON BETA C WRITE( 6,7 ) BETA CALL RNEW( BETA,ANS ) WRITE( 6,19 ) BETA C C CHECK ON X0 - THE INITIAL DISTANCE FROM THE INTERFACE C WRITE( 6,8 ) X0 250 READ( 5,15 ) ANS IF ( ANS.EQ.'Y' ) THEN GO TO 270 ELSE IF ( ANS.NE.'N' ) THEN GO TO 250 ENDIF WRITE( 6,16 ) READ( 5,* ) X0 270 WRITE( 6,8 ) X0 C C ASK FOR THETA - ANGLE OF INCIDENCE TO THE INTERFACE C 300 WRITE( 6,9 ) READ( 5,* ) THETA WRITE( 6,10 ) THETA C C IF THE PROGRAM IS RUNNING IN THE GLANCING ANGLE CONFIGURATION AND C ONE OF THE MEDIA IS VACUUM THEN NEED TO TAKE INTO ACCOUNT THE C REFRACTION OF THE ELECTRON BEAM AT THE INTERFACE. THIS DEPENDS C UPON THE INNER POTENTIAL OF THE MEDIUM. IN THE CASE OF NEITHER C MEDIA BEING A VACUUM, THE REFRACTION OF THE BEAM IS IGNORED. C VINPOT = 0.0D0 IF ( ABS(THETA).GT.1D-10.AND.M(2).EQ.0 ) THEN VINPOT = 18.0 WRITE ( 6,20 ) VINPOT CALL RNEW( VINPOT,ANS ) WRITE ( 6,21 ) VINPOT VINPOT = VINPOT/E/1000.0D0 ENDIF C C ASK IF WANT A GRAPH C 400 WRITE( 6,11 ) READ( 5,15 ) G IF ( G.EQ.'N' ) THEN GO TO 500 ELSE IF ( G.NE.'Y' ) GO TO 400 ENDIF C C A GRAPH IS REQUIRED. ASK FOR HIGHEST ENERGY ON GRAPH. C XXHI = 801.0*ZI WRITE( 6,12 ) XXHI READ( 5,* ) RHE WRITE( 6,13 ) RHE NRI = RHE/ZI IF ( NRI.LE.800 ) GO TO 500 WRITE( 6,17 ) GO TO 400 C C DOUBLE CHECK PARAMETERS C 500 WRITE( 6,14 ) READ( 5,15 ) ANS IF ( ANS.EQ.'N' ) THEN GO TO 200 ELSE IF ( ANS.NE.'Y' ) THEN GO TO 500 ENDIF RETURN 1 FORMAT ( ' PARAMETER CHECK. PRESS Y IF O.K., OR N IF YOU'/ + ' WANT TO ENTER A NEW VALUE.' ) 2 FORMAT ( ' ENERGY OF INCIDENT ELECTRONS (KeV) = E = ', F10.3 ) 3 FORMAT ( ' KYMIN AND KYMAX ARE THE MINIMUM AND MAXIMUM ', + ' COLLECTION SEMIANGLES '/' WITH A SCALE OF ', + ' 1 unit = 0.5/E(eV) mrad. For 100keV electrons ', + ' 200 = 1mrad.'/' KYMIN = ', F10.3 ) 4 FORMAT ( ' KYMAX = ', F10.3 ) 5 FORMAT ( ' DKY0 is the first step length used in the', + ' integration over Ky.'/' ALPHA is used to increase', + ' the step length with Ky during the integration.'/ + ' FIRST STEP DKY0 = ', F10.3 ) 6 FORMAT ( ' ALPHA = STEP/Ky = ', F10.3 ) 7 FORMAT ( ' If you wish to make the calculation', + ' non-relativistically, set BETA=0.'/' BETA = ', F10.3 ) 8 FORMAT ( ' INITIAL DISTANCE OF THE ELECTRON FROM THE ', + ' INTERFACE (Angstroms) = ', F10.3 ) 9 FORMAT ( ' INPUT THETA = ANGLE OF INCIDENCE OF THE ELECTRONS', + ' TO THE INTERFACE (radians)'/' (CAN BE 0):' ) 10 FORMAT ( ' THETA = ', F10.3 ) 11 FORMAT ( ' DO YOU WANT A GRAPH (Y/N)?' ) 12 FORMAT ( ' UP TO WHAT ENERGY (eV,',F6.0,' is maximum)?' ) 13 FORMAT ( ' HIGHEST ENERGY = ', F10.3 ) 14 FORMAT ( ' ARE THE PARAMETERS O.K.? (Y/N):' ) 15 FORMAT ( 1A ) 16 FORMAT ( ' IMPACT PARAMETER =' ) 17 FORMAT ( ' TOO HIGH!' ) 18 FORMAT ( ' FIRST STEP DKY0 = ', F10.3 ) 19 FORMAT ( ' BETA = ', F10.3 ) 20 FORMAT ( ' INNER POTENTIAL OF MEDIUM 1 = ',F6.2,' eV.' ) 21 FORMAT ( ' INNER POTENTIAL OF MEDIUM 1 = ',F6.2 ) 22 FORMAT ( ' KYMIN = ', F10.3 ) END C CCCCCCCCCCCCCCCC C C C C C PCINT C C C C C CCCCCCCCCCCCCCCC SUBROUTINE PCINT CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C This routine calculates the probability of energy loss as a C C function of energy when the electron is at a distance X from C C the interface. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) CHARACTER C, G, ANS COMMON /RAYS/ ER(0:20),EI(0:20),E1R(501),E1I(501),E2R(501), + E2I(501),EP(501),RI(-100:801),IDE(0:20) COMMON /SINGS/ RMAX,E,X0,ALPHA,RKYMAX,DKY0,RHE,THETA,XL2,XL1, + THETAP,T,X,BETA,BETAP,RKYMIN,XSCALE,RL,ZI COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI COMMON /CRACT/ C,G,ANS DO 500 L1=1,NZ IF ( X.LT.0.0D0 ) GO TO 100 RE1R = E1R(L1) RE2R = E2R(L1) RE1I = E1I(L1) RE2I = E2I(L1) GO TO 200 100 RE1R = E2R(L1) RE2R = E1R(L1) RE1I = E2I(L1) RE2I = E1I(L1) 200 A2R = 1.0-RE2R*BETA*BETA A2I = -RE2I*BETA*BETA A1R = 1.0-RE1R*BETA*BETA A1I = -RE1I*BETA*BETA XS = 2.0*ZI*L1*ABS(X)/XSCALE RKY = RKYMIN/L1/ZI RKYM = RKYMAX/L1/ZI DKY = DKY0 S1R = 0.0D0 S1I = 0.0D0 S2R = 0.0D0 S2I = 0.0D0 S3R = 0.0D0 S3I = 0.0D0 C C DO THE INTEGRAL C DO 300 L2=1,10000 RKY = RKY + 0.5*DKY THR = RKY*RKY + A1R TH1 = SQRT(THR*THR+A1I*A1I) TH1R = SQRT(0.5*(THR+TH1)) TH1I = 0.5*A1I/TH1R THR = RKY*RKY + A2R TH2 = SQRT(THR*THR+A2I*A2I) TH2R = SQRT(0.5*(THR+TH2)) TH2I = 0.5*A2I/TH2R EX = EXP(-TH2R*XS) CS = COS(TH2I*XS)*EX SN = -SIN(TH2I*XS)*EX ETHR = RE2R*TH2R-RE2I*TH2I ETHI = RE2R*TH2I+RE2I*TH2R D2R = (TH1R+TH2R)*ETHR-(TH1I+TH2I)*ETHI D2I = (TH1R+TH2R)*ETHI+(TH1I+TH2I)*ETHR BR = RE2R*TH1R-RE2I*TH1I+RE1R*TH2R-RE1I*TH2I BI = RE2R*TH1I+RE2I*TH1R+RE1R*TH2I+RE1I*TH2R D1R = D2R*BR-D2I*BI D1I = D2R*BI+D2I*BR D1 = D1R*D1R+D1I*D1I D2 = D2R*D2R+D2I*D2I ETH = ETHR*ETHR+ETHI*ETHI F1R = 2.0*(THR*(RE1R-RE2R)-A2I*(RE1I-RE2I)) F1I = 2.0*(THR*(RE1I-RE2I)+A2I*(RE1R-RE2R)) RH1R = (F1R*D1R+F1I*D1I)/D1 RH1I = (-F1R*D1I+F1I*D1R)/D1 S1R = S1R+(RH1R*CS-RH1I*SN)*DKY/RL S1I = S1I+(RH1R*SN+RH1I*CS)*DKY/RL F2R = ((TH2R-TH1R)*D2R+(TH2I-TH1I)*D2I)/D2 F2I = ((TH1R-TH2R)*D2I+(TH2I-TH1I)*D2R)/D2 RH2R = A2R*F2R-A2I*F2I RH2I = A2R*F2I+A2I*F2R S2R = S2R-(RH2R*CS-RH2I*SN)*DKY/RL S2I = S2I-(RH2R*SN+RH2I*CS)*DKY/RL S3R = S3R-(A2R*ETHR+A2I*ETHI)*DKY/ETH/RL S3I = S3I-(-A2R*ETHI+A2I*ETHR)*DKY/ETH/RL RKY = RKY + 0.5*DKY IF ( RKY.GT.RKYM ) GO TO 350 DKY = RKY*ALPHA IF ( (RKY+DKY).GT.RKYM ) DKY = RKYM-RKY+0.001 300 CONTINUE 350 IF ( C.EQ.'N' ) GO TO 400 STOT = S1I+S2I+S3I ENERGY = L1*ZI WRITE( 6,1 ) ENERGY,S1I,S2I,S3I,STOT 400 EP(L1) = (S1I+S2I+S3I)*T/10000.0*ZI 500 CONTINUE RETURN 1 FORMAT ( 1X, F11.3, 4( 3X, D14.6 ) ) END SUBROUTINE OPUSFL( IU,L ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C This routine is called when the user wishes to supply C C dielectric data from an external file. The routine reads C C in the filename and opens it on unit IUL and then reads the C C top two lines of information in the file which contains a C C title, the number of sets of data and the energy interval C C (eV) of the data. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) CHARACTER*20 FN1,FN2,FN3 CHARACTER*4 TITLE(18) CHARACTER ANS DIMENSION IU(2) COMMON FN1,FN2,FN3 COMMON /SINGS/ RMAX,E,X0,ALPHA,RKYMAX,DKY0,RHE,THETA,XL2,XL1, + THETAP,T,X,BETA,BETAP,RKYMIN,XSCALE,RL,ZI COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI COMMON /GRPH/ TITLE FN1 = ' ' NZZ = NZ 100 WRITE( 6,1 ) L READ( 5,2 ) FN1 IF ( FN1.EQ.FN2.AND.FN1.NE.' ' ) THEN RMAX = 3.0 + L GO TO 150 ENDIF IF ( FN1.EQ.FN3.AND.FN1.NE.' ' ) THEN IF ( RMAX.GT.1.0 ) RMAX = 9.0 IF ( RMAX.LT.1.0 ) RMAX =-3.0 GO TO 150 ELSE IF ( FN1.EQ.' ' ) THEN GO TO 100 ENDIF IU(L) = 11 + L IUL = IU(L) OPEN( UNIT=IUL, IOSTAT=IOST, ERR=500, FILE=FN1, STATUS='OLD', + ACCESS='SEQUENTIAL',FORM='FORMATTED' ) READ( IUL,3 ) (TITLE(J),J=6*(L-1)+1,6*L) READ( IUL,* ) NZ,ZI IF ( L.EQ.2.AND.IU(1).NE.0.AND.NZ.NE.NZZ) GO TO 300 IF ( NZ.LE.0 ) GO TO 300 IF ( NZ.GT.501 ) GO TO 400 IF ( L.EQ.1 ) FN2 = FN1 150 IF ( L.EQ.2 ) FN3 = FN1 RETURN C C ERROR ON READING FILE C 200 WRITE( 6,4 ) FN1 READ( 5,5 ) ANS IF ( ANS.EQ.'N' ) STOP IF ( ANS.EQ.'Y' ) GO TO 100 GO TO 200 C C ERROR IN USER DEFINED DATA C 300 WRITE( 6,6 ) NZZ,NZ STOP 400 WRITE( 6,7 ) NZ STOP C C ERROR ON OPENING FILE C 500 WRITE( 6,8 ) L, IOST WRITE( 6,456 ) FN1 456 FORMAT( 1X,A20 ) STOP 1 FORMAT ( ' FILENAME OF DATA FOR MEDIUM ',I2,':' ) 2 FORMAT ( A20 ) 3 FORMAT ( 6A4 ) 4 FORMAT ( ' ERROR WHILST READING FILE ',A20/' DO YOU WISH TO', + ' CONTINUE (Y/N)?' ) 5 FORMAT ( A1 ) 6 FORMAT ( ' USER DEFINED DATA IS INCOMPATIBLE'/' NZZ=',I5/ + ' NZ =',I5 ) 7 FORMAT ( ' NZ=',I5,' ARRAYS ARE TOO SMALL ' ) 8 FORMAT ( ' ERROR ON OPENING FILE FOR MEDIUM',I3,' IOSTAT=',I6 ) END SUBROUTINE USDE( IU,L,E2R,E2I ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C This routine reads in the dielectric data from the user- C C supplied file. The real and imaginary parts of the data are C C read into the arrays E2R and E2I. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION E2R(501),E2I(501) COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI DO 100 L1=1,NZ READ( IU,*,ERR=300,END=400 ) BUF1,E2R(L1),E2I(L1) 100 CONTINUE DO 150 L1=1,NZ IF ( E2I(L1).LT.0.001D0 ) E2I(L1) = 0.001D0 150 CONTINUE RETURN 300 WRITE( 6,1 ) IU STOP 400 WRITE( 6,2 ) IU,L STOP 1 FORMAT ( ' ERROR OCCURRED WHILST READING FILE ON UNIT',I4 ) 2 FORMAT ( ' END OF FILE REACHED ON UNIT',I4,' WHILST READING ', + 'DATA FOR MEDIUM ',I4 ) C 3 FORMAT ( 1X,6D13.6 ) END SUBROUTINE RNEW( RNV, ANS ) IMPLICIT REAL*8 (A-H,O-Z) CHARACTER ANS 100 READ( 5,2 ) ANS IF ( ANS.EQ.'Y' ) THEN GO TO 300 ELSE IF ( ANS.NE.'N' ) THEN GO TO 100 ENDIF 200 WRITE( 6,1 ) READ( 5,* ) RNV IF ( RNV.LT.0.0 ) GO TO 200 300 CONTINUE RETURN 1 FORMAT ( ' ENTER NEW VALUE:' ) 2 FORMAT ( 1A ) END SUBROUTINE DISTR ( RI,SIGMA ) IMPLICIT REAL*8 (A-H,O-Z) CHARACTER*1 ANS COMMON /SINGS/ RMAX,E,X0,ALPHA,RKYMAX,DKY0,RHE,THETA,XL2,XL1, + THETAP,T,X,BETA,BETAP,RKYMIN,XSCALE,RL,ZI DIMENSION RI(-100:801) 100 WRITE (6,1) READ (5,*) TSIGMA IF ( TSIGMA.LT.0.0D0 ) GO TO 100 IF ( TSIGMA.GT.10.0D0 ) THEN 120 WRITE (6,2) TSIGMA READ (5,3) ANS IF ( ANS.EQ.'Y' ) GO TO 130 IF ( ANS.NE.'N' ) GO TO 120 GO TO 100 ENDIF 130 SUM = 0.0D0 SIGMA = TSIGMA/2.0D0 PI = 3.14159 FACT = 1.0D0/SIGMA/SQRT( 2.0D0*PI )*ZI EFAC = -0.5D0/SIGMA/SIGMA*ZI*ZI DO 140 L1=0,100 RI(L1) = EXP( EFAC*L1*L1 )*FACT*1D6 RI(-L1) = RI(L1) SUM = SUM + RI(L1) IF ( L1.GT.0 ) SUM = SUM + RI(L1) 140 CONTINUE FACT = 1D6/SUM DO 150 L1=-100,100 RI(L1) = RI(L1)*FACT 150 CONTINUE RETURN 1 FORMAT ( ' INPUT FULL WIDTH AT HALF HEIGHT OF THE INCIDENT BEAM', + ' (eV):' ) 2 FORMAT ( ' FWHH = ',E12.4,' ARE YOU SURE (Y/N)?' ) 3 FORMAT ( 1A ) END SUBROUTINE FOLD( RI,EP ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C This is the folding routine which allows for multiple C C scattering effects. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION EP(501),RI(-100:801) COMMON /SINGS/ RMAX,E,X0,ALPHA,RKYMAX,DKY0,RHE,THETA,XL2,XL1, + THETAP,T,X,BETA,BETAP,RKYMIN,XSCALE,RL,ZI COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI IUP = INT( RHE/ZI + 0.5 ) IF ( IUP.LT.NZ ) IUP = NZ DO 200 L1=1,NZ DO 100 L2=-100,IUP SL = EP(L1)*RI(L2) IF ( (L1+L2).LE.IUP ) RI(L1+L2)=RI(L1+L2)+SL RI(L2) = RI(L2) - SL 100 CONTINUE 200 CONTINUE RETURN END SUBROUTINE GRAPH( RI,SIGMA ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C This routine draws a graph of the calculated spectrum. It C C uses graphic routines from the NGRAPH library available from C C Nestor Zaluzec. C C C C This routine asks the user to specify the lower and upper energy C C range of the spectrum to be plotted (XLOW and XHIGH) and then C C finds the value of the spectrum in this energy range beyond the C C zero loss peak (SUMH). The user is then informed of this maximum C C value and asked to specify the upper limit of the spectrum C C intensity to set the scale of the graph to be plotted. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) COMMON /SINGS/ RMAX,E,X0,ALPHA,RKYMAX,DKY0,RHE,THETA,XL2,XL1, + THETAP,T,X,BETA,BETAP,RKYMIN,XSCALE,RL,ZI COMMON /CRTLST/ AXMIN,DX,AYMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ COMMON /LABEL/ IXTIE(12),IYTIE(12),NILX,NILY DIMENSION RI(-100:801) REAL*4 RONE,XA(801),YA(801),AXMIN,AYMIN,DX,DY,OFFX,OFFY,XL,YL, + XZ,YZ,SYA(801) CHARACTER*2 ANS DATA IXTIE /'EN','ER','GY',' L','OS','S ',' (','eV',') ',' ', + ' ',' '/ DATA IYTIE /'C ','O ','U ','N ','T ','S ',' ',' ',' ',' ', + ' ',' '/ XL = 800.0 YL = 500.0 XZ = 150.0 YZ = 250.0 NILX = 20 NILY = 7 100 WRITE ( 6,* ) ' INPUT ENERGY RANGE OF GRAPH. LOWER ENERGY (eV):' READ ( 5,* ) XLOW IF ( XLOW.LT.0D0 ) GO TO 100 IF ( XLOW.GT.RHE ) GO TO 100 120 WRITE ( 6,* ) ' HIGHER ENERGY (eV):' READ ( 5,* ) XHIGH IF ( XHIGH.LE.XLOW ) GO TO 100 IF ( XHIGH.GT.RHE ) GO TO 120 SUMH = 0.0D0 NLOW = INT(XLOW/ZI+0.5D0) NHIGH = INT(XHIGH/ZI+0.5D0) C C FIND THE MAXIMUM VALUE OF THE SPECTRUM IN THE REQUIRED ENERGY C RANGE BEYOND THE ZERO LOSS PEAK. SIGMA IS THE STANDARD DEVIATION C OF THE ZERO LOSS PEAK. C SUMH = 0.0D0 IF (XLOW.LE.3.0D0*SIGMA.OR.ABS(SIGMA).LE.1D-6) THEN RLAST = RI(NLOW) IF (NLOW.LE.0) RLAST = RI(1) DO 140 L1=NLOW+1,NHIGH NMIN=L1 IF ( RI(L1).GT.RLAST ) GO TO 160 RLAST = RI(L1) 140 CONTINUE ENDIF SUMH = RI(NLOW) GO TO 200 160 DO 170 L1=NMIN,NHIGH IF ( RI(L1).GT.SUMH ) SUMH=RI(L1) 170 CONTINUE 200 WRITE (6,1) SUMH 220 WRITE ( 6,5 ) READ ( 5,* ) YHIGH IF ( YHIGH.LT.0.0D0 ) GO TO 220 C C SET UP ARRAYS TO BE PLOTTED C ENER = XLOW - ZI NP = INT( (XHIGH-XLOW)/ZI + 0.5 ) NP = NP + 1 DO 150 L1=1,NP YA(L1)=RI(NLOW+L1-1) IF ( YA(L1).GT.YHIGH ) YA(L1) = YHIGH SYA(L1) = 0.0 ENER = ENER + ZI XA(L1) = ENER 150 CONTINUE SYA(1) = YHIGH C C PLOT GRAPH - this routine plots a graph using the data in the array C YA as the y value, and array XA as the x values between C XLOW and XHIGH. NP is the number of data values to be plotted. C IZERO = 0 IONE = 1 RONE = 1.0 CALL TERMIN CALL STRTCRT CALL ERASE CALL CRT( XA(1),YA(1),NP,IONE,IZERO,RONE,RONE,XA(1),SYA(1),NP ) C C ASK IF WANT TO WRITE SPECTRUM TO DISC. C CALL PLOT ( 10,150,0 ) CALL WRASK C C CHANGE BACK TO A VT100 TERMINAL C CALL VT100 RETURN 1 FORMAT ( ' HIGHEST INTENSITY OF THE SPECTRUM IS ',F9.1 ) 2 FORMAT ( 12A4 ) 3 FORMAT ( 6A4 ) 4 FORMAT ( ' INPUT REQUIRED INTENSITY RANGE OF GRAPH:'/ + ' LOWEST INTENSITY = ' ) 5 FORMAT ( ' HIGHEST INTENSITY = ' ) 6 FORMAT ( 1A1 ) END SUBROUTINE ALMINA( MDM,IDE,ER,EI,TITLE ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C Dielectric data for alumina. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION IDE(0:20),ER(0:20),EI(0:20) CHARACTER*4 TITLE(18) DIMENSION IRAY(20), RAYR(20), RAYI(20) COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI DATA (IRAY(J),J=1,20)/5,6,7,8,9,10,12,14,15,16,18,20,24, + 28,32,0,0,0,0,0/ DATA (RAYR(J),J=1,20)/2.8,3,3.5,4,5.4,6.3,1.8,1,0.3,0.25, + 0.1,0.2,0.3,0.5,0.6,0,0,0,0,0/ DATA (RAYI(J),J=1,20)/0.001,0.01,0.03,0.1,0.7,3,4,3,2.6,2.4, + 2,1.5,0.6,0.5,0.3,0,0,0,0,0/ TITLE((MDM-1)*6+1) = 'ALUM' TITLE((MDM-1)*6+2) = 'INA ' WRITE( 6,1 ) MDM NPOINT = 15 DO 100 L1=1,20 IDE(L1) = IRAY(L1) ER(L1) = RAYR(L1) EI(L1) = RAYI(L1) 100 CONTINUE RETURN 1 FORMAT ( ' MEDIUM ',I2,' IS ALUMINA' ) END SUBROUTINE ALMNM( MDM,IDE,ER,EI,TITLE ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C Dielectric data for aluminium. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION IDE(0:20),ER(0:20),EI(0:20) CHARACTER*4 TITLE(18) DIMENSION IRAY(20), RAYR(20), RAYI(20) COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI DATA (IRAY(J),J=1,20)/5,6,7,8,9,10,12,14,15,16,18,20,24,28,32, + 0,0,0,0,0/ DATA (RAYR(J),J=1,20)/-7.65,-5.08,-3.5,-2.46,-1.74,-1.23,-0.55, + -0.14,0.004,0.124,0.31,0.44,0.61,0.8,1,0,0,0,0,0/ DATA (RAYI(J),J=1,20)/1.73,1.01,0.64,0.43,0.3,0.22,0.13,0.08, + 0.07,0.05,0.04,0.03,0.02,0.018,0.016,0,0,0,0,0/ TITLE((MDM-1)*6+1) = 'ALUM' TITLE((MDM-1)*6+2) = 'INIU' TITLE((MDM-1)*6+3) = 'M ' WRITE( 6,1 ) MDM NPOINT = 15 DO 100 L1=1,20 IDE(L1) = IRAY(L1) ER(L1) = RAYR(L1) EI(L1) = RAYI(L1) 100 CONTINUE RETURN 1 FORMAT ( ' MEDIUM ',I2,' IS ALUMINIUM' ) END C CCCCCCCCCCCCCCCC C C C C C CARBON C C C C C CCCCCCCCCCCCCCCC SUBROUTINE CARBON( MDM,IDE,ER,EI,TITLE ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C Dielectric data for carbon. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION IDE(0:20),ER(0:20),EI(0:20) CHARACTER*4 TITLE(18) DIMENSION IRAY(20), RAYR(20), RAYI(20) COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI DATA (IRAY(J),J=1,20)/5,6,7,8,9,10,12,14,15,16,18,20,24,28,32, + 0,0,0,0,0/ DATA (RAYR(J),J=1,20)/1.15,1.1,1.5,2.2,2.9,3,2.2,1,0.5,0,0.1, + 0.3,0.5,0.6,0.7,0,0,0,0,0/ DATA (RAYI(J),J=1,20)/3,1.5,0.7,0.4,1,2,5,5,3,2,1.5,1,0.4, + 0.3,0.25,0,0,0,0,0/ TITLE((MDM-1)*6+1) = 'CARB' TITLE((MDM-1)*6+2) = 'ON ' WRITE( 6,1 ) MDM NPOINT = 15 DO 100 L1=1,20 IDE(L1) = IRAY(L1) ER(L1) = RAYR(L1) EI(L1) = RAYI(L1) 100 CONTINUE RETURN 1 FORMAT ( ' MEDIUM ',I2,' IS CARBON' ) END SUBROUTINE DIAMD( MDM,IDE,ER,EI,TITLE ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C Dielectric data for diamond. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION IDE(0:20),ER(0:20),EI(0:20) CHARACTER*4 TITLE(18) DIMENSION IRAY(20), RAYR(20), RAYI(20) COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI DATA (IRAY(J),J=1,20)/2,4,6,8,10,12,14,16,20,24,28,32,0,0,0,0, + 0,0,0,0/ DATA (RAYR(J),J=1,20)/5.5,6,7,10,9,0.001,-5,-3,-1,-0.15,0.04, + 0.2,0,0,0,0,0,0,0,0/ DATA (RAYI(J),J=1,20)/0.001,0.001,0.001,4,9,17,4.5,2,0.9,0.6, + 0.3,0.05,0,0,0,0,0,0,0,0/ TITLE((MDM-1)*6+1) = 'DIAM' TITLE((MDM-1)*6+2) = 'OND ' WRITE( 6,1 ) MDM NPOINT = 12 DO 100 L1=1,20 IDE(L1) = IRAY(L1) ER(L1) = RAYR(L1) EI(L1) = RAYI(L1) 100 CONTINUE RETURN 1 FORMAT ( ' MEDIUM ',I2,' IS DIAMOND' ) END SUBROUTINE DIAMDP( MDM,IDE,ER,EI,TITLE ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C Dielectric data for diamond. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION IDE(0:20),ER(0:20),EI(0:20) CHARACTER*4 TITLE(18) DIMENSION IRAY(20), RAYR(20), RAYI(20) COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI DATA (IRAY(J),J=1,20)/2,4,6,8,10,12,14,16,20,24,28,32,0,0,0,0, + 0,0,0,0/ DATA (RAYR(J),J=1,20)/5.5,6,7,10,9,0,-5,-3,-1,-0.15,0.04,0.2, + 0,0,0,0,0,0,0,0/ DATA (RAYI(J),J=1,20)/0,0,0,4,9,17,4.5,2,0.9,0.6,0.3,0.05,0,0, + 0,0,0,0,0,0/ TITLE((MDM-1)*6+1) = 'DIAM' TITLE((MDM-1)*6+2) = 'OND"' WRITE( 6,1 ) MDM NPOINT = 12 DO 100 L1=1,20 IDE(L1) = IRAY(L1) ER(L1) = RAYR(L1) EI(L1) = RAYI(L1) 100 CONTINUE RETURN 1 FORMAT ( ' MEDIUM ',I2,' IS DIAMOND"' ) END SUBROUTINE GAAS( MDM,IDE,ER,EI,TITLE ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C Dielectric data for gallium arsenide. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION IDE(0:20),ER(0:20),EI(0:20) CHARACTER*4 TITLE(18) DIMENSION IRAY(20), RAYR(20), RAYI(20) COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI DATA (IRAY(J),J=1,20)/2,3,4,5,6,7,8,9,10,12,14,15,17,19,20, + 0,0,0,0,0/ DATA (RAYR(J),J=1,20)/15.0,15.0,9,-6,-2.3,-2.5,-1.2,-0.6,-0.2, + 0.1,0.3,0.4,0.6,0.8,0.9,0,0,0,0,0/ DATA (RAYI(J),J=1,20)/1,17,13,14,6,5,3,2,1.8,1.4,1,0.8,0.6, + 0.4,0.4,0,0,0,0,0/ TITLE((MDM-1)*6+1) = 'GaAs' WRITE( 6,1 ) MDM NPOINT = 15 DO 100 L1=1,20 IDE(L1) = IRAY(L1) ER(L1) = RAYR(L1) EI(L1) = RAYI(L1) 100 CONTINUE RETURN 1 FORMAT ( ' MEDIUM ',I2,' IS GaAs' ) END SUBROUTINE GAP( MDM,IDE,ER,EI,TITLE ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C Dielectric data for gallium phosphide. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION IDE(0:20),ER(0:20),EI(0:20) CHARACTER*4 TITLE(18) DIMENSION IRAY(20), RAYR(20), RAYI(20) COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI DATA (IRAY(J),J=1,20)/1,2,3,4,5,6,7,8,10,13,15,17,19,21,23, + 0,0,0,0,0/ DATA (RAYR(J),J=1,20)/9.6,11.6,19.2,5.5,-0.2,-3.8,-3.8,-3, + -1.4,-0.56,-0.16,0.16,0.36,0.46,0.6,0,0,0,0,0/ DATA (RAYI(J),J=1,20)/0,0,4.5,13.8,15,6.3,4.6,2.3,1.6,0.56, + 0.4,0.24,0.12,0.16,0.16,0,0,0,0,0/ TITLE((MDM-1)*6+1) = 'GaP ' WRITE( 6,1 ) MDM NPOINT = 15 DO 100 L1=1,20 IDE(L1) = IRAY(L1) ER(L1) = RAYR(L1) EI(L1) = RAYI(L1) 100 CONTINUE RETURN 1 FORMAT ( ' MEDIUM ',I2,' IS GaP' ) END SUBROUTINE MGO( MDM,IDE,ER,EI,TITLE ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C Dielectric data for magnesium oxide. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION IDE(0:20),ER(0:20),EI(0:20) CHARACTER*4 TITLE(18) DIMENSION IRAY(20), RAYR(20), RAYI(20) COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI DATA (IRAY(J),J=1,20)/ 5,6,7,8,9,10,11,12,14,15,16,18,20,22,24, + 28,32,0,0,0/ DATA (RAYR(J),J=1,20)/ 3.4,3.9,4.8,4.2,4.1,3.4,1.2,1.8,0.1,0.6, + 0.9,-0.25,-0.25,0.05,0.3,0.5,0.8,0.0,0.0,0.0/ DATA (RAYI(J),J=1,20)/0.001,0.001,0.001,2.4,2.8,3.4,4.2,2.4,2.0, + 1.9,1.9,1.7,0.8,0.3,0.3,0.01,0.005,0.0,0.0,0.0/ TITLE((MDM-1)*6+1) = 'MgO ' WRITE( 6,1 ) MDM NPOINT = 17 DO 100 L1=1,20 IDE(L1) = IRAY(L1) ER(L1) = RAYR(L1) EI(L1) = RAYI(L1) 100 CONTINUE RETURN 1 FORMAT ( ' MEDIUM ',I2,' IS MgO' ) END SUBROUTINE SILCAP( MDM,IDE,ER,EI,TITLE ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C Dielectric data for silica. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION IDE(0:20),ER(0:20),EI(0:20) CHARACTER*4 TITLE(18) DIMENSION IRAY(20), RAYR(20), RAYI(20) COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI DATA (IRAY(J),J=1,20)/5,6,7,8,9,10,12,14,15,16,18,20,24,28,32, + 0,0,0,0,0/ DATA (RAYR(J),J=1,20)/2.3,2.4,2.5,3,4,6.5,1.2,1.4,0.7,0.9,0.2, + 0.2,0.2,0.2,0.3,0,0,0,0,0/ DATA (RAYI(J),J=1,20)/0.01,0.02,0.04,0.06,0.2,2.8,3,2.6,2, + 1.9,1.5,1,0.5,0.3,0.2,0,0,0,0,0/ TITLE((MDM-1)*6+1) = 'SILI' TITLE((MDM-1)*6+2) = 'CA" ' WRITE( 6,1 ) MDM NPOINT = 15 DO 100 L1=1,20 IDE(L1) = IRAY(L1) ER(L1) = RAYR(L1) EI(L1) = RAYI(L1) 100 CONTINUE RETURN 1 FORMAT ( ' MEDIUM ',I2,' IS SILICA"' ) END SUBROUTINE SILICA( MDM,IDE,ER,EI,TITLE ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C Dielectric data for silica. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION IDE(0:20),ER(0:20),EI(0:20) CHARACTER*4 TITLE(18) DIMENSION IRAY(20), RAYR(20), RAYI(20) COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI DATA (IRAY(J),J=1,20)/5,6,7,8,9,10,12,14,15,16,18,20,24,28,32, + 0,0,0,0,0/ DATA (RAYR(J),J=1,20)/2.3,2.4,2.5,2.8,3.5,4.9,1.4,1.2,1.1,1.2, + 0.8,0.6,0.65,0.8,0.9,0,0,0,0,0/ DATA (RAYI(J),J=1,20)/0.01,0.02,0.04,0.06,0.12,2.8,2.1,1.6,1.4, + 1.3,1.2,0.7,0.45,0.2,0.15,0,0,0,0,0/ TITLE((MDM-1)*6+1) = 'SILI' TITLE((MDM-1)*6+2) = 'CA ' WRITE( 6,1 ) MDM NPOINT = 15 DO 100 L1=1,20 IDE(L1) = IRAY(L1) ER(L1) = RAYR(L1) EI(L1) = RAYI(L1) 100 CONTINUE RETURN 1 FORMAT ( ' MEDIUM ',I2,' IS SILICA' ) END SUBROUTINE SILICN( MDM,IDE,ER,EI,TITLE ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C Dielectric data for silicon. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION IDE(0:20),ER(0:20),EI(0:20) CHARACTER*4 TITLE(18) DIMENSION IRAY(20), RAYR(20), RAYI(20) COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI DATA (IRAY(J),J=1,20)/2,3,4,5,6,7,8,9,10,12,14,15,16,18,20, + 24,28,32,0,0/ DATA (RAYR(J),J=1,20)/15,45,-2,-9,-5,-4,-3,-2,-1.5,-1.1,-0.7, + -0.3,-0.1,0.1,0.2,0.5,0.6,0.7,0,0/ DATA (RAYI(J),J=1,20)/0,1,40,10,7,5,3.5,2.3,1.5,0.9,0.5, + 0.3,0.25,0.2,0.15,0.1,0.05,0.02,0,0/ TITLE((MDM-1)*6+1) = 'SILI' TITLE((MDM-1)*6+2) = 'CON ' WRITE( 6,1 ) MDM NPOINT = 18 DO 100 L1=1,20 IDE(L1) = IRAY(L1) ER(L1) = RAYR(L1) EI(L1) = RAYI(L1) 100 CONTINUE RETURN 1 FORMAT ( ' MEDIUM ',I2,' IS SILICON' ) END SUBROUTINE VACUUM( MDM,E2R,E2I,TITLE ) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION E2R(501),E2I(501) CHARACTER*4 TITLE(18) TITLE((MDM-1)*6+1) = 'VACU' TITLE((MDM-1)*6+2) = 'UM ' WRITE( 6,1) MDM DO 100 L1=1,501 E2R(L1) = 1.0 E2I(L1) = 0.00001 100 CONTINUE RETURN 1 FORMAT ( ' MEDIUM ',I2,' IS A VACUUM' ) END SUBROUTINE GAPP( MDM,IDE,ER,EI,TITLE ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C A slightly different version of dielectric data for gallium C C phosphide. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT REAL*8 (A-H,O-Z) DIMENSION IDE(0:20),ER(0:20),EI(0:20) CHARACTER*4 TITLE(18) DIMENSION IRAY(20), RAYR(20), RAYI(20) COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI DATA (IRAY(J),J=1,20)/3,4,5,6,7,8,9,10,11,12,13,14,15,16, + 17,18,19,20,22,24/ DATA (RAYR(J),J=1,20)/14.75,11.7,1,-2.25,-2.5,-1.5,0,-0.05, + -0.25,-0.2,-0.15,-0.05,0,0.2,0.3,0.5,0.65,0.85,0.82,0.8/ DATA (RAYI(J),J=1,20)/0.75,11.5,22,7.7,5.5,3.3,3,3.5,2.25, + 1.75,1.5,1.25,1,0.75,0.6,0.5,0.65,0.5,0.5,0.5/ TITLE((MDM-1)*6+1) = 'GaP"' WRITE( 6,1 ) MDM NPOINT = 20 DO 100 L1=1,20 IDE(L1) = IRAY(L1) ER(L1) = RAYR(L1) EI(L1) = RAYI(L1) 100 CONTINUE RETURN 1 FORMAT ( ' MEDIUM ',I2,' IS GaP"' ) END SUBROUTINE WRASK CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C SUBROUTINE TO WRITE SPECTRUM TO A FILE ON DISC IN THE FORMAT C C REQUIRED TO BE READ BY THE PROGRAM NELS. C C PARAMETERS USED IN THE CALCULATIONS ARE WRITTEN ON THE TOP C C OF THE FILE. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT DOUBLE PRECISION (A-H,O-Z) CHARACTER*20 FN1,FN2,FN3,CZERO CHARACTER*1 ANS CHARACTER*4 TITLE(18) COMMON FN1,FN2,FN3 COMMON /RAYS/ ER(0:20),EI(0:20),E1R(501),E1I(501),E2R(501), + E2I(501),EP(501),RI(-100:801),IDE(0:20) COMMON /SINGS/ RMAX,E,X0,ALPHA,RKYMAX,DKY0,RHE,THETA,XL2,XL1, + THETAP,T,X,BETA,BETAP,RKYMIN,XSCALE,RL,ZI COMMON /INNS/ NPOINT,NS,NMAT,NZ,NRI COMMON /GRPH/ TITLE COMMON /CRACT/ C,G,ANS COMMON /WRTFIL/ IOPEN DIMENSION XS(11) 105 WRITE(6,*) ' WRITE SPECTRUM TO FILE IN NELS DATA FORMAT (Y/N)?' READ(5,2) ANS IF ( ANS.EQ.'N' ) RETURN IF ( ANS.NE.'Y' ) GO TO 105 CZERO = ' ' IF ( IOPEN.EQ.1 ) THEN 100 WRITE(6,1) FN3 WRITE(6,*) ' SHALL I WRITE THE SPECTRUM TO THIS FILE (Y/N)?' READ(5,2) ANS IF ( ANS.NE.'Y'.AND.ANS.NE.'N' ) GO TO 100 ENDIF IF ( ANS.EQ.'N'.OR.IOPEN.EQ.0 ) THEN CLOSE(14) FN3 = CZERO 110 WRITE(6,*) ' INPUT A NEW FILENAME FOR THE SPECTRUM:' READ(5,7) FN3 IF ( FN3.EQ.CZERO ) GO TO 110 OPEN (UNIT=14,FILE=FN3,FORM='FORMATTED',ACCESS='SEQUENTIAL', + ERR=120,IOSTAT=IOST,STATUS='NEW') IOPEN = 1 GO TO 130 120 WRITE(6,3) FN3,IOST GO TO 110 130 CONTINUE ENDIF C C WRITE SPECTRUM TO DISC C WRITE(14,4) (TITLE(J),J=1,18) SA = RKYMAX/E/2.0 WRITE(14,5) SA,THETA,XL1,XL2,BETA NUP = RHE/ZI NOPTS = NUP + 101 WRITE(14,6) NOPTS, (RI(J),J=-99,-91) NOPTS = -90 140 IF ( NOPTS + 10 .GE.NUP ) GO TO 150 WRITE(14,8) (RI(J),J=NOPTS,NOPTS+9) NOPTS = NOPTS + 10 GO TO 140 C C WRITE REST OF SPECTRUM AND PARAMETERS TO DISK C 150 IP = 1 IF ( NOPTS.LE.NUP ) THEN DO 160 L1=NOPTS,NUP XS(IP) = RI(L1) IP = IP + 1 160 CONTINUE DO 165 L1=IP,10 XS(L1)=0.0 165 CONTINUE WRITE(14,8) (XS(J),J=1,10) ENDIF C C OFFSET ENERGY OF FIRST CHANNEL C IP = 1 XS(IP) = -100.0*ZI C C ENERGY PER CHANNEL (eV) C IP = IP + 1 XS(IP) = ZI IP = IP + 1 C C ACCELERATING VOLTAGE (keV) C XS(IP) = E IP = IP + 1 C C INCIDENCE BEAM DIVERGENCE C XS(IP) = 0.0 IP = IP + 1 C C SCATTERING ANGLE (mrad) C XS(IP) = SA IP = IP + 1 C C LIVE TIME C XS(IP) = 0.0 IP = IP + 1 C C DEAD TIME C XS(IP) = 0.0 IP = IP + 1 C C INCIDENT BEAM CURRENT (nA) C XS(IP) = 0.0 IP = IP + 1 C C PROBE DIAMETER (nm) C XS(IP) = 0.0 IP = IP + 1 C C SPECIMEN THICKNESS ALONG BEAM DIRECTION (nm) C XS(IP) = (XL1 - XL2)/10.0 IF ( ABS(THETA).LT.1D-10 ) XS(IP) = XL2/10.0 IP = IP + 1 WRITE(14,9) (XS(J),J=1,10) RETURN 1 FORMAT( ' OUTPUT FILE: ',A20,' ALREADY OPEN' ) 2 FORMAT( 1A1 ) 3 FORMAT( ' ERROR OCCURRED WHILST TRYING TO OPEN FILE ',1A20/ + ' IOSTAT=',I4 ) 4 FORMAT( 18A4 ) 5 FORMAT( ' SCATT ANG='F8.2,'mrad THETA=',F8.6,'rad XL1='F7.1, + ' XL2=',F7.1,' BETA=',F6.4 ) 6 FORMAT( I12,9(F12.0) ) 7 FORMAT( 1A20 ) 8 FORMAT( 10(F12.0) ) 9 FORMAT( 10(F12.4) ) END SUBROUTINE OFILE CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C SUBROUTINE TO ASK FOR THE NAME OF AN OUTPUT FILE AND OPEN IT C C ON UNIT 15 FOR OUTPUT OF SPECTRA AND PROBABILITY DATA. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT DOUBLE PRECISION (A-H,O-Z) CHARACTER*20 CZERO,FN CZERO = ' ' FN = CZERO 100 WRITE(6,*) ' INPUT OUTPUT FILENAME FOR CALCULATED SPECTRA:' READ(5,1) FN IF ( FN.EQ.CZERO ) GO TO 100 OPEN (UNIT=15,FILE=FN,FORM='FORMATTED',ACCESS='SEQUENTIAL', + ERR=120,IOSTAT=IOST,STATUS='NEW') GO TO 130 120 WRITE(6,2) FN,IOST GO TO 100 130 CONTINUE RETURN 1 FORMAT ( 1A20 ) 2 FORMAT ( ' ERROR WHILST TRYING TO OPEN FILE',1A20,' IOSTAT=',I6 ) END