Title :NELS- Version 21 Keywords :EELS Computer :DEC VAX-11/730 - 11/785 or higher series Operating System :VMS Programming Language :Fortran IV Hardware Requirements :Tektronics 4010 Series Graphic Terminal or Emulator Author(s) :Nestor J. Zaluzec Correspondence Address :Electron Microscopy Center, Materials Science, Bldg 212 :Argonne Nat. Lab, Argonne, IL. 60439 USA Abstract: NELS is a generic data analysis computer program for reduction of EELS spectra. The program performs the functions of spectral calibration, KLM edge identification, background fitting & subtraction, edge integration, digital filtering, and general spectral maniputlation. The input spectra are assumed to be in the EMMPDL spectral file format (see RWEMMPDL.ABS). The program can handle simultaneously three spectra of upto 4096 channels in length. A shorter and less elaborate version is available for smaller DEC computers running RT-11. The program NELS requires a computer terminal which responds to the graphics protocol of Tektronics 4010/4014 series terminal (see NGRAPH.ABS) and an atomic inner shell energy level data file (see ATMLVL.ABS). Title :NELS- Version 21 Keywords :EELS Computer :DEC VAX-11/730 - 11/785 or higher series Operating System :VMS Programming Language :Fortran IV Hardware Requirements :Tektronics 4010 Series Graphic Terminal or Emulator Author(s) :Nestor J. Zaluzec Correspondence Address :Electron Microscopy Center, Materials Science, Bldg 212 :Argonne Nat. Lab, Argonne, IL. 60439 USA Documentation: The following is the Compilation and Linking Procedure for NELS V21 for a VAX- based system $fortran/extend nels21 $fortran/extend ngraph $link nels20,ngraph $rename nels20.exe nels.exe $run nels NOTE: you must have current versions of the ATMLVL.TBL file and for the program to implement all options. You should also create the terminal data file TERMIN.DAT for the graphics utilities. for further information on these see EMMPDL entries NGRAPH, ATMLVL. A preliminary set of instructions for the NEDS/NELS series of data (spectral analysis programs follows. Preliminary Instructions for NELS VERSION 21 ---------------------------------------------- The program is started entering the normal RT-11 operating mode. When the (.) appears indicating RT-11 is ready enter the following .RUN NELS or .RUN DY1:NELS where means to enter a carriage RETURN. For Vax-based systems you need only enter RUN NELS after your system prompt. The RT-11 device specification DY1: is only necessary if your copy of NELS is stored on DISK DRIVE 1 instead of the system DISK (DY0:) as most software is normally stored. You should also insure that the file ATMLVL.TBL exists and is on the system DISK. You should refer to the EMMPDL files on this subject for further details. Before running the program you should first inspect the file DK:TERMIN.DAT by typing it out on the crt display. This file is read by NELS and determines the characteristics of your terminal relative .to this program. If the characteristics of your terminal (type and/or baud rate) are different than those listed, then use the normal .VAX/RT-11 EDITOR to change the terminal type (i.e. 0,1,2,3 as defined in the file), and the baud rate. NOTE: donot erase any text from this file as it important in the field specifications. A listing of this file is given below: GRAPHICS TERMINAL TYPE:8 BAUD RATE FOR GRAPHICS:9600 HARDCOPY DEVICE TYPE:5 C -------------------------------------------------------- C TERMINAL DEFINITIONS: C --------------------- C -3 = ASCII TERMINAL: NO GRAPHICS, NO VT-100 compatibility <32COLUMNS C -2 = ASCII TERMINAL: NO GRAPHICS, NO VT-100 compatibility C -1 = ANSI X3.4 ASCII TERMINAL: NO GRAPHICS, VT-100 compatible C 0 = TEKTRONICS 4010-1,4014-1 (WITH HARDWARE CURSOR) C 1 = TEKTRONICS 4006 ,4010 (NO HARDWARE CURSOR) C 2 = LEAR SIEGLER ADM-3,ADM-5 (WITH RETROGRAPHICS RG-512) C 3 = PERITEK VCG-512 BIT MAP COLOR GRAPHICS C 4 = TEKTRONICS 4027 C 5 = HP 7470A/7475A HARDCOPY PLOTTER C 6 = LA100/LA50 HARDCOPY GRAPHICS C 7 = INTECOLOR VT-100 with 4010-1,4014-1 mode C 8 = INTECOLOR VT-100 with 4027 mode C 9 = DEC VT-200 with 4010-1,4014-1 mode C 10 = TEKTRONICS 4105/4107 with VT-100 compatibility mode C 11 = PLESSEY VT-100/TEK 4010-1 C 12 = ESPIRIT VT-100/TEK 4010-1 C 13 = TEKTRONICS 4010-1 HARDCOPY VIDEO TO VERSATIC C 14 = TEKTRONICS 4695 Color Copier C -------------------------------------------------------- C BAUD RATE RANGE FOR GRAPHICS (110 - 19200) C -------------------------------------------------------- C C The current terminal configurations is by the above an INTECOLOR VT-100 with Tektronics 4027 graphics, operating at 9600 BAUD. The program begins by initializing the display, printing out the current version number, date etc. and the hardcopy output file name. For VAX based systems this is the default output device FOR008.DAT. The program next types out the general Menu of available commands the operator need only type in the corresponding two character code to select an option. Hopefully, most of these will be selfexplainitory. A list of the current Menu is given below GR = GRAPH MODIFICATION OF DISPLAY WI = WINDOW SET AND EXPAND TO FILL DISPLAY CU = CURSOR VALUE SM = SMOOTH DATA CA = CALIBRATE SPECTRUM NO = NORMALIZE (REMOVE) GAIN CHANGES BG = BACKGROUND FITTING IN = INTEGRATE EDGE DF = DIGITAL FILTER OF THE SPECTRUM ME = MEMORY AND DATA MANIPULATION RE = READ DATA FILE FROM DISK WR = WRITE DATA TO FILE ON DISK OR = ORIGINAL DATA REDRAWN EX = EXIT PROGRAM THE FOLLOWING MODES ARE AVAILABLE IN THE "ME" OPTION CO = COPY FROM MEM 1 TO 2 OR 3 SS = SPECTRUM-SPECTRUM MATH +,-,X,/ MEM 1 BY MEM 2 SC = SPECTRUM-CONSTANT MATH " " " " MEM 1 BY CONSTANT LO = LOG OF DISPLAY AL = ANTILOG OF DISPLAY OV = OVERLAY MEM1/MEM2 ETC. ZE = SETS ALL VALUES <0. IN CURRENT MEMORY = TO 0. DE = DELETE PORTION OF CURRENT MEMORY EQ = EQUALIZE (SCALE) 2 MEMORIES AT CURSOR LOCATION After listing the options the program automatically requests you to supply a EMMPDL DATA FILE for input. Use the demonstration file NIO.001 if you donot have a suitable spectrum available. Please note that this program only reads data in the EMMPDL file format, if your data file is of a different format (ie. EDAX, KEVEX, ORTEC, TRACOR, PGT, LINK etc.) this program will not accept the data files as is. The program NTRAN will translate some of this data to EMMPDL format and vice-versa and you must run it first. The subroutine RPDL.FOR documents the expected data file format for those interested in writing their own translation routine, however, I would appreciate comments and/or information concerning other data formats inorder that I can keep NTRAN up to date for all users. At this point the program enters its general command sequence the usual series of steps in an analysis would be something like this: 1.) RECALL SPECTRUM (if present data is not suitable) 2.) CALIBRATE, and MODIFY SPECTRUM USING DISPLAY MANIPULATION OPTIONS 3.) NORMALIZE (if gain changes are present between edges) 4.) FIT BACKGROUND CURVE FOR 1ST ELEMENT (lowest energy first!) 5.) INTEGRATE 1ST ELEMENT 6.) FIT BACKGROUND CURVE FOR NEXT ELEMENT 7.) INTEGRATE (repeat lines 5 - 7 for each element) 8.) EXIT PROGRAM General Remarks --------------- This version of NELS exclusively uses memory 1 for data and memory 2 for calculations and/or backgrounds. The maximum spectral size that the current program can support is 4096 channels in 3 memories on a Vax system, or with appropriate modification 1024 channels in 3 memories on a RT-11 system. Remember that some routines modify or even delete portions of the spectrum so always store the spectra on a DISK FILE before analyzing it. GRPH ---- This is the general display control routine which allows the operator to redraw the current display. This should be used when the display gets too messy due to multiple operations or you just wish to erase the entire screen and start over. Its options are few, but must be entered in the following order OVERLAY(Y-N),Y-MAG,SYMBOL(-9 to +1). Here the OVERLAY options allows the user to superimpose ontop of the current display a modification of the present display . Thus by specifying a Y-MAG the user can expand the display to look at a small feature of interest. SYMBOL refers to the various plotting symbols +1= line, 0= dots, etc. (-9 = solid fill and is good for photography). A typical answer sequence would therefore be: Y,3.,+1 or N,1.,-9 Remember the commas here are essential to seperate command inputs!!. VALUE OF CURSOR EMail: Zaluzec at ANLMST on BITNET ------------------------- The following is NIO test data referenced in the instructions ------------------------- NICKEL OXIDE (NIO) THIN & CALIBRATED 1020., 93., 51., 75., 60., 79., 62., 64., 141., 66., 64., 76., 54., 99., 58., 63., 98., 69., 62., 84., 132., 83., 61., 75., 68., 75., 61., 49., 76., 76., 73., 68., 66., 78., 58., 78., 81., 66., 67., 78., 82., 72., 94., 64., 93., 84., 77., 94., 95., 78., 61., 70., 78., 91., 69., 76., 95., 63., 79., 96., 73., 82., 63., 68., 101., 84., 75., 68., 82., 95., 82., 95., 67., 94., 68., 67., 83., 66., 61., 89., 96., 114., 82., 88., 79., 79., 112., 100., 98., 111., 89., 86., 94., 107., 82., 104., 103., 111., 88., 90., 90., 96., 81., 96., 122., 128., 90., 115., 113., 147., 119., 134., 135., 151., 161., 155., 161., 197., 197., 243., 279., 267., 330., 371., 415., 578., 615., 719., 958., 1202., 1532., 2034., 2658., 3250., 4374., 5908., 7733., 10446., 13964., 19649., 27850., 40443., 70906.,287377.,907273.,999999.,999999.,378770., 66816., 22640., 19922., 18842., 18193., 17625., 17627., 18512., 20334., 22886., 24198., 27896., 32936., 38326., 42489., 45744., 45637., 45694., 42882., 39301., 35152., 31583., 29351., 28482., 28329., 28231., 27694., 27198., 27029., 25240., 23441., 21900., 21761., 20889., 19835., 19086., 18714., 18160., 17165., 15978., 15562., 14447., 13456., 12783., 12715., 12453., 11348., 11020., 10826., 10562., 10939., 11753., 12310., 11939., 12096., 11224., 10748., 10746., 10011., 9657., 9138., 8772., 7979., 7886., 7974., 7469., 7350., 6893., 6750., 6689., 6406., 6630., 6251., 6156., 5988., 5570., 5426., 5285., 5220., 4784., 4774., 4711., 4548., 4602., 4369., 4196., 4170., 4039., 3939., 3861., 3809., 3644., 3583., 3468., 3446., 3326., 3296., 3245., 2997., 3223., 2869., 3161., 2924., 2977., 2852., 2748., 2856., 2681., 2573., 2367., 2314., 2378., 2200., 2131., 2149., 2140., 1977., 1892., 2087., 2029., 2534., 1768., 1661., 1866., 1781., 1756., 1720., 1675., 1608., 1617., 1508., 1466., 1435., 1472., 1510., 1398., 1368., 1329., 1337., 1252., 1453., 1239., 1300., 1183., 1143., 1201., 1099., 1156., 1118., 1022., 1108., 982., 1125., 1025., 982., 944., 1050., 975., 928., 838., 823., 760., 765., 886., 813., 747., 763., 766., 748., 730., 757., 722., 649., 615., 738., 715., 693., 622., 680., 737., 584., 641., 564., 655., 689., 604., 592., 528., 607., 604., 522., 603., 478., 595., 490., 567., 545., 415., 482., 465., 426., 425., 431., 457., 421., 452., 404., 426., 417., 402., 453., 351., 402., 390., 414., 336., 387., 393., 393., 362., 342., 320., 408., 344., 338., 376., 369., 325., 324., 388., 334., 286., 368., 313., 294., 276., 287., 278., 309., 298., 376., 276., 272., 285., 289., 303., 285., 262., 246., 241., 282., 276., 233., 293., 282., 240., 239., 249., 221., 281., 256., 271., 149., 233., 222., 257., 216., 267., 218., 223., 247., 257., 232., 247., 219., 218., 212., 155., 140., 210., 223., 241., 189., 202., 245., 221., 165., 251., 175., 216., 194., 149., 171., 190., 164., 186., 216., 196., 185., 224., 150., 187., 206., 189., 133., 487., 571., 646., 693., 687., 1008., 1647., 2813., 6044., 40099.,277728.,277409.,272166.,270886.,266481.,263601.,261307., 259766.,256997.,255017.,252812.,251045.,248411.,244405.,242740.,238405.,237245., 234104.,232296.,229025.,225395.,224023.,221182.,218490.,215402.,213702.,211487., 207637.,207735.,204466.,202871.,202065.,199026.,196004.,195179.,193225.,189142., 190629.,185958.,183435.,182445.,181003.,177714.,177066.,176459.,173183.,172498., 168982.,168526.,166561.,166676.,163205.,161700.,159002.,157138.,156360.,155075., 153980.,152122.,150120.,148595.,146894.,146346.,144359.,143181.,140798.,140292., 140132.,135705.,137252.,134358.,133585.,132797.,130784.,129705.,129194.,127651., 125425.,123154.,123534.,121702.,120843.,120157.,118150.,118032.,116287.,115672., 115285.,113973.,113388.,111495.,110880.,109568.,109774.,106461.,106323.,104967., 104643.,103502.,104092.,102644.,104190.,105582.,106395.,109830.,114841.,129101., 137277.,142550.,149411.,158227.,168096.,174291.,172829.,166554.,157611.,148944., 144233.,139998.,136330.,131842.,129541.,129004.,128481.,129033.,130659.,132068., 134279.,139090.,140299.,141070.,139709.,137719.,134209.,131433.,128674.,126290., 124120.,123793.,122970.,121338.,120432.,119968.,117555.,116829.,115500.,113855., 112896.,111465.,110772.,110957.,112426.,111254.,111707.,111340.,111586.,111292., 110646.,109792.,108029.,107390.,105776.,105663.,103972.,103007.,101590.,100175., 99562., 98693., 98671., 98117., 97191., 97018., 97311., 96194., 95550., 95144., 93372., 93893., 93201., 92543., 92113., 90136., 89461., 88138., 87355., 85956., 85598., 86478., 84765., 82918., 83029., 82903., 82383., 80994., 80362., 80657., 79170., 78562., 77901., 77659., 76640., 76562., 74963., 74868., 74225., 73558., 72662., 72417., 71988., 72266., 70731., 70142., 70512., 69499., 70128., 69245., 69133., 69000., 68097., 67476., 68128., 67207., 66850., 66931., 65032., 65035., 64151., 64225., 64499., 63902., 64004., 62786., 62028., 61914., 61245., 60335., 60764., 59707., 59732., 59747., 58824., 57525., 58046., 57436., 57468., 56501., 55913., 55589., 54643., 54997., 54668., 54088., 53744., 53649., 52389., 51945., 51164., 51927., 50908., 50968., 50256., 49614., 49768., 49726., 49569., 48483., 48473., 48055., 48228., 47425., 47197., 46796., 46509., 46731., 45584., 45705., 44891., 45651., 45364., 44624., 43876., 44809., 44017., 43414., 44280., 43760., 43098., 42721., 41847., 42834., 42514., 42229., 41980., 41429., 41665., 41478., 40280., 39655., 40277., 39400., 39810., 39155., 38942., 38524., 38987., 38578., 38149., 38412., 38098., 37363., 36986., 36513., 37225., 36471., 36750., 35670., 35654., 35111., 34972., 34597., 35546., 34397., 33960., 33984., 34213., 33654., 32605., 33079., 32964., 32997., 32004., 31906., 32380., 31624., 32007., 31770., 31195., 30987., 31307., 30791., 30367., 30066., 29673., 30650., 29888., 29718., 29515., 29614., 29248., 29056., 29206., 28654., 29277., 29832., 30293., 31739., 33143., 35386., 38005., 42283., 48606., 56194., 69613., 79228., 78678., 70758., 61475., 55560., 48418., 46756., 46864., 46747., 47284., 48928., 52790., 55932., 58923., 59484., 57711., 55943., 56622., 56076., 54291., 55565., 54753., 55316., 56510., 56593., 58280., 58112., 59354., 59508., 59830., 59774., 59720., 60352., 61062., 60197., 60383., 62077., 61434., 61409., 62046., 61559., 60863., 60725., 59335., 59226., 58273., 58599., 57555., 58101., 58111., 57519., 57933., 58333., 58251., 57303., 57649., 58043., 57903., 57674., 57288., 57416., 57397., 57375., 56610., 56347., 55843., 56699., 57092., 56753., 57061., 56104., 56540., 55981., 55924., 56330., 56196., 55254., 55222., 56086., 55215., 54807., 54411., 54094., 53300., 53501., 53106., 52811., 52449., 52493., 52201., 51868., 51787., 52627., 52250., 51317., 51700., 50665., 51154., 50692., 50475., 50328., 49430., 50273., 49697., 49659., 49063., 48921., 49142., 48431., 48374., 48587., 47717., 48105., 47322., 47748., 47022., 46769., 47170., 46016., 46376., 46418., 46073., 46342., 45866., 45809., 45854., 45620., 46447., 45839., 45895., 45773., 45814., 46322., 46128., 46276., 45886., 45706., 45632., 45236., 45691., 46823., 45675., 45593., 45620., 45590., 45132., 45595., 44343., 44604., 44320., 44520., 44621., 43953., 43822., 42680., 43177., 42882., 42749., 42090., 41447., 41790., 41044., 41576., 41097., 40515., 40769., 40758., 40750., 40255., 40484., 40040., 40121., 39719., 39147., 39235., 39314., 39923., 39422., 39004., 38943., 38509., 12805., 12619., 12776., 12568., 12480., 12535., 12879., 12543., 12161., 12394., 12572., 12042., -180.37, 1.2520, 120.00, 1.000, 2.600, 0.000, 0.000, 0.000, 0.000, 0.000, ------------------------ End of NIO Test Data ------------------------ Title :NELS- Version 20 Keywords :EELS Computer :DEC VAX-11/730 - 11/785 or higher series Operating System :VMS Programming Language :Fortran IV Hardware Requirements :Tektronics 4010 Series Graphic Terminal or Emulator Author(s) :Nestor J. Zaluzec Correspondence Address :Electron Microscopy Center, Materials Science, Bldg 212 :Argonne Nat. Lab, Argonne, IL. 60439 USA Abstract: NELS is a generic data analysis computer program for reduction of EELS spectra. The program performs the functions of spectral calibration, KLM edge identification, background fitting & subtraction, edge integration, digital filtering, and general spectral maniputlation. The input spectra are assumed to be in the EMMPDL spectral file format (see RWEMMPDL.ABS). The program can handle simultaneously three spectra of upto 4096 channels in length. A shorter and less elaborate version is available for smaller DEC computers running RT-11. The program NELS requires a computer terminal which responds to the graphics protocol of Tektronics 4010/4014 series terminal (see NGRAPH.ABS) and an atomic inner shell energy level data file (see ATMLVL.ABS). PROGRAM NELS C C C ELS DATA ANALYSIS PROGRAM C FOR TEKTRONICS 4010 TERMINALS C C NESTOR J. ZALUZEC C VERSION 19 C C ORIGINAL 1977 UNIV. OF ILL. C C MODIFIED 01/09/79 AT ORNL C C MODIFIED 10/02/83 AT ANL C NOW COMPATABILE WITH LSI ADM-5A TERMINAL C EQUIPPED WITH RETRO-GRAPHICS RG-512 C VERSION 14 MAXIMUM BAUD RATE 9600 C NOTE: TO CONVERT CURRENT PROGRAM TO 4010-1 C EQUIVALENT REPLACE ALL CALLS TO THE C SUBROUTINE CRSSPT(X,Y) BY CRSSHR(X,Y) C C MODIFIED 06/29/84 AT ANL C VERSION 15 NOW HARDCOPY OPTION AVAILABLE FOR C HP7470 PLOTTER C C MODIFIED 07/21/84 AT ANL C VERSION 16 ADDITION OF DIGITAL FILTER C C MODIFIED 10/04/84 AT ANL C VERSION 17 ADDITION OF 4027 AND VT-200 MODES C C MODIFIED 10/06/84 AT ANL C VERSION 18 ADDITION OF DEFAULT DATA FILE OUTPUT C TO WD:NELOUT.DAT C C C C MODIFIED 4/26/85 AT ANL C NEW VERSION OF CRT SUBROUTINE (CLEAN UP) C VERSION 19 MAINLY FOR VAX OPERATION C ADDITION OF COMMON /LABEL/ FOR NGRAPH LIBRARY C C C VERSION 20 MULTIPLE WINDOW OPTION ADDED FOR POLYNOMINAL C BACKGROUND FITTING ROUTINES C 5/15/85 C C VERSION 21 VAX VERSION ARRAY LIMITS INCREASED TO 4110 C AND UPDATED NGRAPH LIBRARY 18/JUN/87 NJZ C C ------------ C DEFINITIONS: C ------------ C C AX,BX = ARRAYS HOLDING VALUES FOR X-AXIS C AY,BY = ARRAYS HOLDING VALUES FOR Y-AXIS (# OF COUNTS) C AX,AY = ARRAYS FOR CURRENT SUBPLOT DISPLAY C BX,BY = ARRAYS FOR CURRENT DATA=MEMORY #1 C CCY = ARRAY FOR BACKGROUND STORAGE=MEMORY#2 C DDY = ARRAY FOR ?=MEMORY#3 C C OFFSET= CALIBRATION CONSTANT HORIZONTAL AXIS C EVCH = CALIBRATION CONSTANT HORIZONTAL AXIS C C NNPT = NUMBER OF POINTS IN ORIGINAL SPECTRUM (USUALLY 1020) C NPT = NUMBER OF POINTS IN CURRENT SUBPLOT DISPLAY C C IBFLG = BACKGROUND FLAG C = 0 WHEN BACKGROUND HAS NOT BEEN SUBTRACTED C = 1 WHEN BACKGROUND HAS BEEN SUBTRACTED C IBGND = TYPE OF BACKGROUND FIT C = 0 - NO FIT YET C = 1 - A*E**R MODEL C = 2 - POLYNOMINAL MODEL C = 3 - LOGPOLYNOMINAL MODEL C C IZERO = STARTING CHANNEL OF CURRENT SUBPLOT I.E. AX(1)=BX(IZERO) C C NEL= # OF ELEMENTS INTEGRATED C SYM= ARRAY CONTAINING ATOMIC SYMBOLS OF INTEGRATED ELEMENTS C AEDID= ARRAY CONTAINING EDGE ID`S " " " C PK= ARRAY CONTAINING EDGE INTEGRALS C BGD= ARRAY CONTAINING BGND INTEGRALS C ES = ARRAY CONTAINING WINDOW WIDTHS IN EV FOR ABOVE INTEGRALS C C C INAME= ARRAY CONTAINING NAME OF CURRENT DATA FILE C EO= ACCELERATING VOLTAGE OF ABOVE DATA (KV) C CONV= INCIDENT BEAM CONVERGENCE (MR) C DIV = SCATTERING ANGLE (MR) C C A,R= COEFFICIENTS OF BGND FIT A*E**R MODEL C COE= COEFFICIENTS OF POLYNOMINAL &/OR LOGPOLYNOMINAL FITS C C FLAG= SUBPLOT CONTROL FLAG C > 0 REPLOT CURRENT DISPLAY WINDOW C = 0 RESET ENTIRE WINDOW TO 1-1020 LIMITS C =-1 REPLOT NEW SUBPLOT C =-2 SUPERIMPOSE BGND FIT C =-3 DECONVOLUTION COMPLETED C C AREAZL= ZERO LOSS INTEGRAL FOR DECONVOLUTION SCALING C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL,AREAZL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /NAME/ INAME(8),EO,CONV,DIV C C C COMMON BLOCKS for NGRAPH C C C COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ COMMON /TERMIN/ ITERM,IBAUD,IPLOT,ICOLOR,IMODE,IHARD COMMON /LABEL/ ILX(12),ILY(12),NILX,NILY C C DEFAULT LABELS FOR X & Y AXIS ARE DEFINED HERE C C C ILX= LABEL OF XAXIS STORED AS TWO CHARACTERS PER BYTE C ILY= LABEL OF YAXIS STORED AS ONE CHARACTER PER BYTE C NILX= TOTAL NUMBER OF CHARACTERS IN XAXIS LABEL INCLUDING SPACES C NILY=TOTAL NUMBER OF CHARACTERS IN YAXIS LABEL INCLUDING SPACES DATA ILX /'En','er','gy',' (','eV',') ',' ',' ',' ',' ',' ',' '/ DATA ILY /' ','C ','o ','u ','n ','t ','s ',' ',' ',' ',' ',' '/ DATA XL,YL,XZ,YZ /800.,500.,150.,250./ NILX=12 NILY=7 ICOLOR=1 C C C C C C DEFINE MENU CHARACTERS C C DATA HC,EINT,G,C,W,RE,S,BG,YES,SN,E/'HC','IN','GR','CA','WI','RE', C 'SM','BG','Y','NO','EX'/ DATA WR,CU,OR,AMEM,DE /'WR','CU','OR','ME','DE'/ C C C C C INITIALIZE PARAMETERS, WRITE TITLES, & GET OUTPUT FILE NAMES C CALL TITLE CALL TERMIN ITORG=ITERM !STARTUP VALUE OF ITERM C CALL MENU !LIST MENU FOR NEWUSERS C C C C C READ IN DATA FROM EMMPDL FORMAT FILE C CALL RPDL(BY,NNPT) NPT=NNPT CALL DELAY(3) C C C READ IN OPERATING CONDITIONS IF UNKNOWN C CALL CONDIT(BY) C C C C USE CALIBRATION CONSTANTS TO SCALE HORIZONTAL AXIS C AND FILL UP ARRAYS WITH ORIGINAL DATA FOR BACKUP C IF(BY(NNPT+2).NE.0) EVCH=BY(NNPT+2) IF(BY(NNPT+1).NE.0) OFFSET=BY(NNPT+1) IF(BY(NNPT+2).EQ.0.) BY(NNPT+2)=1. !SAFETY FACTOR NEVER LET EVCH=0. C C JUMP TO 999 TO DRAW DATA AND THEN ENTER GENERAL COMMAND LOOP C GO TO 999 C C C ENTER GENERAL COMMAND INPUT SEQUENCE C 11 CONTINUE CALL PLOT(0.,180.,0) IRPLT=0 YM=1. C C C REQUEST GRAPHICS &/OR DATA MANIPULATION FROM OPERATOR C C CALL DELAY(5) WRITE(7,12) 12 FORMAT($,' GRPH,HC,WIN,CUR,DF,SM,CAL,NOR,BGD,INT,READ, C WRIT,ORIG,MEM,DECON,EX: ') READ(5,13) RES 13 FORMAT(1A2) C C C READ RES(PONSE) AND CHECK WITH MENU OF AVAILABLE OPTIONS C TYPE LIST IF NULL COMMAND C C C IF(RES.EQ.RE) CALL NELOUT !WRITE OUT RESULTS IF NEL>0 IF(RES.EQ.RE) CALL RPDL(BY,NNPT) IF(RES.EQ.RE) NPT=NNPT IF(RES.EQ.RE) CALL CONDIT(BY) C C IF(RES.EQ.WR) CALL WRTPDL(BY,IL,IU,NPT,NNPT) C C IF((RES.EQ.WR).OR.(RES.EQ.RE)) RES=OR !REDRAW SPECTRUM C C IF(RES.EQ.BG) GO TO 50 IF(RES.EQ.EINT) GO TO 51 IF(RES.EQ.SN) GO TO 52 IF(RES.EQ.CU) GO TO 53 IF(RES.EQ.AMEM) GO TO 54 IF(RES.EQ.DE) GO TO 55 IF(RES.EQ.C) GO TO 15 IF(RES.EQ.W) GO TO 30 C C C SPECIAL SECTION OF MENU FOR WINDOWS C COPIES PARTIAL ARRAY FROM BX,BY INTO AX,AY C C 18 IF (FLAG.EQ.0) IL=1 IF(FLAG.EQ.0) IU=NNPT IF(FLAG.GE.0.5) GO TO 17 IF (FLAG.GE.-1.)FLAG=1. NPT=0 DO 17 I=IL,IU J=I-IL+1 AX(J)=BX(I) AY(J)=BY(I) NPT=NPT+1 17 CONTINUE C C STORE CURRENT IL VALUE IN IZERO C IZERO=IL C C C IF((RES.EQ.BG).AND.(FLAG.EQ.-2)) GO TO 510 IF (RES.EQ.G) GO TO 14 IF(RES.EQ.HC) ITORG=ITERM IF(RES.EQ.HC) ITERM=IHARD IF(RES.EQ.HC) GO TO 14 C IF(RES.EQ.OR) FLAG=0. IF(RES.EQ.OR) GO TO 999 IF (RES.EQ.S) GO TO 40 IF(RES.EQ.'DF') GO TO 999 IF(RES.EQ.AMEM) GO TO 11 IF(RES.EQ.C) GO TO 33 IF (RES.EQ.BG) GO TO 33 IF (RES.EQ.W) GO TO 33 IF(RES.EQ.SN) GO TO 33 IF(RES.EQ.EINT) GO TO 33 IF(RES.EQ.E) GO TO 9999 YM=YLAST C C C TYPE MENU IF NULL OR INVALID INPUT FOR COMMAND C C CALL MENU C C RETURN FOR NEW COMMAND C IRPLT=0 GO TO 11 C C CALIBRATE ROUTINE C C 15 CALL CALIB FLAG=-1. NMODE=1 IRPLT=0 GO TO 18 C C A*E**R C C BACKGROUND FITTING ROUTINE C C 50 IF(IBFLG.EQ.0) GO TO 520 !NO BACKGROUND SUBTRACTION C C CHECK FOR DECONVOLUTION IF SO PREVIOUS BGND C WAS ERASED & FFT IS STORED IN MEM#2 C IN THIS CASE DONOT ATTEMPT TO ADD DATA BACK C IF(FLAG.EQ.-3.) GO TO 520 !FLAG=-3 MEANS DECONVOLUTION C C CHECK FOR PREVIOUS BGND SUBTRACTION C IF SO ADD IT BACK!! & REDRAW C DO 500 I=1,NNPT BY(I)=BY(I) + CCY(I) 500 CCY(I)=0. C C RESET IBFLG=0 SINCE ADD BGND BACK C ALSO SET FLAG=-2 FOR REDRAWING PURPOSES C IN BGND OVERLAY SITUATIONS C IBFLG=0 FLAG=-2. GO TO 18 510 AY(NPT)=0. CALL CRT( AX,AY,NPT,1,0,1.,1.,AX,AY,NPT) 520 CONTINUE CALL BGDMOD NMODE=1 FLAG=-1. IRPLT=0 IF(IBFLG.EQ.0) GO TO 11 !NO SUBTRACTION DONOT REDRAW GO TO 18 C C INTEGRATION ROUTINE C C C 51 CALL DINT IF(FLAG.NE.-3) FLAG=-1. NMODE=1 IRPLT=0 GO TO 18 C C NORMALIZE ROUTINE C C C 52 CALL NORM FLAG=-1. NMODE=1. IRPLT=0. GO TO 18 C C C WINDOW ROUTINE C C 30 CALL WINDOW(IL,IU) NMODE=1 IRPLT=0 IF (FLAG.NE.-3) FLAG=-1. GO TO 18 C C C CURSOR ROUTINE C 53 CALL CURSOR GO TO 11 C C C MULTIPLE SCATTERING DECONVOLUTION C USING LEAPMAN/SWYT METHOD C 55 CALL DECONL IF(FLAG.GE.-2) GO TO 11 C C REZERO DATA ARRAY C DO 56 I=1,NNPT CCY(I)=0. 56 CONTINUE C C COPY DECONVOLUTED DATA INTO MEM#2 C DO 57 I=1,512 CCY(I+IROI-1)=AY(I) 57 CONTINUE C C COPY BACK CALIB. CONSTANTS INTO MEM#2 C DO 58 I=NNPT+1,NNPT+10 58 CCY(I)=BY(I) GO TO 11 C C C C MEMORY & DATA MANIPULATION ROUTINES C C C 54 CALL MEMORY GO TO 999 C C SMOOTH ROUTINE C 40 CALL SMOOTH NMODE=1 IRPLT=1 YM=YLAST GO TO 33 C C DIGITAL FILTER ROUTINE C PROGRAM SHOULD HAVE RESTORED ORIGINAL FULL SCREEN WIDTH C 400 WRITE(7,401) 401 FORMAT($,' Select Type of Digital Filter ? [x,y]=1, [x,y,z]=2 [1]:') READ(5,402) NMODE !USE NMODE TO SAVE CREATING XTRA VAR. 402 FORMAT(I2) NTM=1 !DEFAULT TO [X,Y,Z] IF (NMODE.EQ.2) NTM=2 CALL LDIGFIL(NTM) NMODE=1 IRPLT=O YM=1. GO TO 33 C C C GRAPHICS REDRAW ROUTINE C C 14 CONTINUE IF(ITERM.EQ.IHARD) CALL HCSCRN !COPY CURRENT SCREEN IF(ITERM.GE.11.AND.ITERM.LE.12) GO TO 888 CALL DELAY(1) WRITE(7,6) 6 FORMAT($,' OVERLAY(Y-N),Y-MAG,SYM(-9 TO +1),ICOLOR(-1,0 TO 7)') READ(5,66) RES,YM,ANMODE,ACOLOR 66 FORMAT(1A1,1X,3(F10.0)) NMODE=INT(ANMODE) ICOLOR=INT(ACOLOR) IRPLT=0 IF (RES.EQ.YES) IRPLT=1 C C MAKE SURE LAST POINT OF CURRENT PLOT HAS ZERO VALUE FOR PLOTTING C TO AVOID PROBLEMS WITH SCALING C BUT DONOT PLOT THIS POINT ON GRAPH C 33 AY(NPT)=0. CALL CRT( AX,AY,NPT-1,NMODE,IRPLT,XM,YM,AX,AY,NPT) YLAST=YM 888 IF(ITERM.NE.ITORG) ITERM=ITORG !GET OUT OF HARDCOPY MODE GO TO 11 C C C RECALCULATE BX ARRAY & REFILL AX,AY ARRAYS C C C 999 FLAG=0. OFFSET=BY(NNPT+1) IF(BY(NNPT+2).EQ.0.) BY(NNPT+2)=1. !SAFETY FACTOR NEVER LET EVCH=0 EVCH=BY(NNPT+2) DO 16 I=1,NNPT BX(I)=OFFSET + EVCH*I AY(I)=BY(I) 16 AX(I)=BX(I) C C C STORE SPECTRAL CONSTANTS C DO 4 I=NNPT+1,NNPT+10 4 AY(I)=BY(I) C C IF MEMORY CHANGE DON'T REPLOT THE WHOLE THING C IF(RES.EQ.AMEM) GO TO 18 NPT=NNPT NMODE=1 IRPLT=0 C C IF DIGITAL FILTER GO TO 400 & FILTER ENTIRE SPECTRUM C IF(RES.EQ.'DF') GO TO 400 GO TO 33 9999 CONTINUE C C MUST RESTART TO FOOL SYSTEM AND ERASE GRAPHICS AT END OF PROGRAM C CALL STRTCRT CALL ERASE CALL VT100 C C CREATE OUTPUT DATA FILE FOR USE WITH NELQNT C CALL NELOUT STOP END C-----------------------------------------------------------------------C SUBROUTINE TITLE C-----------------------------------------------------------------------C C C INITIALIZE PARAMETERS & PRINT TITLES ETC C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ COMMON /TERMIN/ ITERM,IBAUD,IPLOT,ICOLOR,IMODE,IHARD C C C C INITIALIZE SOME PARAMETERS C C IZERO=0 NEL=0 IBFLG=0 IBGND=0 EVCH=1. OFFSET=0. CFLAG=0. A=1. R=0. IL=1 IU=4110 NPT=4110 NNPT=NPT NMODE=0 IRPLT=0 XM=1. YM=1. IMODE=-1 C C INITIALIZE ALL DATA ARRAY'S C DO 992 I=1,NNPT BY(I)=0. CCY(I)=0. DDY(I)=0. 992 CONTINUE C C C WRITE TITLES C C C CALL ERASE CALL DELAY(10) 9991 WRITE(7,1) 1 FORMAT(//////, C 18X,'******************************',/, C 18X,'* NELS *',/, C 18X,'* ELS DATA ANALYSIS PROGRAM *',/, C 18X,'* 8706220021-NJZ *',/, C 18X,'******************************',//) CALL DELAY(10) C C C C GIVE OPERATOR OPTION FOR OUTPUT NL:,TT:,LP:,DEV:NAME.EXT C ON VAX YOU DONT HAVE THIS OPTION EVERYTHING GOES TO A FOR008.DAT FILE C WRITE(7,2) 2 FORMAT(/,' Hardcopy Print File = FOR008.DAT') C CALL ASSIGN(8,'DUMMY',-6) WRITE(8,1) RETURN END C-----------------------------------------------------------------------C SUBROUTINE DECONL C-----------------------------------------------------------------------C C C MULTIPLE SCATTERING DECONVOLUTION ROUTINE C LOW-LOSS FROM HIGH-LOSS C LEAPMAN METHOD C REPLACE ZERO LOSS PEAK BY DELTA FUNCTION OF EQUAL AREA C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),AHL(515),ALL(515) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL,AREAZL COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ DIMENSION SX(2),SY(2) COMPLEX ALL,AHL C C C C ----------- C DEFINITIONS C ----------- C C ILLS= LOW LOSS START CHANNEL C ILLF= LOW LOSS FINISH CHANNEL C IZLS= ZERO LOSS START CHANNEL C IZLF= ZERO LOSS FINISH CHANNEL C IZLC= ZERO LOSS CENTROID ( MAXIMUM VALUE) C IHLS= HIGH LOSS START CHANNEL C IHLF= HIGH LOSS FINISH CHANNEL C ALL= COMPLEX ARRAY CONTAINING LOW LOSS FFT C AHL= COMPLEX ARRAY CONTAINING HIGH LOSS FFT C THESE ARE STORED IN CCY & DDY MEMORY AREAS TO SAVE COMMON SPACE C C LET OPERATOR CHOOSE OPTION C THIS WAY HE CAN SWITCH BACK & FORTH TO SET NEW C WINDOWS TO OPTIMIZE DISPLAY FOR LOW LOSS & HIGH LOSS REGIONS C C C C CALL ADM5E(9) !ERASE 9 LINES C C C CALL PLOT(0.,136.,0) CALL DELAY(5) 1 WRITE(7,2) CALL DELAY(5) 2 FORMAT($,' SELECT OPTION:LOW LOSS,HIGH LOSS,DECONVOL,BACK ') READ(5,3) RES 3 FORMAT(1A1) IF(RES.EQ.'L') GO TO 10 IF(RES.EQ.'H') GO TO 75 IF(RES.EQ.'D') GO TO 105 IF(RES.EQ.'B') RETURN GO TO 1 C C C OPERATOR HAS SELECTED LOW LOSS FFT C C C 10 CALL DELAY(5) WRITE(7,20) CALL DELAY(5) 20 FORMAT($,' SET WINDOW AROUND LOW LOSS') CALL DELAY(5) CALL WINDOW(ILLS,ILLF) CALL DELAY(5) C C COMPUTE AREA UNDER LOW LOSS REGION C WIDTH=(ILLF-ILLS)*EVCH AREAZL=0. DO 21 I=ILLS,ILLF 21 AREAZL=AREAZL +BY(I) WRITE (8,22) AREAZL,WIDTH 22 FORMAT(' LOW LOSS INTEGRAL = ',1PE10.4,' WIDTH (EV) = ',1F6.1) CALL PLOT(0.,70.,0) CALL DELAY(5) CALL DELAY(5) WRITE(7,30) CALL DELAY(5) 30 FORMAT($,' SET WINDOW AROUND ZERO LOSS') CALL WINDOW(IZLS,IZLF) CALL DELAY(5) CALL PLOT(0.,48.,0) CALL DELAY(5) C C COMPUTE AREA UNDER ZERO LOSS C ZMAX=0. WIDTH=(IZLF-IZLS)*EVCH AREAZL=0. DO 40 I=IZLS,IZLF C C WHILE COMPUTING AREA ALSO FIND CHANNEL OF MAXIMUM C IF(BY(I).GT.ZMAX) IZLC=I IF(BY(I).GT.ZMAX) ZMAX=BY(I) 40 AREAZL=AREAZL + BY(I) WRITE(8,41) AREAZL,WIDTH C WRITE(7,41) AREAZL 41 FORMAT(' ZERO LOSS INTEGRAL = ',1PE10.4,' WIDTH (EV) = ',1F6.1) C C FIND MINIMUM BETWEEN ZERO LOSS & 1ST PLASMON C IMIN=IZLC DMIN=BY(IMIN) DO 42 I=1,ILLF-IZLC J=IMIN+I-1 IF(BY(J).LT.DMIN) DMIN=BY(J) IF(BY(J).GT.DMIN) GO TO 43 42 CONTINUE 43 IZLF=J C C EXTRAPOLATE 1ST PLASMON LINEARLY TO ZERO ON C LOW LOSS SIDE TO ZERO LOSS PEAK C STORE DATA IN ALL ARRAY (ALOWLOSS) C SLOPE=-BY(IZLF)/FLOAT(1-(IZLF-IZLC)) SINTER= BY(IZLF)-SLOPE*(IZLF-IZLC) DO 50 I=1,IZLF-IZLC !IZLF-IZLC=# CHANNELS IZLC TO IZLF TEMP=SLOPE*I +SINTER 50 ALL(I)=CMPLX(TEMP,0.) !STORE AS REAL PART OF COMPLEX # C C REPLACE ZERO LOSS BY DELTA FUNCTION = AREA C ALL(1)=CMPLX(AREAZL,0.) C C FILL IN REST OF DATA ARRAY C DO 60 I=IZLF-IZLC+1,ILLF-IZLF TEMP=BY(IZLC+I) IF(TEMP.LT.0.) TEMP=0. C C MAKE SURE NO NEGATIVE COUNTS EXIST!! C 60 ALL(I)=CMPLX(TEMP,0.) C C EXTRAPOLATE REMAINDER OF ARRAY TO ZERO C SLOPE=BY(ILLF)/FLOAT(ILLF-IZLC-512) SINTER=BY(ILLF)-SLOPE*(ILLF-IZLC) DO 70 I=ILLF-IZLF+1, 512 TEMP=SLOPE*I+SINTER 70 ALL(I)=CMPLX(TEMP,0.) C C DRAW MODIFIED DATA C D CALL DECONP(RES,IZLC,ILLS,ILLF,IHLS,IHLF) CALL PLOT(0.,48.,0) CALL DELAY(5) WRITE(7,71) CALL DELAY(5) 71 FORMAT(' CALCULATING FFT') CALL CFFT(1,ALL,9) !USE M.DISKO'S ROUTINE IROI=3 CALL DELAY(5) WRITE(7,72) IROI CALL DELAY(5) 72 FORMAT(' FFT IN MEM #',I1) C C RETURN TO DISPLAY SO OPERATOR CAN SET HIGH LOSS REGION C RETURN C C OPERATOR HAS CHOOZEN HIGH LOSS FFT C DELINEATE HIGH LOSS EDGE WITH WINDOW C 75 CALL DELAY(5) WRITE(7,80) CALL DELAY(5) 80 FORMAT($,' SET WINDOW AROUND HIGH LOSS') CALL WINDOW(IHLS,IHLF) CALL PLOT(0.,70.,0) CALL DELAY(5) C C STORE DATA IN AHL ARRAY (AHIGHLOSS) C DO 90 I=1,IHLF-IHLS TEMP=BY(IHLS+I-1) IF(TEMP.LT.0.) TEMP=0. 90 AHL(I)=CMPLX(TEMP,0.) C C EXTRAPOLATE TO ZERO AT CHANNEL 512 C SLOPE= BY(IHLF)/FLOAT(IHLF-IHLS-512) SINTER=BY(IHLF)-SLOPE*(IHLF-IHLS) DO 100 I=IHLF-IHLS+1,512 TEMP=SLOPE*I +SINTER 100 AHL(I)=CMPLX(TEMP,0.) C C DRAW MODIFIED DATA C D CALL DECONP(RES,IZLC,ILLS,ILLF,IHLS,IHLF) CALL PLOT(0.,48.,0) C C CALCULATE FOURIER COEFFICIENTS C CALL DELAY(5) WRITE(7,71) CALL DELAY(5) CALL CFFT(1,AHL,9) IROI=2 CALL DELAY(5) WRITE(7,72) IROI CALL DELAY(5) C C RETURN TO DISPLAY MODE C RETURN C C DECONVOLUTE LOW LOSS FROM HIGH LOSS C BY DIVISION OF FOURIER COEFFICIENTS C 105 DO 130 I=1,512 130 AHL(I)=AHL(I)/ALL(I) C C RESULT IN AHL ARRAY C C C CALCULATE INVERSE TRANSFORM C CALL DELAY(5) WRITE(7,71) CALL DELAY(5) CALL CFFT(-1,AHL,9) C C TELL OPERATOR WHERE DATA IS STORED C IROI=2 CALL DELAY(5) WRITE(7,72) IROI CALL DELAY(5) C C DRAW DECONVOLUTED DATA ON TOP OF CURRENT DISPLAY C THIS ASSUMES HIGH LOSS WAS LAST DISPLAY!! C C DRAW MODIFIED DATA ON TOP OF OLD DATA C C 131 CALL DECONP(RES,IZLC,ILLS,ILLF,IHLS,IHLF) C C SET FLAG=-3 TO INDICATE DECONVOLUTION COMPLETED C C STORE IHLS IN IROI C IROI=IHLS FLAG=-3. RETURN END C-----------------------------------------------------------------------C SUBROUTINE RRTNJZ(XS) C-----------------------------------------------------------------------C C----------------------------------------------------------------------- C SUBROUTINE TO READ AN RT-11 DATA FILE IN NJZ FORMAT C INTO AN ARRAY XS(4110) C FIRST TWO LINES OF TEXT ARE FOR IDENTIFICATION 72A1 FORMAT C XS(1-4110)=DATA 10 VALUES/LINE @F7.0 WITH NNPT+ LINES C XS(NNPT+1-4110)=CALIBRATION CONSTANTS OF SPECTRAL DATA C NNPT+1= OFFSET ENERGY OF FIRST CHANNEL C NNPT+2= EV/CHANNEL C NNPT+3= ACCELERATING VOLTAGE (KEV) C NNPT+4= INCIDENT BEAM DIVERGENCE (MR) C NNPT+5= SCATTERING ANGLE (MR) C NNPT+6-4110= RESERVED FOR FUTURE USE C C C C MODIFIED 24-OCT-84 TO CHECK DATA FORMAT C C C SOME OLDER DATA FILES RUN INTO AN ERROR WHEN USING C 1F7.0 FORMAT IF THE DATA VALUES EXCEED 999,999 AND C RESULT IN AN ERROR IN THE DECIMAL POINT POSITION C THE CORRECTION IN NTRAN PROGRAM WAS TO ADD A SPACE IN FRONT C OF EACH LINE BEFORE WRITING IT TO A DISK FILE TO CHECK FOR C THIS IN RRTNJZ WE LOOK FOR THE POSITION OF THE FIRST DECIMAL C POINT ON THE THIRD LINE IF IT IS IN THE OLD FORMAT POSITION C WE USE A DIFFERENT FORMAT STATEMENT FOR THE DATA READ C C MODIFIED 8-NOV-84 TO ADD FILE NAME & DAY/TIME C----------------------------------------------------------------------------- C C C COMMON BLOCKS FOR NELS C C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) C COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT C COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) C COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL C COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ COMMON /TERMIN/ ITERM,IBAUD,IPLOT,ICOLOR,IMODE,IHARD COMMON /NAME/ INAME(8),EO,CONV,DIV C C C C DIMENSION XS(1024),IXL(72) BYTE DAY(9),TIM(8) C C CALL ERASE !ERASE SCREEN C C READ IN RT-11 FILE NAME C 204 WRITE(7,1) 1 FORMAT(' ENTER RT-11/NJZ DATA FILE NAME (DEV:NAME.EXT):',/) READ(5,205) (INAME(J),J=1,7) 205 FORMAT(7A2) INAME(8)=0 OPEN (UNIT=11,NAME=INAME,TYPE='OLD',DISPOSE='SAVE', C ACCESS='SEQUENTIAL') C C C WRITE IT OUT C C CALL DATE(DAY) CALL TIME(TIM) WRITE(8,206) DAY,TIM,(INAME(J),J=1,7) 206 FORMAT(1X,9A1,' at ',8A1,/,' Input Data File = ',7A2) C C READ IN TWO LINES OF TEXT C WRITE(7,21) WRITE(8,21) 21 FORMAT(' SPECTRUM IDENTIFICATION',/,23('-')) C C FILL XS(1021-4110) WITH ZERO'S C XS(NNPT+1)=0. XS(NNPT+2)=0. !DEFAULT TO 0EV/CH IF NO CALIBRATIONS DO 22 I=NNPT+3,NNPT+10 22 XS(I)=0.0 READ(11,2) IXL 2 FORMAT(72A1) WRITE(7,3) IXL WRITE(8,3) IXL 3 FORMAT(1X,72A1) READ(11,2) IXL WRITE(7,3) IXL WRITE(8,3) IXL C C CHECK FOR FORMAT OF FILE. C IDEC = 7 !IDEC SHOULD HAVE COLUMN NUMBER WHERE C !DECIMAL POINT APPEARS IN OLD FORMAT IOLD = 0 INEW = 1 READ(11,2) IXL IFORM = IOLD !OLD FORM OF FILE IF (IXL(IDEC) .NE. '.') IFORM = INEW !NEW FORM OF FILE C C NOW THAT WE KNOW THE TYPE OF DATA FILE REWIND C REREAD TITLES AND SKIP THEN CONTINUE C REWIND 11 READ(11,2) IXL READ(11,2) IXL C C READ IN DATA C C C WATCH OUT FOR OVERFLOWS DATA FORMAT => MAX VALUE =999,999 C XMAX=999999. DO 4 I=1,102 IF=10*I-9 IL=10*I C C MAX DATA ARRAY IN-CASE OF ERROR C DO 45 J=IF,IL 45 XS(J)=XMAX IF (IFORM.EQ.IOLD) READ(11,5,END=4,ERR=4) (XS(J), J=IF,IL) IF (IFORM.EQ.INEW) READ(11,8,END=4,ERR=4) (XS(J), J=IF,IL) 5 FORMAT(10F7.0) 8 FORMAT(1X,10F7.0) 4 CONTINUE C C READ IN CALIBRATION CONSTANTS C IF THEY EXIST C NOTE: SOME VERY OLD NJZ DATA FILES DO NOT HAVE CALIBRATIONS C THIS LOOP ALLOWS BOTH OLD DATA AND NEW DATA TO BE ACCESSED C IF (IFORM.EQ.IOLD) READ(11,51,END =52) (XS(J), J =1021,1025) IF (IFORM.EQ.INEW) READ(11,54,END =52) (XS(J), J =1021,1025) 51 FORMAT(10F7.3) 54 FORMAT(1X,10F7.3) 52 CONTINUE IF((XS(1021).EQ.0).AND.(XS(1022).EQ.0)) WRITE(7,53) IF((XS(1021).EQ.0).AND.(XS(1022).EQ.0)) WRITE(8,53) 53 FORMAT(' UNCALIBRATED SPECTRUM:') WRITE(7,7) XS(1021),XS(1022),XS(1023),XS(1024),XS(1025) 7 FORMAT(//, C ' OFFSET = ',1F9.3,' EV',/, C ' EV/CHANNEL = ',1F9.3,/, C ' ACCELERATING VOLTAGE = ',1F9.1,' KEV',/, C ' INCID. BEAM DIVERGENCE = ',1F9.3,' MR',/, C ' SCATTERING ANGLE = ',1F9.3,' MR',/) 6 WRITE(7,600) 600 FORMAT (' HIT TO CONTINUE') READ(5,5) X CALL CLOSE(11) CALL ERASE !ERASE SCREEN RETURN END C------------------------------------------------------------------ SUBROUTINE WRTPDL(BY,IL,IU,NPT,NNPT) C------------------------------------------------------------------ C C DRIVER ROUTINE TO WRTIE DATA IN EMMPDL DATA FORMAT C BUT ALLOWING ENTIRE DATA ARRAY OR ONLY THE SUBSET C DEFINED BY THE CURRENT WINDOW C C DIMENSION BY(4110),CY(4110) WRITE (7,100) 100 FORMAT ($,' ENTIRE SPECTRUM OR CURRENT WINDOW {Entire,Current} [E]: ') READ (5,110) ANS 110 FORMAT (1A1) ILL=1 IUL=NNPT IF (ANS.NE.'C') GO TO 120 ILL=IL IUL=IU 120 DO 112 I=ILL,IUL 112 CY(I+1-ILL)=BY(I) NPNT=IUL-ILL+1 DO 114 I=1,10 114 CY(NPNT+I) = BY(NNPT+I) EVCH=BY(NNPT+2) OFFSET=BY(NNPT+1) CY(NPNT+1) = (EVCH*(IL-1)+OFFSET) - EVCH ! NEW OFFSET CALL WPDL(CY,NPNT) RETURN END C------------------------------------------------------------------ SUBROUTINE WRTNJZ(XS) C------------------------------------------------------------------ C C C SUBROUTINE TO WRITE AN RT-11 DATA FILE IN NJZ FORMAT C FROM AN ARRAY XS(4110) C FIRST TWO LINES OF TEXT ARE FOR IDENTIFICATION 72A1 FORMAT C XS(1-4110)=DATA 10 VALUES/LINE @F7.0 WITH 102 LINES C XS(1021-4110)=CALIBRATION CONSTANTS OF SPECTRAL DATA C 1021= OFFSET ENERGY OF FIRST CHANNEL C 1022= EV/CHANNEL C 1023= ACCELERATING VOLTAGE (KEV) C 1024= INCIDENT BEAM DIVERGENCE (MR) C 1025= SCATTERING ANGLE (MR) C 1026-4110= RESERVED FOR FUTURE USE C C---------------------------------------------------------------------- C C DIMENSION XS(1024),IXL(72) INTEGER IB(8) C C ERASE SCREEN C CALL ERASE !ERASE SCREEN C C READ IN RT-11 FILE NAME C 204 WRITE(7,1) 1 FORMAT(' ENTER RT-11/NJZ DATA FILE NAME (DEV:NAME.EXT):',/) C CALL ASSIGN(11,'DUMMY',-6,'NEW','NC',1) READ(5,205) (IB(J),J=1,7) 205 FORMAT(7A2) IB(8)=0 !FILE NAME=IB MUST BE NULL TERMINATED OPEN(UNIT=11,NAME=IB,TYPE='NEW',DISP='SAVE', C FORM='FORMATTED',ACCESS='SEQUENTIAL',ERR=200) GO TO 201 C C ERROR CONDITION FILE DOESNT EXIST ? TYPO? OR QUIT C 200 WRITE(7,202) (IB(J),J=1,7) 202 FORMAT(' THE FILE --->',7A2,'<--- CANNOT BE OPENED ',/, C $,' DO YOU WISH TO REENTER THE FILE NAME?') READ(5,203) ANS 203 FORMAT(1A1) IF(ANS.EQ.'N') RETURN GO TO 204 201 CONTINUE C C READ IN TWO LINES OF TEXT C WRITE(7,21) 21 FORMAT(' ENTER SPECTRUM IDENTIFICATION (2 LINES OF TEXT)') READ(5,2) IXL 2 FORMAT(72A1) WRITE(11,3) IXL 3 FORMAT(1X,72A1) READ(5,2) IXL WRITE(11,3) IXL C C WRITE IN DATA C C C CHECK DATA TO INSURE THAT MAX VALUE<999,999. (1F7.0) C IF GREATER THEN SCALE EVERYTHING DOWN!! C DMAX=0. DO 31 I=1,NNPT 31 DMAX=AMAX1(DMAX,XS(I)) IF(DMAX.LT.1.E06) GO TO 33 DO 32 I=1,NNPT 32 XS(I)=XS(I)*999999./DMAX 33 CONTINUE DO 4 I=1,102 IF=10*I-9 IL=10*I WRITE(11,5) (XS(J), J=IF,IL) C C C This format statement changed by RJH on 6-18-84 to eliminate problem C with losing first character on output line. C5 FORMAT(10F7.0) C C 5 FORMAT(1X,10F7.0) 4 CONTINUE C C WRITE IN CALIBRATION CONSTANTS C IF THEY EXIST C WRITE(11,51) (XS(J), J =1021,4110) 51 FORMAT(1X,F7.2,F7.4,F7.2,7F7.3) WRITE(7,7) XS(1021),XS(1022),XS(1023),XS(1024),XS(1025) 7 FORMAT(//, C ' OFFSET = ',1F9.3,' EV',/, C ' EV/CHANNEL = ',1F9.3,/, C ' ACCELERATING VOLTAGE = ',1F9.1,' KEV',/, C ' INCID. BEAM DIVERGENCE = ',1F9.3,' MR',/, C ' SCATTERING ANGLE = ',1F9.3,' MR',/) CALL DELAY(2) 6 WRITE(7,600) 600 FORMAT (' HIT TO CONTINUE') READ(5,203) ANS CALL CLOSE (11) CALL ERASE !ERASE SCREEN RETURN END C-----------------------------------------------------------------------C SUBROUTINE MEMORY C-----------------------------------------------------------------------C C C SUBROUTINE TO CHANGE CURRENT SPECTRUM MEMORY CHARACTERISTICS C C COPY BETWEEN MEMORIES 1->2 OR 1->3 C SPECTRUM-SPECTRUM MATH (SSM) C SPECTRUM-CONSTANT MATH (SCM) C LOG/ANTILOG DISPLAY (LOG/ALOG) C MEMORY COMPARE 1/2, 1/3 C FFT AND INVERSE FFT OF MEM#1 (FFT=FAST FOURIER TRANSFORM) C DELETE PORTIONS OF CURRENTLY DISPLAYED SPECTRUM C EQUALIZE (SCALE) SPECTRA TO SAME VALUE C ZERO ALL NEGATIVE VALUES STORED IN CURRENT SPECTRUM C C NOTE ARRAY DEFINITIONS: BY=1, CCY=2, DDY=3 C C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ COMMON /TERMIN/ ITERM,IBAUD,IPLOT,ICOLOR,IMODE,IHARD C C C C COMPLEX AFFT(515) EQUIVALENCE (AFFT,DDY) C C CALL ADM5E(6) !ERASE 6 LINES C C C 1 WRITE(7,10) 10 FORMAT($,' SELECT OPTION:COPY,SSM,SCM,LOG,ALOG,OVRLY,EQUAL,DEL C ,ZERO:') READ(5,20) RES 20 FORMAT(1A2) IF(RES.EQ.'CO') CALL COPY IF(RES.EQ.'SS') CALL SSM IF(RES.EQ.'SC') CALL SCM IF(RES.EQ.'ZE') GO TO 30 IF(RES.EQ.'LO') GO TO 30 IF(RES.EQ.'AL') GO TO 30 IF(RES.EQ.'OV') CALL OVRLAY(RES) IF(RES.EQ.'DE') CALL DELETE IF(RES.EQ.'EQ') CALL EQUAL IF(RES.EQ.'FF') IDIR = 1 IF(RES.EQ.'IF') IDIR =-1 IF((RES.EQ.'FF').OR.(RES.EQ.'IF')) GO TO 60 GO TO 50 C C CALCULATE LOG OR ANTI-LOG AS SELECTED C OR SET ALL NEGATIVE VALUES TO ZERO C 30 DO 40 I=1,NNPT IF((RES.EQ.'ZE').AND.(BY(I).LT.0.)) BY(I)=0. IF(RES.EQ.'LO') BY(I)=ALOG10(BY(I)) IF(RES.EQ.'AL') BY(I)=10.**BY(I) 40 CONTINUE 50 RETURN C C FOURIER TRANSFORMS SELECTED C NOTE MEM#2 & #3 USED & DESTROYED C 60 IF(NPT.GT.512) WRITE(7,100) NPT IF(NPT.GT.512) RETURN 100 FORMAT(' # PTS = ',I3,' REDUCE WINDOW TO 512') DO 70 I=1,NPT C C IF OPERATOR REQUESTS IFT,THEN ASSUME HE HAS COPIED C DATA BACK FROM MEM#2 CCY INTO MEM#1 BY & THAT C REAL PART IS IN 1ST 1/2 AND IMAG PART IS IN 2ND 1/2 C C IF OPERATOR REQUESTS FFT, THEN ASSUME IMAG PART IS ZERO C AND ONLY TAKE 512 CHANNELS AS ALL REAL DATA C C IF(RES.EQ.'IF') AFFT(I)=CMPLX(AY(I),AY(I+512)) IF(RES.EQ.'FF') AFFT(I)=CMPLX(AY(I),0.) 70 CONTINUE CALL CFFT(IDIR,AFFT,9) DO 80 I=1,NPT CCY(I)=REAL(AFFT(I)) CCY(I+512)=AIMAG(AFFT(I)) 80 CONTINUE C C CHANGE CALIBRATION CONSTANTS TO CHANNEL NUMBER C CCY(NNPT+1)=1. CCY(NNPT+2)=1. WRITE(7,90) 90 FORMAT(' REAL TRANS = 1ST 1/2 OF MEM#2, IMAG = 2ND 1/2') CALL CRT( AX,CCY,NPT,1,0,1.,1.,AX,CCY,NPT) DO 110 I=1,NPT 110 AY(I)=CCY(I+512) CALL CRT( AX,AY,NPT,1,1,1.,1.,AX,AY,NPT) RETURN END C-----------------------------------------------------------------------C SUBROUTINE OVRLAY(ANS) C-----------------------------------------------------------------------C C C SUBROUTINE TO OVERLAY CONTENTS OF MEM#1 WITH MEM#2 OR #3 C USING THE CURRENTLY SELECTED WINDOW & SCALE FACTORS C C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ COMMON /TERMIN/ ITERM,IBAUD,IPLOT,ICOLOR,IMODE,IHARD C C C C DIMENSION SX(2),SY(2) C C GET CURRENT SCALE FACTORS C SX(1)=XMIN SX(2)=XMIN + XL/DX SY(1)=YMIN SY(2)=YMIN + YL/DY IF((ANS.EQ.'1/2').OR.(ANS.EQ.'1/3')) GO TO 21 1 WRITE (7,10) 10 FORMAT($,' SELECT OVERLAY:1/2,1/3:') READ(5,20) ANS 20 FORMAT(1A3) C C GET CURRENT ROI LIMITS C 21 IL=IZERO IU=IL + NPT DO 30 I=1,NPT J=I + IZERO -1 IF(ANS.EQ.'1/2') AY(I)=CCY(J) IF(ANS.EQ.'1/3') AY(I)=DDY(J) 30 CONTINUE C C DRAW 2ND SPECTRA IN OVERLAY AND CHANGE COLOR C ICOLD=ICOLOR ICOLOR=ICOLOR+1 IF(ICOLOR.GT.7) ICOLOR=1 CALL CRT( AX,AY,NPT,1,1,1.,1.,SX,SY,2) ICOLOR=ICOLD C C NORMAL EXIT C RETURN END C-----------------------------------------------------------------------C SUBROUTINE EQUAL C-----------------------------------------------------------------------C C C THIS ROUTINE EQUALIZES (SCALES) TWO SPECTRA C SO THAT THEIR VALUES ARE EQUAL AT THE OPERATOR SELECTED C CHANNEL C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ COMMON /TERMIN/ ITERM,IBAUD,IPLOT,ICOLOR,IMODE,IHARD C C C C C LET OPERATOR CHOOSE TWO SPECTRA C ANS=' ' !SET ANS(WER) TO BLANK TO QUIRIE OPERAT CALL OVRLAY(ANS)!SAVE ANS FOR REUSE LATER CALL PLOT(0.,92.,0) 10 WRITE(7,20) 20 FORMAT($,' SET CRSSPT TO CHANNEL') CALL DELAY(3) CALL CRSSPT(X,Y) CALL TRANXY(X,Y,E,C) C C COMPUTE CHANNEL NUMBER C IEQUAL=ICHAN(E) C C COMPUTE VALUES AT EACH CHANNEL & FIND LARGEST C TEMP1=BY(IEQUAL) TEMP2=CCY(IEQUAL) IF(ANS.EQ.'1/3') TEMP2=DDY(IEQUAL) D WRITE(7,21) TEMP1,TEMP2 D21 FORMAT(1X,2PE10.3) C C SCALE LARGER TO SMALLER SO NOISE IS NOT MULTIPLIED C DO 30 I=1,NNPT IF(TEMP1.GE.TEMP2) BY(I)=BY(I)*TEMP2/TEMP1 IF((TEMP2.GT.TEMP1).AND.(ANS.EQ.'1/3')) DDY(I)=DDY(I)*TEMP1/TEMP2 IF((TEMP2.GT.TEMP1).AND.(ANS.EQ.'1/2')) CCY(I)=CCY(I)*TEMP1/TEMP2 30 CONTINUE C C REDRAW C CALL OVRLAY(ANS) RETURN END C-----------------------------------------------------------------------C SUBROUTINE NORM C-----------------------------------------------------------------------C C C C SPECTRUM NORMALIZE (GAIN CHANGE REMOVAL) ROUTINE C BY FITTING A*E**R MODEL EITHER SIDE OF GAIN CHANGE POINT C C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C C C CALL ADM5E(6) !ERASE 6 LINES C C C C SET WINDOW CENTERED ABOUT GAIN CHANGE PT C CALL WINDOW (IL,IU) C C SET CURSOR TO CENTROID C CALL CRSSPT(X,Y) CALL TRANXY(X,Y,EC,C) C C GET CHANNEL OF CENTROID C IC=ICHAN(EC) 10 CALL PLOT(0.,120.,0) WRITE(7,6) 6 FORMAT($,' CALCULATE GAIN?') READ(5,61) ANS 61 FORMAT(1A1) 11 IF(ANS.EQ.'Y') GO TO 8 C C OPERATOR WANTS TO ENTER GAIN VALUE C WRITE(7,9) 9 FORMAT($,' ENTER GAIN: ') C C NOTE A GAIN OF < 1. MULTIPLIES SPECTRUM C READ (7,2) GAIN 2 FORMAT(2F10.0) IF (GAIN.EQ.0.) GO TO 10 GO TO 13 8 CALL BGDFIT(IL,IC-2) AL=A*(OFFSET+EVCH*IC)**R CALL BGDFIT(IC+2,IU) AU=A*(OFFSET+EVCH*IC)**R GAIN=AU/AL 13 WRITE(7,14) EC,GAIN WRITE(8,14) EC,GAIN 14 FORMAT(5X,' GAIN CHANGE AT ',/,5X,1F7.2,' EV = ',1F8.3) DO 12 I=1,IC 12 BY(I)=BY(I)*GAIN CALL BGDFIT(IL,IC-4) IL=IC-4 IU=IC+4 DO 15 JJ=IL,IU 15 BY(JJ)=A*(OFFSET+EVCH*JJ)**R RETURN END C C C C C C-----------------------------------------------------------------------C SUBROUTINE BGDMOD C-----------------------------------------------------------------------C C C C MODIFIED 5-15-85 FOR MULTIPLE WINDOWS FOR POLNOMINAL ROUTINES C C SUBROUTINE TO FIT A*E**R BACKGROUND MODEL TO SINGLE WINDOW C SUBROUTINE TO FIT POLYNOMIAL BACKGROUND MODEL TO WINDOWS C SUBROUTINE TO FIT LOG-POLYNOMIAL BACKGROUND MODEL TO WINDOWS C SELECTED BY USER C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ DOUBLE PRECISION SUM(5),X,Y C C C DIMENSION SX(2),SY(2) C C CALCULATE SCALING ARRAY'S FOR DISPLAY C SX(1)=XMIN SX(2)=XMIN + XL/DX SY(1)=YMIN SY(2)=YMIN + YL/DY D5 FORMAT(1X,4(1PE10.3,1X)) C C CHOOSE TYPE OF BACKGROUND FIT C CALL PLOT(0.,180.,0) WRITE(7,10) 10 FORMAT($,' SELECT BGD FIT:A*E**R,POLYNOMIAL,LOG-POLYNOMIAL ') READ(5,20) ANS 20 FORMAT(1A1) IBGND=1 !A*E**R BGD IF(ANS.EQ.'P') IBGND=2 !POLYNOMIAL BGD IF(ANS.EQ.'L') IBGND=3 !LOG-POLYNOMIAL BGD C C READ DATA FROM BX,& BY ARRAYS C C INITIALIZE WINDOW LIMITS TO ZERO SO THAT C SUBROUTINE BGDFIT WILL AS USER FOR LIMITS C IBL=0 IBU=0 C C BY SELECTING WINDOWS IN THE CURRENT DISPLAY C THIS IS NOW DONE IN THE BGDFIT & POLFIT ROUTINES C C DO LINEAR LEAST SQUARES FIT OF ROI IF A*E**R C OTHERWISE DO 6TH ORDER POLYNOMIAL FIT C C IF(IBGND.EQ.1) CALL BGDFIT (IBL,IBU) IF(IBGND.NE.1) CALL POLFIT (IBL,IBU) C C FINISHED FIT WRITE OUT COEFFICIENTS C CALL DELAY(3) C C C CALCULATE START CHANNEL OF BGND FIT IN AX UNITS C STARTING CHANNEL IN AX = EMIN IN BX - EMIN IN AX / EVCH C IAL=(BX(IBL)-AX(1))/EVCH BLAST=BY(IBL-1) DO 300 I=1,NPT-1 AY(I)=0. IF(I.LT.IAL) GO TO 300 !SKIP CALCULATION IF < IAL IF(IBGND.EQ.1) AY(I)=A*AX(I)**R !A*E**R BGD IF(IBGND.EQ.1) GO TO 312 !SKIP REST IF A*E**R BGD EN=BX(IZERO+I-1)/BX(IBL) !SCALE ALL ENERGY TO 1ST POINT IF(IBGND.EQ.3) EN=ALOG(BX(IZERO+I-1)) !EXCEPT IN LOG MODE BG=0. !TEMPORARY BG CALCULATION DO 310 IJ=1,7 !LOOP THROUGH COEFFICIENTS IEXP=1-IJ !EXPONENT IF(IBGND.EQ.3) IEXP=-IEXP !NEGATIVE EXP IF LOG MODE BG=BG+COE(IJ)*EN**(IEXP) !CALCULATE EACH POWER & SUM 310 CONTINUE C IF(BG.GT.BLAST) BG=BLAST !CONSTRAIN TO DECREASING FNCT IF(BG.LT.0) BG=0. !CONSTRAIN TO POSITITVE # BLAST=BG !LAST VALUE OF BGND CALCULATED IF(IBGND.EQ.3) BG=EXP(BG) !LOG MODE MUST EXP TO GET VALUE AY(I)=BG 312 CONTINUE 300 CONTINUE C C DRAW THIS ON EXISTING GRAPH C CALL CRT(AX,AY,NPT-1,1,1,1.,1.,SX,SY,2) C C NOTE ASSUMES THAT NPT IN THIS MODE IS ALWAYS <4110 TO AVOID NPTS RECAL C CALL PLOT(0.,120.,0) CALL DELAY(3) WRITE(7,400) 400 FORMAT($,' CHANGE BGFIT?') READ(5,401) ANS 401 FORMAT(1A1) IF(ANS.EQ.'Y') RETURN WRITE(7,402) 402 FORMAT($,' SUBTRACT BGND?') READ(5,401) ANS IF(ANS.EQ.'Y') IBFLG=1 !BACKGROUND SUBTRACTED FROM DATA IF(ANS.EQ.'N') IBFLG=0 !NO BACKGROUND SUBTRACTION C C OPERATOR HAS CHOZEN TO SUBTRACT BACKGROUND C C SET IBFLG=1 TO INDICATE SUBTRACTION C SET IBGND=1 TO INDICATE A*E**R MODEL C SET IBGND=2 TO INDICATE POLYNOMINAL MODEL C SET IBGND=3 TO INDICATE LOG-POLY MODEL C BLAST=BY(IBL-1) DO 500 I=1,NNPT IF(I.LT.IBL) GO TO 501 IF(IBGND.EQ.1) CCY(I)=A*BX(I)**R IF(IBGND.EQ.1) GO TO 501 EN=BX(I)/BX(IBL) !SCALE ALL ENERGY TO 1ST POINT IF(IBGND.EQ.3) EN=ALOG(BX(I)) !EXCEPT IN LOG MODE BG=0. !TEMPORARY BG CALCULATION DO 510 IJ=1,7 !LOOP THROUGH COEFFICIENTS C IEXP=IJ-7 !EXPONENT IEXP=1-IJ !EXPONENT IF(IBGND.EQ.3) IEXP=-IEXP !NEGATIVE EXP IF LOG MODE BG=BG+COE(IJ)*EN**(IEXP) !CALCULATE EACH POWER & SUM 510 CONTINUE IF(BG.GT.BLAST) BG=BLAST !CONSTRAIN TO DECREASING FNCT IF(BG.LT.0) BG=0. !CONSTRAIN TO POSITITVE # BLAST=BG !LAST VALUE OF BGND CALCULATED IF(IBGND.EQ.3) BG=EXP(BG) !LOG MODE MUST EXP TO GET VALUE CCY(I)=BG 501 IF(ANS.EQ.'Y') BY(I)=BY(I)-CCY(I) 500 CONTINUE C C EXIT ROUTINE C RETURN END C-----------------------------------------------------------------------C SUBROUTINE BGDFIT (IBL,IBU) C-----------------------------------------------------------------------C C C MODIFIED 5-15-85 WINDOW ROUTINE MOVED INTO THIS SUBROUTINE C C SUBROUTINE TO FIT A*E**R BACKGROUND MODEL TO WINDOW C SELECTED BY USER C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ DOUBLE PRECISION SUM(5),X,Y C C LET USER SELECT THE WINDOW FOR THE FIT C IF(IBL.EQ.IBU) CALL WINDOW (IBL,IBU) C C DO LINEAR LEAST SQUARES FIT OF ROI C CALL LSTSQR(1,,,,,,SUM) !USE BRESSLER'S ROUTINE HERE NUMPTS=IBU-IBL+1 DO 200 I = IBL,IBU X=DLOG(DBLE(BX(I))) Y=DLOG(DBLE(BY(I))) CALL LSTSQR(2,X,Y,,,,SUM) 200 CONTINUE CALL LSTSQR(3,,,NUMPTS,R,A,SUM) A=EXP(A) D WRITE(7,300) BX(IBL),BX(IBU),A,R WRITE(8,300) BX(IBL),BX(IBU),A,R 300 FORMAT(' BGND REGION ',1F6.1,' TO ',1F6.1,' EV',/, C ' A = ',1PE10.3,' R = ',0PF6.2) RETURN END C-----------------------------------------------------------------------C SUBROUTINE POLFIT (IBL,IBU) C-----------------------------------------------------------------------C C C THIS ROUTINE FITS UP TO A 6TH ORDER POLYNOMIAL EQUATION C TO A SET OF DATA DEFINED IN several windows C THIS IS A MODIFICATION OF THE ORIGINAL POLFIT ROUTINE C WHICH STARTED USING ONLY A SINGLE WINDOW REGION C C MODIFIED 5-15-85 NJZ C C USING STANDARD DEC ROUTINES C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) C COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C COMMON /TERMIN/ ITERM,IBAUD,IPLOT,ICOLOR,IMODE,IHARD C C C DEFINE ARRAYS NEEDED BY DEC SS ROUTINES C REAL DI(36),D(28) REAL B(6),SB(6),T(6),E(6),ANS(10) DIMENSION XBAR(7),STD(7),SUMSQ(7),ISAVE(7) C C DEFINE WINDOW ARRAYS LOW (WL) AND HIGH (WH) FOR X-AXIS POINTS C MAXIMUM OF 10 WINDOWS ALLOWED C INTEGER WL(10),WH(10) C C BEGIN DATA ANALYSIS & ZERO WINDOWS C 4 DO 5 I=1,10 WL(I)=0. 5 WH(I)=0. C C NOW ASK USER TO INPUT WINDOW C STOP WHEN THERE ARE 10 WINDOWS OR C USER HAS NOT ENTERED A CHANGE FROM THE LAST ONE C NBGPT=0. !NUMBER OF POINTS IN WINDOW NWIND=1 !NUMBER OF WINDOWS DO 6 I=1,10 IF(I.GT.1) WRITE(7,600) IF(I.GT.1) READ (5,601) AANS 600 FORMAT(' ANOTHER WINDOW [N]?') 601 FORMAT(1A1) IF((I.GT.1).AND.(AANS.NE.'Y')) GO TO 8 7 CALL WINDOW(WL(I),WH(I)) !WL & WH ARE CHANNEL NUMBERS IF ((I.EQ.1).AND.(WL(I).EQ.WH(I))) GO TO 7 !ERROR NO WINDOW NBGPT=WH(I)- WL(I) + NBGPT 6 NWIND=NWIND + 1 8 NWIND=NWIND-1 IBL=WL(1) IBU=WH(NWIND) IF (NBGPT.GT.100) WRITE(7,9) !99 POINTS = MAXIMUM 9 FORMAT(' ERROR- TOO MANY POINTS REENTER ALL') IF (NBGPT.GT.100) GO TO 4 !REENTER ALL WINDOWS NJ=0. !COUNTER TO ZERO DO 10 I=1,NWIND !GET DATA FROM WINDOW DO 10 II=WL(I),WH(I) NJ=NJ+1 EX(NJ)=BX(II)/BX(WL(1)) !SCALE DATA TO STOP OVERFLOWS 10 CONTINUE C C ZERO COEFFICIENT ARRAY FOR SAFETY C DO 20 I=1,7 B(I)=0. 20 COE(I)=0. C C C EFIRST=BX(WL(1)) !ENERGY VALUE OF 1ST POINT EX(100)=EFIRST !SAVE FIRST VALUE IN EX(100) FOR DINT ROUTINE INDEX=NJ !# OF DATA POINTS C ------------------------------------------------------------ C Put independent and dependent data in AY array C ------------------------------------------------------------ 105 N = INDEX !Number of observations M = 6 !Highest degree specified L = N * M DO 110 I = 1,N J = L + I ENERGY = EX(I) * EFIRST !Actual energy AY(I) = 1.0 / EX(I) !Inverse powers of X IF( IBGND .EQ. 3 ) AY(I) = ALOG(ENERGY) !No inverse for Log IF( AY(I) .EQ. 0.0 ) AY(I) = 0.01 AY(J)=BY(ICHAN(ENERGY)) !DATA IF( IBGND .EQ. 3 ) AY(J) = ALOG(AY(J)) !Log poly fit 110 CONTINUE C ------------------------------------------------------------ C Generate the other powers etc. of X C ------------------------------------------------------------ CALL GDATA(N,M,AY,XBAR,STD,D,SUMSQ) MM = M + 1 SUM = 0.0 IDEGRE = 0 C ------------------------------------------------------------ C Try all powers of X until no statistical improvement C ------------------------------------------------------------ DO 200 I = 1,M ISAVE(I) = I CALL ORDER(MM,D,MM,I,ISAVE,DI,E) CALL MINV(DI,I,DET,B,T) CALL MULTR(N,I,XBAR,STD,SUMSQ,DI,E,ISAVE,B,SB,T,ANS) IF( ANS(7) .LT. 0.0 ) GO TO 1000 SUMIP = ANS(4) - SUM IF( SUMIP .LE. 0.0 ) GO TO 1000 SUM = ANS(4) COE(1)=ANS(1) ! AFTER DEC POLYRG ROUTINE DO 160 J = 1,I !Save coefficients here COE(J+1)=B(J) !ALA DEC 160 CONTINUE C COE(I+1) = ANS(1) IDEGRE = I !Highest degree used in fit 200 CONTINUE C ------------------------------------------------------------ C List result of fit and allow Printing also C ------------------------------------------------------------ 1000 EMIN = EX(1) * EFIRST EMAX = EX(INDEX) * EFIRST D WRITE(7,500) EMIN,EMAX,(COE(JJ),JJ=7,1,-1) WRITE(8,500) EMIN,EMAX,(COE(JJ),JJ=7,1,-1) 500 FORMAT(' BGD REGION ',1F6.1,' TO ',1F6.1,' EV',/ C ' -6 -5 -4 -3 -2 -1 ',/, C ' AX + BX + CX + DX + EX + FX + G',/, C 10X,'A = ',1PE10.3,/, C 10X,'B = ',1PE10.3,/, C 10X,'C = ',1PE10.3,/, C 10X,'D = ',1PE10.3,/, C 10X,'E = ',1PE10.3,/, C 10X,'F = ',1PE10.3,/, C 10X,'G = ',1PE10.3,/) IF(IBGND.EQ.2) WRITE(8,510) D IF(IBGND.EQ.2) WRITE(7,510) IF(IBGND.EQ.3) WRITE(8,520) D IF(IBGND.EQ.3) WRITE(7,520) 510 FORMAT(' POLYNOMIAL FIT') 520 FORMAT(' LOG-POLYNOMIAL FIT') RETURN END C-----------------------------------------------------------------------C SUBROUTINE DINT C-----------------------------------------------------------------------C C C EDGE INTEGRATION ROUTINE FOR NELS C C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ COMMON /NAME/ INAME(8),EO,CONV,DIV C C C C C CALL ADM5E(6) !ERASE 6 LINES BENEATH GRAPH C C C START COUNTING # OF EDGES INTEGRATED (NEL) C C DFLAG=0. NEL=NEL+1 IF(NEL.GT.20) RETURN 3 CALL PLOT(0.,158.,0) CALL DELAY(5) WRITE(7,1) NEL 1 FORMAT($,' ENTER ATOMIC SYMBOL FOR ELEMENT ',I2,':') READ(5,4) SYM(NEL) WRITE(7,2) 2 FORMAT($,' ENTER EDGE SYM. (K-L-M) :') READ(5,4) AEDID(NEL) 4 FORMAT(1A2,1X,1A2) WRITE(7,5) 5 FORMAT($,' ENTER INTEGRATION WINDOW (EV):') READ(5,6) WIND 6 FORMAT(1F10.0) C C STORE WINDOW WIDTH C ES(NEL)=WIND C C LOCATE EDGE ONSET C WRITE(7,61) 61 FORMAT(' SET CURSOR TO EDGE ONSET') C C GET EDGE ONSET USING CRSSHR C CALL DELAY(7) CALL CRSSPT(X,Y) CALL TRANXY(X,Y,E,C) IWB=ICHAN(E) ESTART=OFFSET+EVCH*IWB EFIN=ESTART+WIND IWF=ICHAN(EFIN) BGD(NEL)=0. !BGD = ZERO C C IF FLAG=-3 THEN DECONVOLUTION C SET ALL BGD SUMS TO ZERO C IF(FLAG.EQ.-3) GO TO 71 !SKIP BGD SUM IF DECON EXISTS C C IF CCY(I)=0. THEN ASSUME HAVE SUBTRACTED BACKGROUND C USE LAST CALCULATED COEFFICIENTS TO DETERMINE C BGND VALUES FOR STATISITICAL ANALYSIS C C DO 7 I=IWB,IWF DFLAG=0. IF(CCY(I).EQ.0.) DFLAG=1. B=CCY(I) E=OFFSET+EVCH*I C C TEST WHAT TYPE OF BGND FIT C THEN USE CORRECT FORMULAE C IF IBGND=1 THEN A*E**R MODEL C IF((B.EQ.0.).AND.(IBGND.EQ.1)) B=A*E**R C C POLYNOMIAL FIT IF IBGND=2 C LOG-POLYNOMINAL FIT IF IBGND=3 C IF((B.EQ.0).AND.(IBGND.GE.2)) GO TO 8 GO TO 9 8 DO 10 J=1,7 IEXP=1-J IF(IBGND.EQ.3) IEXP=-IEXP !LOG POLYNOMINAL FIT IF(IBGND.EQ.3) E=ALOG(E) IF(IBGND.EQ.2) E=E/EX(100) !POLYNOMIAL FIT B=B+COE(J)*E**(IEXP) 10 CONTINUE IF(IBGND.EQ.3) B=EXP(B) !LOG POLYNOMINAL FIT 9 BGD(NEL)=BGD(NEL)+B C C C TOTALIZE SPECTRUM COUNTS C 7 CONTINUE C 71 TOT=SUMMEM(IWB,IWF) C C CHECK FOR PREVIOUS BGND SUBTRACTION C C IF(DFLAG.EQ.1) TOT=TOT+BGD(NEL) C C CALCULATE NET PEAK INTENSITY C C PK(NEL)=TOT-BGD(NEL) PB=0. IF(FLAG.NE.-3) PB=PK(NEL)/BGD(NEL) PT=PK(NEL)/TOT CALL VT100 CALL DELAY(10) IF (ITERM.EQ.8) CALL ENI4027 WRITE(7,13) NEL,ESTART,EFIN,WIND WRITE(8,13) NEL,ESTART,EFIN,WIND 13 FORMAT(//,' INTEGRATION OF ELEMENT ',I2,/, C 32('-'),/' FROM ',1F7.2 C ' TO ',1F7.2,' EV',/, C ' ENERGY WINDOW= ',1F7.2,/) WRITE(7,12) SYM(NEL),AEDID(NEL),PB,PT, CTOT,BGD(NEL),PK(NEL) WRITE(8,12) SYM(NEL),AEDID(NEL),PB,PT, CTOT,BGD(NEL),PK(NEL) 12 FORMAT(1X,/, C ' ELEMENT =',8X,1A2,/, C ' EDGE =',9X,1A2,/, C ' PEAK/BGND =',1F10.3,/, C ' PEAK/TOT =',1F10.3,/, C ' TOTAL INT =',1PE10.3,/, C ' BGND INT =',1PE10.3,/, C ' PEAK INT =',1PE10.3,/) WRITE(7,14) 14 FORMAT($,' CHANGE INT. WINDOW?') READ(5,15) ANS 15 FORMAT(1A1) CALL ERASE IF(ANS.EQ.'N') RETURN BGD(NEL)=0. PK(NEL)=0. B=0. C CALL CRT( AX,AY,NPT,1,0,1.,1.,AX,AY,NPT) GO TO 3 16 CONTINUE RETURN END C C C C INTEGRATION FUNCTION FOR NELS C C C C-----------------------------------------------------------------------C FUNCTION SUMMEM(IWB,IWF) C-----------------------------------------------------------------------C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) SUMMEM=0. DO 1 I=IWB,IWF 1 SUMMEM=SUMMEM + BY(I) RETURN END C-----------------------------------------------------------------------C SUBROUTINE CONDIT(BY) C-----------------------------------------------------------------------C C C THIS ROUTINE REQUESTS OPERATING CONDITIONS OF MICROSCOPE C FOR STORAGE IN APPROPRIATE LOCATION OF DATA ARRAY C C C C COMMON BLOCKS FOR NELS C C C COMMON AX(4110),AY(4110),BX(4110),BY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT C COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) C COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL C COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ COMMON /NAME/ INAME(8),EO,CONV,DIV DIMENSION BY(4110) C C C C IF(BY(NNPT+3).EQ.0) WRITE(7,1) IF(BY(NNPT+3).EQ.0) READ (5,2) BY(NNPT+3) 1 FORMAT($,' ENTER INCIDENT BEAM ENERGY (KEV) : ') 2 FORMAT(1F10.0) IF(BY(NNPT+4).EQ.0) WRITE (7,3) 3 FORMAT($,' ENTER INCIDENT BEAM DIVERGENCE (MR): ') IF(BY(NNPT+4).EQ.0) READ(5,2) BY(NNPT+4) IF(BY(NNPT+5).EQ.0) WRITE(7,4) 4 FORMAT($,' ENTER SCATTERING ANGLE (MR) : ') IF(BY(NNPT+5).EQ.0) READ(5,2) BY(NNPT+5) C C NOW STORE NEW VALUES IN COMMON C EO=BY(NNPT+3) CONV=BY(NNPT+4) DIV=BY(NNPT+5) RETURN END C-----------------------------------------------------------------------C SUBROUTINE CURSOR C-----------------------------------------------------------------------C C C C SUBROUTINE TO WRITE VALUE OF DATA POINT SELECTED C BY CURSOR ON THE DISPLAY TERMINAL C OR TO CALCULATE THE AREA IN A USER SELECTED WINDOW C C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C C CALL ADM5E(7) !ERASE 7 LINES C C C C C CALL PLOT(0.,114.,0) CALL DELAY(3) C C SELECT OPTION C WRITE(7,10) 10 FORMAT($,' SELECT OPTION:VALUE AT CURSOR,AREA OF WINDOW:') READ(5,20) ANS 20 FORMAT(1A1) IF(ANS.EQ.'A') GO TO 30 C C USER WANTS VALUE OF CURSOR C WRITE (7,1) CALL DELAY(5) 1 FORMAT(' SET CURSOR TO DATA POINT:') CALL CRSSPT(X,Y) CALL TRANXY(X,Y,EN,CT) CALL PLOT(490.,92.,0) CALL DELAY(3) WRITE(7,2) EN,BY(ICHAN(EN)) CALL DELAY(3) 2 FORMAT(' EN=',1F9.3,' CTS=',1PE10.3) RETURN C C USER WANTS AREA IN WINDOW C 30 CALL AREA RETURN END C-----------------------------------------------------------------------C SUBROUTINE AREA C-----------------------------------------------------------------------C C C THIS SUBROUTINE COMPUTES THE AREA IN THE ROI C DEFINED BY THE USER (I.E. SUMS # OF COUNTS IN A WINDOW) C AND DOESNOT PERFORM ANY BACKGROUND CALCULATIONS C IT IS MAINLY USED FOR ZERO LOSS INTEGRATIONS C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C C C C C C SELECT WINDOW C WRITE (7,100) 100 FORMAT(' Area defined by Cursor Limits (L), or Width (W): [L] ') READ (5,110) ANS 110 FORMAT (1A1) IF (ANS.NE.'W') CALL WINDOW (IL,IU) IF (ANS.NE.'W') GO TO 120 C C SET START POINT AND WIDTH C WRITE (7,1130) 1130 FORMAT($,' ENTER START ENERGY (0=CURSOR)') READ (5,112) EN IF(EN.NE.0) GO TO 1111 WRITE (7,113) 113 FORMAT (' SET CURSOR TO STARTING POINT') CALL DELAY(5) CALL CRSSPT(X,Y) CALL TRANXY(X,Y,EN,CNT) CALL PLOT(490.,92.,0) 1111 WRITE (7,111) 111 FORMAT($,' Enter Width in EV :') READ (5,112) WIDEEN 112 FORMAT (1F10.0) IL=ICHAN(EN) IU=IL + WIDEEN/EVCH C C SET WAREA=0 & SUM COUNTS C 120 WAREAABS=0. WAREA=0. AMAXW=0. AMINW=1.E10 C C NOW FIND SIMPLE INTEGRAT WAREA C ABSOLUTE AREA WAREAABS C MIN AMINW C MAX AMAXW C C DO 10 I=IL,IU IF(BY(I).LE.AMINW) AMINW=BY(I) IF(BY(I).GE.AMAXW) AMAXW=BY(I) WAREAABS=WAREAABS + ABS(BY(I)) 10 WAREA=WAREA + BY(I) C C COMPUTE ENERGY LIMITS FOR OUTPUT C EL=OFFSET + EVCH*IL EU=OFFSET + EVCH*IU C C WRITE IT OUT C CALL PLOT(0.,49.,0) CALL DELAY(3) WAREA=WAREA*EVCH WAREAABS=WAREAABS*EVCH DELW = EU-EL DEL=AMAXW-AMINW WRITE(7,20) EL,EU,DELW,WAREA,WAREAABS,AMINW,AMAXW,DEL WRITE(8,20) EL,EU,DELW,WAREA,WAREAABS,AMINW,AMAXW,DEL 20 FORMAT(1X, C ' WINDOW = ',1F6.1,' TO ',1F6.1,' (EV) = ',1F7.2,/, C ' SIMPLE AREA = ',1PE12.4,' COUNTS*EVCH',/, C ' ABSOLUTE AREA = ',1PE12.4,' COUNTS*EVCH',/, C ' MIN,MAX,DEL = ',3 (1PE12.4,1X),' COUNTS') RETURN END C-----------------------------------------------------------------------C SUBROUTINE CALIB C-----------------------------------------------------------------------C C C C CALIBRATION SUBROUTINE FOR HORIZONTAL AXIS C C ALLOWS ENTER,LINEAR,SHIFT,RESCALE OPTIONS C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C C C CALL ADM5E(6) !ERASE 6 LINES C C C C C C 10 CONTINUE C10 CALL PLOT(0.,100.,0) WRITE(7,1) 1 FORMAT($,' SELECT OPTION:ENTER,LINEAR,RESCALE,SHIFT,TYPE,BACK:') READ(5,2) RES 2 FORMAT(1A1) IF(RES.EQ.'T') GO TO 150 IF(RES.EQ.'E') GO TO 250 IF(RES.EQ.'L') GO TO 200 IF(RES.EQ.'S') GO TO 5000 IF(RES.EQ.'R') GO TO 4000 IF(RES.EQ.'B') GO TO 999 GO TO 10 C C C LINEAR CALIBRATION SECTION C C 200 CONTINUE WRITE(7,2000) CALL DELAY(3) 2000 FORMAT($,' CURSOR TO LOW ENERGY :') C C TURN ON CURSOR & CALCULATE CHANNEL NUMBER OF ARRAY C CALL CRSSPT(X,Y) CALL TRANXY(X,Y,EL,Y) IL=ICHAN(EL) CALL PLOT(512.,80.,0) CALL DELAY(3) WRITE(7,3000) CALL DELAY(3) 3000 FORMAT($,' ENTER VALUE:') READ(5,2001)EL 2001 FORMAT(F12.0) WRITE(7,2002) CALL DELAY(3) 2002 FORMAT($,' CURSOR TO HIGH ENERGY ') CALL CRSSPT(X,Y) CALL TRANXY(X,Y,EH,Y) IH=ICHAN(EH) CALL PLOT(512.,60.,0) CALL DELAY(3) WRITE(7,3000) CALL DELAY(3) READ(5,2001)EH EVCH = (EH-EL)/(IH-IL) OFFSET = EL - IL*EVCH 210 WRITE(7,2003) OFFSET,EVCH CALL DELAY(3) 2003 FORMAT(' OFFSET = ',F10.4,' EVCH = ',F10.4) GO TO 999 C C C ENTER CALIBRATION DATA HERE C C C C OPERATOR WANTS TO TYPE IN VALUES C 150 WRITE(7,160) 160 FORMAT($,' ENTER OFFSET(EV):') READ(5,2001) OFFSET WRITE(7,2501) READ(5,2001) EVCH GO TO 210 C C OPERATOR WANTS TO INPUT EV/CH & USE CURSOR TO CAL. OFFSET C 250 WRITE(7,2500) CALL DELAY(3) 2500 FORMAT($,' CURSOR TO KNOWN POSITION ') CALL CRSSPT(X,Y) CALL TRANXY(X,Y,EK,Y) CALL PLOT(512.,80.,0) CALL DELAY(3) WRITE(7,3000) CALL DELAY(3) IK=ICHAN(EK) READ(5,2001) EK WRITE(7,2501) CALL DELAY(3) 2501 FORMAT($,' ENTER EV/CH :') READ(5,2001)EVCH OFFSET=EK-IK*EVCH GO TO 210 C C C CALCULATION OF RESCALE OF SPECTRUM C C 4000 CONTINUE C WRITE(7,4001) CALL DELAY(3) 4001 FORMAT($,' IS CURRENT CALIB. CORRECT?') READ(5,4002) 4002 FORMAT(1A1) IF(ANS.EQ.'N') GO TO 10 WRITE(7,4003) CALL DELAY(4) 4003 FORMAT(' SET CURSOR TO LOC. OF RESCALED SPECTRUM ') CALL CRSSPT(X,Y) CALL TRANXY(X,Y,EL,Y) INEWCH=ICHAN(EL) CALL PLOT(0.,80.,0) WRITE(7,3000) CALL DELAY(3) READ(5,4004) EN 4004 FORMAT(1F10.0) WRITE(7,2501) READ(5,4004) EVNEW CALL DELAY(3) OFFNEW=EN-EVNEW*INEWCH CALL RESCAL(OFFNEW,EVNEW) GO TO 210 C C C CALCULATE SHIFTED SPECTRUM C C 5000 CONTINUE CALL SHIFT GO TO 210 C C C EXIT ROUTINE AND STORE RESULTS C C 999 CONTINUE DO 6000 I=1,NNPT 6000 BX(I)=OFFSET+EVCH*I BY(NNPT+1)=OFFSET BY(NNPT+2)=EVCH AY(NNPT+1)=BY(NNPT+1) AY(NNPT+2)=BY(NNPT+2) RETURN END C-----------------------------------------------------------------------C SUBROUTINE SMOOTH C-----------------------------------------------------------------------C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C C C C C C DOES A 3 POINT BINOMIAL SMOOTH OF C DATA IN CURRENTLY DISPLAYED MEMORY C C C C C C SKIP FIRST TWO & LAST TWO DATA POINTS C DO 1 I=2,NPT-1 C C TEMPORAIRLY STORE SMOOTHED DATA IN AX THEN RECALCULATE C AX(I)=0.5*BY(I) +0.25*BY(I-1) +0.25*BY(I+1) 1 CONTINUE DO 2 I=2,NNPT BY(I)=AX(I) 2 AX(I)=OFFSET + EVCH*I C C NOTE FIRST CALCULATE ENTIRE ARRAY C CAN CAUSE ERRORS DUE TO USING ROUNDED DATA C RATHER THAN RAW DATA TO COMPUTE SMOOTH C C CALLING SMOOTH N TIMES IS EQUIVALENT TO A 2N+1 C DATA POINT SMOOTHING C RETURN END C-----------------------------------------------------------------------C SUBROUTINE DIGFIL(N) C-----------------------------------------------------------------------C C C 8407200001-NJZ C 8408070002-NJZ ADD TRIANGULAR FILTER N=6 & GENERAL CLEAN UP C 8408070003-NJZ ADD 4TH DERIVITATIVE TOPHAT FILTER AS N=7 C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C C C THIS ROUTINE PERFORMS A DIGITAL FILTER OF THE CURRENTLY C DISPLAYED MEMORY THE CURRENT VERSION USES THE FOLLOWING C ALOGRITHM TO FILTER THE SPECTRUM: C C FOR N= 1 TO 5 A TOPHAT FILTER IS USED C FOR N=6 TRIANGULAR FILTER WITH NEGATIVE LOBES C FOR N=7 MODIFIED TOPHAT 2 POSITIVE ONE NEGATIVE LOBES C C IF N=1 THEN DF(I)= -1*A(I-1) + 2*A(I) -1*A(I+1) OR C THIS IS A HIGH RESOLUTION FILTER BUT WILL BE C OF LITTLE USE IN A NOISEY SPECTRUM C IF N=2 THEN DF(I) = EQUIVALENT OF N=1 BUT 3 CHANNELS WIDE C THIS IS A MEDIUM RESOLUTION FILTER AND IN C A NOISEY SPECTRUM IMPROVES THE NOISE SUPPRESSION C DF(I)=-1*(A(I-4)+A(I-3)+A(I-2)) +2*(A(I-1)+A(I)+A(I+1)) C -1*(A(I+2)+A(I+3)+A(I+4)) C IF N=3 THEN DF(I) = C -1.*(BY(I-7)+BY(I-6)+BY(I-5)+BY(I-4)+BY(I-3)+BY(I-2)) C +4.*(BY(I-1)+BY(I)+BY(I+1)) C -1.*(BY(I+7)+BY(I+6)+BY(I+5)+BY(I+4)+BY(I+3)+BY(I+2)) C THIS IS A LOWER RESOLUTION FILTER WHICH IMPROVES NOISE C SUPPRESSION OVER THAT OF 2 BUT ALSO DECREASES THE RESOLUTION C IF N=4 THEN DF(I) = C -1.*(BY(I-10)+BY(I-9)+BY(I-8)+BY(I-7)+BY(I-6)+BY(I-5) C +BY(I-4)+BY(I-3)+BY(I-2))+6.*(BY(I-1)+BY(I)+BY(I+1)) C -1.*(BY(I+10)+BY(I+9)+BY(I+8)+BY(I+7)+BY(I+6)+BY(I+5) C +BY(I+4)+BY(I+3)+BY(I+2)) C THIS WILL BE THE LOWEST RESOLUTION FILTER C IF N>5 THEN USER IS ASKED TO ENTER HIS OWN TOPHAT FILTER FUNCTION C IN TERMS OF THE WIDTH OF THE NEGATIVE LOBES (NLOWER) C AND THE WIDTH OF THE POSITIVE LOBE (NUPPER) C NOTE: THE POSITIVE LOBE MUST HAVE AN ODD NUMBER OF CHANNELS C C C WHERE BY(I) IS THE VALUE STORED IN MEMORY LOCATION I OF ARRAY "BY" C C C C SKIP FIRST & LAST DATA POINTS AS APPROPRIATE C NNEG=1 !RECTANGULAR TOPHAT FILTER N=1 DEFAULT NPOS=1 IF(N.EQ.2) NNEG=3 IF(N.EQ.2) NPOS=3 IF(N.EQ.3) NNEG=6 IF(N.EQ.3) NPOS=3 IF(N.EQ.4) NNEG=9 IF(N.EQ.4) NPOS=3 ANEG=FLOAT(NNEG) APOS=FLOAT(NPOS) IF(N.LE.4) GO TO 60 C C USER HAS ELECTED TO ENTER A NONSTANDARD TOPHAT FILTER WIDTH C C READ IN NEG AND POS LOBE WIDTHS THEN COMPUTE ZERO AREA FUNCTION C C IF(N.LE.6) WRITE(7,20) IF(N.EQ.7) WRITE(7,21) 20 FORMAT($,' ENTER WIDTH OF POS. & NEG. LOBES IN CHANNELS C (POS. MUST BE ODD !):') 21 FORMAT($,' ENTER WIDTH OF CENTRAL, NEG, & POS LOBES: ') IF(N.LE.6) READ (5,30) APOS,ANEG IF(N.EQ.7) READ (5,30) APOS,ANEG,APOS2 30 FORMAT(4F10.0) NNEG=INT(ANEG) NPOS=INT(APOS) NPOS2=INT(APOS2) C C COMPUTE ZERO AREA TERMS I.E. AREA OF NEG= -1(AREA OF POS) C 60 IF(ANEG.EQ.0) ANEG=1. !NO ZEROS ALLOWED FOR AFACT IF(APOS.EQ.0) APOS=1. !IBID ATR IF(N.LE.6) AFACT=2.*ANEG/APOS IF(N.EQ.7) AFACT=(2.*ANEG-2.*APOS2)/APOS ATR=ANEG/(APOS+1) C C DEFINE START AND END POINT TO DO FILTERING C I.E. MAKE SURE WE ARE IN LIMITS OF DATA ARRAYS C ISTRT=NNEG+(NPOS+1)/2 IEND=NPT-NNEG-NPOS/2 IF(N.EQ.7) ISTRT=ISTR+NPOS2 IF(N.EQ.7) IEND=IEND-NPOS2 C C C LET'S DO IT C C 10 DO 1 I=ISTRT,IEND C C TEMPORARILY STORE SMOOTHED DATA IN AX THEN RECALCULATE C C NEGATIVE LOBES FIRST C TEMP1=0. IF(NNEG.EQ.0) GO TO 41 DO 40 J=1,NNEG TEMP1=TEMP1 -1.*BY(I-J-NPOS/2) -1.*BY(I+J+NPOS/2) 40 CONTINUE C C NOW THE POSITIVE LOBES C 41 TEMP2=0. DO 50 J=-NPOS/2,NPOS/2 C C CHECK FOR TOPHAT OR TRIANGLE FUNCTION C IF(N.LE.5)TEMP2=TEMP2 + AFACT*BY(I-J) IF(N.EQ.7)TEMP2=TEMP2 + AFACT*BY(I-J) IF(N.EQ.6)TEMP2=TEMP2 + ATR*((NPOS+1)/2 -IABS(J))*BY(I-J) 50 CONTINUE C C ADDITIONAL LOOP FOR N=7 CASE C TEMP3=0. IF(N.LE.6) GO TO 52 IF(NPOS2.EQ.0) GO TO 52 DO 51 J=1,NPOS2 TEMP3=TEMP3 +1.*BY(I-J-NNEG-NPOS/2) +1.*BY(I+J+NNEG+NPOS/2) 51 CONTINUE 52 AX(I)=TEMP1 + TEMP2 + TEMP3 C C DONE WITH THIS POINT CONTINUE THE LOOP C 1 CONTINUE C C C FINISHED THE FILTERING NOW C STORE ORIGINAL SPECTRUM IN MEMORY #2 'CCY' C STORE DIGITAL FILTERED SPECTRUM IN MEMORY#1 'BY' C COPY DIGITAL FILTERED SPECTRUM INTO DISPLAY MEMORY 'AY' C COPY X AXIS CALIBRATION INTO 'AX' C DO 2 I=1,NNPT CCY(I)=BY(I) BY(I)=AX(I) AY(I)=AX(I) 2 AX(I)=OFFSET + EVCH*I DO 200 I=NNPT+1,NNPT+10 200 CCY(I)=BY(I) !COPY CALIBRATION CONSTANTS TOO C C LET EVERYBODY KNOW WHAT'S HAPPENING TO THE DATA C C READ(5,13) RES 13 FORMAT(1A2) WRITE(7,1000) 1000 FORMAT(' DIG. FIL. IN MEM#1,ORIG. DATA IN MEM#2: TO CONTINUE') RETURN END C-----------------------------------------------------------------------C SUBROUTINE MENU C-----------------------------------------------------------------------C C C SUBROUTINE TO LIST OUT MENU DEFINITIONS C CALL VT100 WRITE(7,1) 1 FORMAT(/, C ' GR = GRAPH MODIFICATION OF DISPLAY',/, C ' HC = HARCOPY OF DISPLAY ON HP7470 PLOTTER',/, C ' WI = WINDOW SET AND EXPAND TO FILL DISPLAY',/, C ' CU = CURSOR VALUE',/, C ' SM = SMOOTH DATA',/, C ' CA = CALIBRATE SPECTRUM',/, C ' NO = NORMALIZE (REMOVE) GAIN CHANGES',/, C ' BG = BACKGROUND FITTING',/, C ' IN = INTEGRATE EDGE',/, C ' DE = DECONVOLUTION',/, C ' ME = MEMORY & DATA MANIPULATION',/, C ' RE = READ DATA FILE FROM DISK',/, C ' WR = WRITE DATA TO FILE ON DISK',/, C ' OR = ORIGINAL DATA REDRAWN',/, C ' EX = EXIT PROGRAM',/) WRITE(7,1000) 1000 FORMAT (' ENTER TO CONTINUE') READ(5,13) RES 13 FORMAT(1A2) WRITE(7,2) 2 FORMAT(/, C ' THE FOLLOWING MODES ARE AVAILABLE IN THE "ME" OPTION',/, C ' CO = COPY FROM MEM #1 TO 2 OR 3',/, C ' DF = DIGITAL FILTER OF MEM #1',/, C ' SSM= SPECTRUM-SPECTRUM MATH +,-,X,/ MEM 1 BY MEM 2',/, C ' SCM= SPECTRUM-CONSTANT MATH " " " " MEM 1 BY CONSTANT',/, C ' LOG= LOG OF DISPLAY',/, C ' AL = ANTILOG OF DISPLAY',/, C ' OV = OVERLAY MEM1/MEM2 ETC.',/, C ' ZE = SETS ALL VALUES <0. IN CURRENT MEMORY = TO 0.',/, C ' DEL= DELETE PORTION OF CURRENT MEMORY',/, C ' EQ = EQUALIZE (SCALE) 2 MEMORIES AT CURSOR LOCATION',//) WRITE(7,1000) READ(5,13) RES RETURN END C-----------------------------------------------------------------------C SUBROUTINE NELOUT C-----------------------------------------------------------------------C C C C C THIS SUBROUTINE CREATES AN OUTPUT DATA FILE C SUITABLE FOR INPUT FOR THE NELQNT DATA ANALYSIS PROGRAM C CONTAINING ALL THE RELEVENT OUTPUT FROM THE LAST C ANALYSIS PERFORMED BY THE DATA ANALYSIS PROGRAM NELS C C OUTPUT FILE IS WRITTEN TO DEVICE ASSIGNED TO THE LOGICAL C NAME WD: (WRITE DEVICE) WITH THE NAME WD:NELOUT.DAT C C ------------ C DEFINITIONS: C ------------ C C C C NEL= # OF ELEMENTS INTEGRATED C SYM= ARRAY CONTAINING ATOMIC SYMBOLS OF INTEGRATED ELEMENTS C AEDID= ARRAY CONTAINING EDGE ID`S " " " C PK= ARRAY CONTAINING EDGE INTEGRALS C BGD= ARRAY CONTAINING BGND INTEGRALS C ES = ARRAY CONTAINING WINDOW WIDTHS FOR INTEGRALS IN EV C C C COMMON BLOCKS FOR NELS C C C COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT C COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL,AREAZL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /TERMIN/ ITERM,IBAUD,IPLOT,ICOLOR,IMODE,IHARD COMMON /NAME/ INAME(8),EO,CONV,DIV BYTE DAY(9),TIM(8) C C C OPEN THE DEFAULT DATA FILE ONLY IF SOME ELEMENTS WERE ANALYZED C C IF(NEL.LE.0) RETURN C C THERE WAS AT LEAST ONE ELEMENT INTEGRATED C REWIND 8 !FINISH ALL HARDCOPY OUTPUT BEFORE OPENING FILE CLOSE (UNIT=8,DISP='SAVE') !CLOSE HARDCOPY DEVICE C C C OPEN(UNIT=30,TYPE='UNKNOWN',ACCESS='SEQUENTIAL', C NAME='WD:NELOUT.DAT',DISP='SAVE',FORM='FORMATTED',ERR=100) C C C REWIND 30 !REWIND FILE FOR NEW INPUT CALL DATE (DAY) CALL TIME (TIM) IEO=INT(EO) WRITE(30,10) DAY,TIM WRITE(7,10) DAY,TIM 10 FORMAT(' NELS/8504260019-NJZ Data Analysis on ',9A1,' at ',8A1) WRITE(30,20) (INAME(J),J=1,7), NEL, IEO, CONV, DIV WRITE(7,20) (INAME(J),J=1,7), NEL, IEO, CONV, DIV 20 FORMAT( C ' INPUT DATA FILE NAME = ',7A2,/, C ' NUMBER OF ELEMENTS ANALYZED = ',I4,/, C ' ACCELERATING VOLTAGE (KEV) = ',I4,/, C ' INCIDENT BEAM CONVERGENCE (MR)= ',1F6.3,/, C ' SCATTERING ANGLE (MR) = ',1F6.3,/, C ' ----------------------------------------------',/, C ' ELMNT LINE EDGE INT. BGND INT. WINDOW (EV)',/, C ' ----------------------------------------------') IF (NEL.EQ.0) GO TO 60 DO 30 I=1,NEL WRITE(30,40) SYM(I),AEDID(I),PK(I),BGD(I),ES(I) WRITE(7,40) SYM(I),AEDID(I),PK(I),BGD(I),ES(I) 40 FORMAT(1X,1A2,3X,1A2,1X,1PE12.4,1X,1PE12.4,2X,0PF6.2) 30 CONTINUE 60 CLOSE (UNIT=30, DISP='SAVE') RETURN 100 WRITE(7,50) 50 FORMAT(' Error in attempt to open NELOUT.DAT file') RETURN END C-----------------------------------------------------------------------C SUBROUTINE CFFT(IDIR,A,M) C-----------------------------------------------------------------------C C C THIS ROUTINE IS FROM M. DISKO AT A.S.U. C C C COMPLEX FOURIER TRANSFORM ROUTINE IN ONE DIMENSION C IDIR=1 : FORWARD TRANSFORM C IDIR=-1: INVERSE TRANSFORM C THE NUMBER OF POINTS = 2**M C COMPLEX A(1),U,W,T DATA PI /3.14159265/ IDIR=IDIR/IABS(IDIR) N=2**M NV2=N/2 DIR=IDIR NM1=N-1 J=1 DO 7 I=1, NM1 IF(I.GE.J) GOTO 5 T=A(J) A(J)=A(I) A(I)=T 5 K=NV2 6 IF(K.GE.J) GOTO 7 J=J-K K=K/2 GOTO 6 7 J=J+K DO 20 L=1,M LE=2**L LE1=LE/2 U=(1.0,0.) W=CMPLX(COS(PI/LE1),DIR*SIN(PI/LE1)) DO 20 J=1, LE1 DO 10 I=J, N, LE IP=I+LE1 T=A(IP)*U A(IP)=A(I)-T 10 A(I)=A(I)+T 20 U=U*W IF(IDIR.EQ.-1) RETURN ANORM=1./FLOAT(N) DO 30 I=1, N 30 A(I)=A(I)*ANORM RETURN END C-----------------------------------------------------------------------C SUBROUTINE DECONP(RES,IZLC,ILLS,ILLF,IHLS,IHLF) C-----------------------------------------------------------------------C C PLOTTING SECTION FOR DECONL ROUTINE C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),AHL(515),ALL(515) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL,AREAZL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ DIMENSION SX(2),SY(2) COMPLEX ALL,AHL C CALCULATE OLD SCALE FACTORS C C 131 SX(1)=XMIN SX(2)=XMIN + XL/DX SY(1)=YMIN SY(2)=YMIN + YL/DY C C FILL UP AX & AY WITH NEW DATA C C C FIND MAXIMUM VALUE IN ALL(I) ARRAY FOR SCALING C DMAX=0. DO 135 I=1,512 IF((RES.EQ.'H').OR.(RES.EQ.'D')) DMAX=AMAX1(DMAX,CABS(AHL(I))) IF(RES.EQ.'L') DMAX=AMAX1(DMAX,CABS(ALL(I))) 135 CONTINUE NP=512 IF(RES.EQ.'D') NP=IHLF-IHLS IF(RES.EQ.'L') IST=IZLC IF((RES.EQ.'H').OR.(RES.EQ.'D')) IST=IHLS DO 140 I=1,512 AX(I)=BX(IST+I-1) C C MULTIPLY BY SCALE FACTOR SO DECONVOLUTED DATA C IS VISIBLE ON DISPLAY C SFACT=SY(2)/DMAX IF((RES.EQ.'D').OR.(RES.EQ.'H')) AY(I)=SFACT*CABS(AHL(I)) IF(RES.EQ.'L') AY(I)=SFACT*CABS(ALL(I)) 140 CONTINUE C C TELL OPERATOR SCALE FACTOR C CALL DELAY(5) IF(RES.EQ.'D') WRITE (7,150) SFACT IF(RES.EQ.'D') WRITE (8,150) SFACT CALL DELAY(5) 150 FORMAT(' FFT DECON SCALING FACTOR = ',1PE10.4) C C DRAW IT!! C CALL CRT( AX,AY,NP,1,1,1.,1.,SX,SY,2) C C DIVIDE OUT SCALE FACTOR C MULTIPLY BY AREA OF ZERO LOSS FOR CONSTANT SCALING C DO 160 I=1,512 160 AY(I)=AY(I)*AREAZL/SFACT RETURN END C-----------------------------------------------------------------------C SUBROUTINE DELETE C-----------------------------------------------------------------------C C C C THIS SUBROUTINE ALLOWS THE OPERATOR TO SELECTIVELY C DELETE DATA FROM THE CURRENTLY DISPLAYED SPECTRUM C C I.E. SET DATA VALUES TO ZERO C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ COMMON /TERMIN/ ITERM,IBAUD,IPLOT,ICOLOR,IMODE,IHARD C C C C C ALLOW USER TO SET WINDOW AROUND REGION TO BE DELETED C CALL WINDOW(IL,IU) C C SET DATA VALUES TO ZERO C DO 10 I=IL,IU BY(I)=0. 10 CONTINUE RETURN END C-----------------------------------------------------------------------C FUNCTION ICHAN(E) C-----------------------------------------------------------------------C C C COVERTS ENERGY TO CHANNEL NUMBER C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C C C C ICHAN=(E-OFFSET)/EVCH RETURN END C-----------------------------------------------------------------------C SUBROUTINE SHIFT C-----------------------------------------------------------------------C C C SUBROUTINE TO SHIFT A SPECTRUM ALONG ENERGY AXIS C WITH ROLL AROUND C C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C C WRITE(7,5000) 5000 FORMAT($,' SHIFT BY ENERGY OR CHANNEL?:') READ (5,5005) ANS 5005 FORMAT(1A1) WRITE(7,5001) 5001 FORMAT($,' ENTER SHIFT :') READ(5,4004) EN 4004 FORMAT(1F10.0) IF(ANS.EQ.'C') NCHN=EN IF(ANS.EQ.'C') EN=NCHN*EVCH IF(ANS.EQ.'C') GO TO 5006 VAL=(EN)/EVCH NCHN = VAL + SIGN (0.5,VAL) 5006 WRITE(7,5002) EN,NCHN 5002 FORMAT(' ENERGY SHIFT = ',1F8.3, ' EV', C ' NUMBER OF CHANNELS = ',I5) DO 5003 I=1,NNPT J=I-NCHN IF(J.LE.0) J=4110-J IF(J.GT.4110) J=J-4110 BX(J)=BY(I) 5003 CONTINUE OFFSET=OFFSET+EN DO 5004 I=1,NNPT BY(I)=BX(I) 5004 BX(I)=OFFSET + EVCH*I 999 RETURN END C-----------------------------------------------------------------------C SUBROUTINE SSM C-----------------------------------------------------------------------C C C SPECTRUM CONSTANT MATH ROUTINE C OPERATES ONLY ON MEMORY #1 (BY ARRAY) C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C C C C 10 WRITE (7,20) 20 FORMAT($,'MEM#1MEM#2 = MEM#1: SELECT OPERATION: +,-,X,/ :') READ(5,30) ANS 30 FORMAT(1A1) DO 60 I=1,NNPT VAL=CCY(I) C C DON'T ALLOW MULTIPLICATION OR DIVISION BY ZERO C IF((VAL.EQ.0.).AND.((ANS.EQ.'X').OR.(ANS.EQ.'/'))) VAL=1. IF(ANS.EQ.'+') BY(I)=BY(I) + VAL IF(ANS.EQ.'-') BY(I)=BY(I) - VAL IF(ANS.EQ.'X') BY(I)=BY(I)*VAL IF(ANS.EQ.'/') BY(I)=BY(I)/VAL 60 CONTINUE RETURN END C-----------------------------------------------------------------------C SUBROUTINE SCM C-----------------------------------------------------------------------C C C SPECTRUM CONSTANT MATH ROUTINE C OPERATES ONLY ON MEMORY #1 (BY ARRAY) C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C C C C 10 WRITE (7,20) 20 FORMAT($,' SELECT MEM#1 OPERATION: +,-,X,/ :') READ(5,30) ANS 30 FORMAT(1A1) WRITE(7,40) 40 FORMAT($,' ENTER CONSTANT:') READ(5,50) VAL 50 FORMAT(1F10.0) IF((VAL.EQ.0.).AND.((ANS.EQ.'X').OR.(ANS.EQ.'/'))) VAL=1. C C DON'T ALLOW MULTIPLICATION OR DIVISION BY ZERO C DO 60 I=1,NNPT IF(ANS.EQ.'+') BY(I)=BY(I) + VAL IF(ANS.EQ.'-') BY(I)=BY(I) - VAL IF(ANS.EQ.'X') BY(I)=BY(I)*VAL IF(ANS.EQ.'/') BY(I)=BY(I)/VAL 60 CONTINUE RETURN END C-----------------------------------------------------------------------C SUBROUTINE RESCAL(OFFNEW,EVNEW) C-----------------------------------------------------------------------C C C THIS SUBROUTINE USES A CUBIC SPLINE FITTING ROUTINE C TO INTERPOLATE POINTS INORDER TO EXPAND OR CONTRACT C A GIVEN SPECTRUM TO A DESIRED EV SPACING C C C C COMMON BLOCKS FOR NELS C C C C NOTE COMMON ARRAYS ALTERED HERE C COMMON AX(4110),AY(4110),BX(4110),BY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C C C C COMMON E(512),DI(512),COEFF(512),DDB(514),BX(4110) C,BY(4110),CCY(4110),DDY(4110) C C C CALCULATION OF RESCALE OF SPECTRUM C C C INEW=1 C C SUBTRACT 10 CHANNELS FROM CHOSEN START C SO THAT ERRORS WON'T OCCUR AT START POSITION C C ISTRT=(OFFNEW-OFFSET)/EVCH - 10 IF(ISTRT.LE.0) ISTRT=0 C C CALCULATE ABSOLUTE MAX OF ENERGY SCALES C EABS=OFFSET+EVCH*4110 DO 4006 I=1,512 E(I)=OFFSET+EVCH*(I+ISTRT) DI(I)=BY(I+ISTRT) 4006 CONTINUE CALL SPLIN1(512) EMAX=OFFSET+EVCH*(512+ISTRT) DO 4007 I=1,512 DD=0. EN=OFFNEW +EVNEW*I IF(EN.LE.OFFSET) GO TO 4007 IF(EN.GE.EMAX) GO TO 4007 DD=SPLIN2(512,EN) BX(INEW)=DD INEW=INEW+1 4007 CONTINUE C C STORE RESCALED DATA IN BX UNTIL FINISHED C THEN CALL BACK TO PRIMARY MEMORY C C C C FINISHED FIRST 512 CHANNEL SEGMENT C NOW START ON REMAINING 512 CHANNELS C ELAST= ENERGY OF LAST CHANNEL RESCALED C ILAST= CHANNEL NUMBER OF LAST CHANNEL RESCALED C BOTH OF THESE ARE IN OLD UNITS (EFFECTIVELY) ALTHOUGH C ELAST IS CALCULATED USING NEW SCALE FACTORS C C ALSO CHECK TO SEE THAT WE HAVEN'T EXCEEDED LIMITS OF C SPECTRUM ( EMAX= MAX ENERGY OF OLD SPECTRUM IN CURRENT C 512 CHANNEL CHUNK C C ELAST=OFFNEW+EVNEW*512 IF(ELAST.GT.EMAX) ELAST=EMAX 4010 ILAST=(ELAST-OFFSET)/EVCH 4011 EMAX=OFFSET+EVCH*(ILAST+512) DO 4008 I=1,512 E(I)=OFFSET+EVCH*(I+ILAST) DI(I)=0. IF((I+ILAST).GE.4110) GO TO 4008 DI(I)=BY(I+ILAST) 4008 CONTINUE CALL SPLIN1(512) DO 4009 I=1,512 DD=0. EN=EVNEW*I + ELAST IF(EN.LE.OFFSET) GO TO 4009 IF(EN.GT.EMAX) GO TO 4012 DD=SPLIN2(512,EN) BX(INEW)=DD INEW=INEW+1 IF(EN.GT.EABS) GO TO 4014 IF(INEW.GT.4110) GO TO 4014 4009 CONTINUE C C RECOMPUTE LAST ENERGY C RECOMPUTE NEW LAST CHANNEL (ILAST) C REDRAW NEW DATA C CHECK TO SEE IF MAX # OF CHANNELS IS CONVERTED (4110) C C C IF NOT DO NEXT 512 CHUNK OF DATA C C 4012 ELAST=EN-EVNEW ILAST=(ELAST-OFFSET)/EVCH IF(ILAST.LE.1010) GO TO 4011 C C STORE NEW CALIBRATIONS & RESTORE ARRAY'S BX,BY WITH RESCALED VALUES C C 4014 OFFSET=OFFNEW EVCH=EVNEW DO 4013 I=1,NNPT BY(I)=BX(I) BX(I)=OFFSET + I*EVCH 4013 CONTINUE 999 RETURN END C-----------------------------------------------------------------------C SUBROUTINE COPY C-----------------------------------------------------------------------C C C SUBROUTINE TO COPY BETWEEN MEMORIES 1,2,3 C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C C C 10 WRITE(7,20) 20 FORMAT($,' SELECT OPTION:1>2,1>3,2>3:') READ(5,30) ANS 30 FORMAT(1A3) C C COPY SPECTRA BETWEEN MEMORIES C DO 40 I=1,NNPT C C SET DEFAULTS FOR COPY C BN=BY(I) CN=CCY(I) DN=DDY(I) IF(ANS.EQ.'1>2') BN=CCY(I) IF(ANS.EQ.'1>2') CN=BY(I) IF(ANS.EQ.'1>3') BN=DDY(I) IF(ANS.EQ.'1>3') DN=BY(I) IF(ANS.EQ.'2>3') DN=CCY(I) IF(ANS.EQ.'2>3') CN=DDY(I) BY(I)=BN CCY(I)=CN DDY(I)=DN 40 CONTINUE OFFSET=BY(NNPT+1) EVCH=BY(NNPT+2) RETURN END C-----------------------------------------------------------------------C SUBROUTINE WINDOW(IL,IU) C-----------------------------------------------------------------------C C C SUBROUTINE TO LET USER SET A WINDOW IN CURRENT DISPLAY C USER SELECTS MODE OF WINDOW SELECTION BY INPUTING VALUES C OR USING CURSOR. ROUTINE RETURNS CHANNEL VALUES FOR C IL= LOWER CHANNEL C IU= UPPER CHANNEL C OF SELECTED WINDOW C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C C C SELECT TYPE OF INPUT VALUE OR CURSOR C C CALL PLOT(0.,136.,0) WRITE(7,31) 31 FORMAT($,1X,'ENTER LIMITS:EMIN,EMAX (CRSSPT=) ') READ(5,32) EMIN,EMAX 32 FORMAT(2F10.0) IF((EMIN.EQ.0).AND.(EMAX.EQ.0.)) GO TO 35 C C USER HAS SELECTED VALUE METHOD C IL=INT((EMIN-OFFSET)/EVCH) IU=INT((EMAX-OFFSET)/EVCH) CALL DELAY(2) RETURN C C USER HAS SELECTED CRSSHR INPUT C 35 CALL PLOT(0.,92.,0) CALL DELAY(5) C WRITE(7,36) CALL DELAY(5) 36 FORMAT(' SET CURSOR TO LOWER LIMIT') CALL CRSSPT(X1,Y1) CALL TRANXY(X1,Y1,XPT1,YPT1) CALL PLOT(392.,92.,0) CALL DELAY(5) C WRITE(7,37) CALL DELAY(5) 37 FORMAT(' SET CURSOR TO UPPER LIMIT') CALL CRSSPT(X2,Y2) CALL TRANXY(X2,Y2,XPT2,YPT2) IL=INT((XPT1-OFFSET)/EVCH) IU=INT((XPT2-OFFSET)/EVCH) RETURN END C C .................................................................. C C-----------------------------------------------------------------------C C SUBROUTINE GDATA C-----------------------------------------------------------------------C C C PURPOSE C GENERATE INDEPENDENT VARIABLES UP TO THE M-TH POWER (THE C HIGHEST DEGREE POLYNOMIAL SPECIFIED) AND COMPUTE MEANS, C STANDARD DEVIATIONS, AND CORRELATION COEFFICIENTS. THIS C SUBROUTINE IS NORMALLY CALLED BEFORE SUBROUTINES ORDER, C MINV AND MULTR IN THE PERFORMANCE OF A POLYNOMIAL C REGRESSION. C C USAGE C CALL GDATA (N,M,X,XBAR,STD,D,SUMSQ) C C DESCRIPTION OF PARAMETERS C N - NUMBER OF OBSERVATIONS. C M - THE HIGHEST DEGREE POLYNOMIAL TO BE FITTED. C X - INPUT MATRIX (N BY M+1) . WHEN THE SUBROUTINE IS C CALLED, DATA FOR THE INDEPENDENT VARIABLE ARE C STORED IN THE FIRST COLUMN OF MATRIX X, AND DATA FOR C THE DEPENDENT VARIABLE ARE STORED IN THE LAST C COLUMN OF THE MATRIX. UPON RETURNING TO THE C CALLING ROUTINE, GENERATED POWERS OF THE INDEPENDENT C VARIABLE ARE STORED IN COLUMNS 2 THROUGH M. C XBAR - OUTPUT VECTOR OF LENGTH M+1 CONTAINING MEANS OF C INDEPENDENT AND DEPENDENT VARIABLES. C STD - OUTPUT VECTOR OF LENGTH M+1 CONTAINING STANDARD C DEVIATIONS OF INDEPENDENT AND DEPENDENT VARIABLES. C D - OUTPUT MATRIX (ONLY UPPER TRIANGULAR PORTION OF THE C SYMMETRIC MATRIX OF M+1 BY M+1) CONTAINING CORRELA- C TION COEFFICIENTS. (STORAGE MODE OF 1) C SUMSQ - OUTPUT VECTOR OF LENGTH M+1 CONTAINING SUMS OF C PRODUCTS OF DEVIATIONS FROM MEANS OF INDEPENDENT C AND DEPENDENT VARIABLES. C C REMARKS C N MUST BE GREATER THAN M+1. C IF M IS EQUAL TO 5 OR GREATER, SINGLE PRECISION MAY NOT BE C SUFFICIENT TO GIVE SATISFACTORY COMPUTATIONAL RESULTS. C C SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED C NONE C C METHOD C REFER TO B. OSTLE, 'STATISTICS IN RESEARCH', THE IOWA STATE C COLLEGE PRESS, 1954, CHAPTER 6. C C .................................................................. C SUBROUTINE GDATA (N,M,X,XBAR,STD,D,SUMSQ) DIMENSION X(1),XBAR(1),STD(1),D(1),SUMSQ(1) C C ............................................................... C C IF A DOUBLE PRECISION VERSION OF THIS ROUTINE IS DESIRED, THE C C IN COLUMN 1 SHOULD BE REMOVED FROM THE DOUBLE PRECISION C STATEMENT WHICH FOLLOWS. C C DOUBLE PRECISION X,XBAR,STD,D,SUMSQ,T1,T2 C C THE C MUST ALSO BE REMOVED FROM DOUBLE PRECISION STATEMENTS C APPEARING IN OTHER ROUTINES USED IN CONJUNCTION WITH THIS C ROUTINE. C C THE DOUBLE PRECISION VERSION OF THIS SUBROUTINE MUST ALSO C CONTAIN DOUBLE PRECISION FORTRAN FUNCTIONS. SQRT AND ABS IN C STATEMENT 180 MUST BE CHANGED TO DSQRT AND DABS. C C ............................................................... C C GENERATE INDEPENDENT VARIABLES C IF(M-1) 105, 105, 90 90 L1=0 DO 100 I=2,M L1=L1+N DO 100 J=1,N L=L1+J K=L-N 100 X(L)=X(K)*X(J) C C CALCULATE MEANS C 105 MM=M+1 DF=N L=0 DO 115 I=1,MM XBAR(I)=0.0 DO 110 J=1,N L=L+1 110 XBAR(I)=XBAR(I)+X(L) 115 XBAR(I)=XBAR(I)/DF C DO 130 I=1,MM 130 STD(I)=0.0 C C CALCULATE SUMS OF CROSS-PRODUCTS OF DEVIATIONS C L=((MM+1)*MM)/2 DO 150 I=1,L 150 D(I)=0.0 DO 170 K=1,N L=0 DO 170 J=1,MM L2=N*(J-1)+K T2=X(L2)-XBAR(J) STD(J)=STD(J)+T2 DO 170 I=1,J L1=N*(I-1)+K T1=X(L1)-XBAR(I) L=L+1 170 D(L)=D(L)+T1*T2 L=0 DO 175 J=1,MM DO 175 I=1,J L=L+1 175 D(L)=D(L)-STD(I)*STD(J)/DF L=0 DO 180 I=1,MM L=L+I SUMSQ(I)=D(L) 180 STD(I)= SQRT( ABS(D(L))) C C CALCULATE CORRELATION COEFFICIENTS C L=0 DO 190 J=1,MM DO 190 I=1,J L=L+1 190 D(L)=D(L)/(STD(I)*STD(J)) C C CALCULATE STANDARD DEVIATIONS C DF=SQRT(DF-1.0) DO 200 I=1,MM 200 STD(I)=STD(I)/DF RETURN END C-----------------------------------------------------------------------C SUBROUTINE LSTSQR(ICOM,X,Y,NUMPTS,SLOPE,YINTER,SUM) C-----------------------------------------------------------------------C C C$ FILENAME: ELLTSQ.FOR C$ PROJECT: PV9100/60 C$ SUB-PROJECT: EELS C$ VERSION: 1.0 C$ ENVIRONMENT: 2.2 C$ CREATED: 15-OCT-81 BY W. Bresler C$ EDIT LOG: C$ NO. DATE BY REASON C$ 1 25-AUG-82 WB Use double precision for better accuracy C 2 15-FEB-83 NJZ FOR USE ON GENERAL RT-11 ROUTINE WITH NO CLRA FNCT C$ DOCUMENTATION: C C ********************************************************************** C * C * This routine will perform least squares fitting C * to a linear function. C * C * Parameters: C * Inputs: C * ICOM - the command parameter: C * 1- zeros the summations C * 2- adds a point to the summations C * 3- computes the slope and y-intercept C * 4- computes the correlation coeff. C * NUMPTS - the number of points in the fit. C * needed only when ICOM = 3,4 C * X - the value of x in an (x,y) pair C * Y - the value of y in an (x,y) pair C * SLOPE - the returned value input to compute C * the corr. coeff. C * SUM - the array used to hold the summations C * C * Outputs: C * SLOPE - the slope of the fit; returned when C * ICOM = 3. C * YINTER - the y-intercept of the fit returned C * when ICOM = 3; C * the corr. coeff of the fit returned C * when ICOM = 4 C * C * Errors: C * NUMPTS = 1 C * Result - SLOPE = 0 ; YINTER = y C * Action - none necessary C * C ********************************************************************** C C### DEFINITIONS C DOUBLE PRECISION SUM(5),X,Y,DENOM,SIGMAX,SIGMAY C C### CODE C C Zero Sigma's Fit Corr coeff GO TO (100, 200, 300, 400)ICOM C ------------------------------------------------------------ C ZERO ALL SUMMATIONS C ------------------------------------------------------------ 100 DO 101 I=1,5 101 SUM(I)=0. GO TO 500 !Return point C ------------------------------------------------------------ C ADD A POINT TO THE FIT C ------------------------------------------------------------ 200 SUM(1) = SUM(1) + X !SIGMA(X'S) SUM(2) = SUM(2) + Y !SIGMA(Y'S) SUM(3) = SUM(3) + X * X !SIGMA(X*X) SUM(4) = SUM(4) + X * Y !SIGMA(X*Y) SUM(5) = SUM(5) + Y * Y !SIGMA(Y*Y) GO TO 500 C ------------------------------------------------------------ C COMPUTE THE FIT PARAMETERS C ------------------------------------------------------------ 300 IF ( NUMPTS .GT. 1 ) GOTO 310 SLOPE = 0 ! Single point fit YINTER = SUM(2) GO TO 500 310 DENOM = NUMPTS * SUM(3) - SUM(1) * SUM(1) ! Multiple point fit SLOPE = (SUM(4) * NUMPTS - SUM(1) * SUM(2))/DENOM YINTER = (SUM(3) * SUM(2) - SUM(1) * SUM(4))/DENOM GO TO 500 C ------------------------------------------------------------ C COMPUTE THE CORRELATION COEFFICIENT C ------------------------------------------------------------ 400 SIGMAX = DSQRT( SUM(3) / NUMPTS - ( SUM(1) / NUMPTS ) ** 2 ) SIGMAY = DSQRT( SUM(5) / NUMPTS - ( SUM(2) / NUMPTS ) ** 2 ) YINTER = SLOPE * SIGMAX / SIGMAY !RETURNS VALUE IN YINTER GO TO 500 C ------------------------------------------------------------ C NORMAL RETURN PATH C ------------------------------------------------------------ 500 RETURN END C C .................................................................. C C-----------------------------------------------------------------------C C SUBROUTINE MINV C-----------------------------------------------------------------------C C C PURPOSE C INVERT A MATRIX C C USAGE C CALL MINV(A,N,D,L,M) C C DESCRIPTION OF PARAMETERS C A - INPUT MATRIX, DESTROYED IN COMPUTATION AND REPLACED BY C RESULTANT INVERSE. C N - ORDER OF MATRIX A C D - RESULTANT DETERMINANT C L - WORK VECTOR OF LENGTH N C M - WORK VECTOR OF LENGTH N C C REMARKS C MATRIX A MUST BE A GENERAL MATRIX C C SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED C NONE C C METHOD C THE STANDARD GAUSS-JORDAN METHOD IS USED. THE DETERMINANT C IS ALSO CALCULATED. A DETERMINANT OF ZERO INDICATES THAT C THE MATRIX IS SINGULAR. C C .................................................................. C SUBROUTINE MINV(A,N,D,L,M) DIMENSION A(1),L(1),M(1) C C ............................................................... C C IF A DOUBLE PRECISION VERSION OF THIS ROUTINE IS DESIRED, THE C C IN COLUMN 1 SHOULD BE REMOVED FROM THE DOUBLE PRECISION C STATEMENT WHICH FOLLOWS. C C DOUBLE PRECISION A,D,BIGA,HOLD,DABS C C THE C MUST ALSO BE REMOVED FROM DOUBLE PRECISION STATEMENTS C APPEARING IN OTHER ROUTINES USED IN CONJUNCTION WITH THIS C ROUTINE. C C THE DOUBLE PRECISION VERSION OF THIS SUBROUTINE MUST ALSO C CONTAIN DOUBLE PRECISION FORTRAN FUNCTIONS. ABS IN STATEMENT C 10 MUST BE CHANGED TO DABS. C C ............................................................... C C SEARCH FOR LARGEST ELEMENT C D=1.0 NK=-N DO 80 K=1,N NK=NK+N L(K)=K M(K)=K KK=NK+K BIGA=A(KK) DO 20 J=K,N IZ=N*(J-1) DO 20 I=K,N IJ=IZ+I 10 IF( ABS(BIGA)- ABS(A(IJ))) 15,20,20 15 BIGA=A(IJ) L(K)=I M(K)=J 20 CONTINUE C C INTERCHANGE ROWS C J=L(K) IF(J-K) 35,35,25 25 KI=K-N DO 30 I=1,N KI=KI+N HOLD=-A(KI) JI=KI-K+J A(KI)=A(JI) 30 A(JI) =HOLD C C INTERCHANGE COLUMNS C 35 I=M(K) IF(I-K) 45,45,38 38 JP=N*(I-1) DO 40 J=1,N JK=NK+J JI=JP+J HOLD=-A(JK) A(JK)=A(JI) 40 A(JI) =HOLD C C DIVIDE COLUMN BY MINUS PIVOT (VALUE OF PIVOT ELEMENT IS C CONTAINED IN BIGA) C 45 IF(BIGA) 48,46,48 46 D=0.0 RETURN 48 DO 55 I=1,N IF(I-K) 50,55,50 50 IK=NK+I A(IK)=A(IK)/(-BIGA) 55 CONTINUE C C REDUCE MATRIX C DO 65 I=1,N IK=NK+I HOLD=A(IK) IJ=I-N DO 65 J=1,N IJ=IJ+N IF(I-K) 60,65,60 60 IF(J-K) 62,65,62 62 KJ=IJ-I+K A(IJ)=HOLD*A(KJ)+A(IJ) 65 CONTINUE C C DIVIDE ROW BY PIVOT C KJ=K-N DO 75 J=1,N KJ=KJ+N IF(J-K) 70,75,70 70 A(KJ)=A(KJ)/BIGA 75 CONTINUE C C PRODUCT OF PIVOTS C D=D*BIGA C C REPLACE PIVOT BY RECIPROCAL C A(KK)=1.0/BIGA 80 CONTINUE C C FINAL ROW AND COLUMN INTERCHANGE C K=N 100 K=(K-1) IF(K) 150,150,105 105 I=L(K) IF(I-K) 120,120,108 108 JQ=N*(K-1) JR=N*(I-1) DO 110 J=1,N JK=JQ+J HOLD=A(JK) JI=JR+J A(JK)=-A(JI) 110 A(JI) =HOLD 120 J=M(K) IF(J-K) 100,100,125 125 KI=K-N DO 130 I=1,N KI=KI+N HOLD=A(KI) JI=KI-K+J A(KI)=-A(JI) 130 A(JI) =HOLD GO TO 100 150 RETURN END C C .................................................................. C C-----------------------------------------------------------------------C C SUBROUTINE MULTR C-----------------------------------------------------------------------C C C PURPOSE C PERFORM A MULTIPLE LINEAR REGRESSION ANALYSIS FOR A C DEPENDENT VARIABLE AND A SET OF INDEPENDENT VARIABLES. THIS C SUBROUTINE IS NORMALLY USED IN THE PERFORMANCE OF MULTIPLE C AND POLYNOMIAL REGRESSION ANALYSES. C C USAGE C CALL MULTR (N,K,XBAR,STD,D,RX,RY,ISAVE,B,SB,T,ANS) C C DESCRIPTION OF PARAMETERS C N - NUMBER OF OBSERVATIONS. C K - NUMBER OF INDEPENDENT VARIABLES IN THIS REGRESSION. C XBAR - INPUT VECTOR OF LENGTH M CONTAINING MEANS OF ALL C VARIABLES. M IS NUMBER OF VARIABLES IN OBSERVATIONS. C STD - INPUT VECTOR OF LENGTH M CONTAINING STANDARD DEVI- C ATIONS OF ALL VARIABLES. C D - INPUT VECTOR OF LENGTH M CONTAINING THE DIAGONAL OF C THE MATRIX OF SUMS OF CROSS-PRODUCTS OF DEVIATIONS C FROM MEANS FOR ALL VARIABLES. C RX - INPUT MATRIX (K X K) CONTAINING THE INVERSE OF C INTERCORRELATIONS AMONG INDEPENDENT VARIABLES. C RY - INPUT VECTOR OF LENGTH K CONTAINING INTERCORRELA- C TIONS OF INDEPENDENT VARIABLES WITH DEPENDENT C VARIABLE. C ISAVE - INPUT VECTOR OF LENGTH K+1 CONTAINING SUBSCRIPTS OF C INDEPENDENT VARIABLES IN ASCENDING ORDER. THE C SUBSCRIPT OF THE DEPENDENT VARIABLE IS STORED IN C THE LAST, K+1, POSITION. C B - OUTPUT VECTOR OF LENGTH K CONTAINING REGRESSION C COEFFICIENTS. C SB - OUTPUT VECTOR OF LENGTH K CONTAINING STANDARD C DEVIATIONS OF REGRESSION COEFFICIENTS. C T - OUTPUT VECTOR OF LENGTH K CONTAINING T-VALUES. C ANS - OUTPUT VECTOR OF LENGTH 10 CONTAINING THE FOLLOWING C INFORMATION.. C ANS(1) INTERCEPT C ANS(2) MULTIPLE CORRELATION COEFFICIENT C ANS(3) STANDARD ERROR OF ESTIMATE C ANS(4) SUM OF SQUARES ATTRIBUTABLE TO REGRES- C SION (SSAR) C ANS(5) DEGREES OF FREEDOM ASSOCIATED WITH SSAR C ANS(6) MEAN SQUARE OF SSAR C ANS(7) SUM OF SQUARES OF DEVIATIONS FROM REGRES- C SION (SSDR) C ANS(8) DEGREES OF FREEDOM ASSOCIATED WITH SSDR C ANS(9) MEAN SQUARE OF SSDR C ANS(10) F-VALUE C C REMARKS C N MUST BE GREATER THAN K+1. C C SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED C NONE C C METHOD C THE GAUSS-JORDAN METHOD IS USED IN THE SOLUTION OF THE C NORMAL EQUATIONS. REFER TO W. W. COOLEY AND P. R. LOHNES, C 'MULTIVARIATE PROCEDURES FOR THE BEHAVIORAL SCIENCES', C JOHN WILEY AND SONS, 1962, CHAPTER 3, AND B. OSTLE, C 'STATISTICS IN RESEARCH', THE IOWA STATE COLLEGE PRESS, C 1954, CHAPTER 8. C C .................................................................. C SUBROUTINE MULTR (N,K,XBAR,STD,D,RX,RY,ISAVE,B,SB,T,ANS) DIMENSION XBAR(1),STD(1),D(1),RX(1),RY(1),ISAVE(1),B(1),SB(1), 1 T(1),ANS(1) C C ............................................................... C C IF A DOUBLE PRECISION VERSION OF THIS ROUTINE IS DESIRED, THE C C IN COLUMN 1 SHOULD BE REMOVED FROM THE DOUBLE PRECISION C STATEMENT WHICH FOLLOWS. C C DOUBLE PRECISION XBAR,STD,D,RX,RY,B,SB,T,ANS,RM,BO,SSAR,SSDR,SY, C 1 FN,FK,SSARM,SSDRM,F,DSQRT,DABS C C THE C MUST ALSO BE REMOVED FROM DOUBLE PRECISION STATEMENTS C APPEARING IN OTHER ROUTINES USED IN CONJUNCTION WITH THIS C ROUTINE. C C THE DOUBLE PRECISION VERSION OF THIS SUBROUTINE MUST ALSO C CONTAIN DOUBLE PRECISION FORTRAN FUNCTIONS. SQRT AND ABS IN C STATEMENTS 122, 125, AND 135 MUST BE CHANGED TO DSQRT AND DABS. C C ............................................................... C MM=K+1 C C BETA WEIGHTS C DO 100 J=1,K 100 B(J)=0.0 DO 110 J=1,K L1=K*(J-1) DO 110 I=1,K L=L1+I 110 B(J)=B(J)+RY(I)*RX(L) RM=0.0 BO=0.0 L1=ISAVE(MM) C C COEFFICIENT OF DETERMINATION C DO 120 I=1,K RM=RM+B(I)*RY(I) C C REGRESSION COEFFICIENTS C L=ISAVE(I) B(I)=B(I)*(STD(L1)/STD(L)) C C INTERCEPT C 120 BO=BO+B(I)*XBAR(L) BO=XBAR(L1)-BO C C SUM OF SQUARES ATTRIBUTABLE TO REGRESSION C SSAR=RM*D(L1) C C MULTIPLE CORRELATION COEFFICIENT C 122 RM= SQRT( ABS(RM)) C C SUM OF SQUARES OF DEVIATIONS FROM REGRESSION C SSDR=D(L1)-SSAR C C VARIANCE OF ESTIMATE C FN=N-K-1 SY=SSDR/FN C C STANDARD DEVIATIONS OF REGRESSION COEFFICIENTS C DO 130 J=1,K L1=K*(J-1)+J L=ISAVE(J) 125 SB(J)= SQRT( ABS((RX(L1)/D(L))*SY)) C C COMPUTED T-VALUES C 130 T(J)=B(J)/SB(J) C C STANDARD ERROR OF ESTIMATE C 135 SY= SQRT( ABS(SY)) C C F VALUE C FK=K SSARM=SSAR/FK SSDRM=SSDR/FN F=SSARM/SSDRM C ANS(1)=BO ANS(2)=RM ANS(3)=SY ANS(4)=SSAR ANS(5)=FK ANS(6)=SSARM ANS(7)=SSDR ANS(8)=FN ANS(9)=SSDRM ANS(10)=F RETURN END C C .................................................................. C C-----------------------------------------------------------------------C C SUBROUTINE ORDER C-----------------------------------------------------------------------C C C PURPOSE C CONSTRUCT FROM A LARGER MATRIX OF CORRELATION COEFFICIENTS C A SUBSET MATRIX OF INTERCORRELATIONS AMONG INDEPENDENT C VARIABLES AND A VECTOR OF INTERCORRELATIONS OF INDEPENDENT C VARIABLES WITH DEPENDENT VARIABLE. THIS SUBROUTINE IS C NORMALLY USED IN THE PERFORMANCE OF MULTIPLE AND POLYNOMIAL C REGRESSION ANALYSES. C C USAGE C CALL ORDER (M,R,NDEP,K,ISAVE,RX,RY) C C DESCRIPTION OF PARAMETERS C M - NUMBER OF VARIABLES AND ORDER OF MATRIX R. C R - INPUT MATRIX CONTAINING CORRELATION COEFFICIENTS. C THIS SUBROUTINE EXPECTS ONLY UPPER TRIANGULAR C PORTION OF THE SYMMETRIC MATRIX TO BE STORED (BY C COLUMN) IN R. (STORAGE MODE OF 1) C NDEP - THE SUBSCRIPT NUMBER OF THE DEPENDENT VARIABLE. C K - NUMBER OF INDEPENDENT VARIABLES TO BE INCLUDED C IN THE FORTHCOMING REGRESSION. K MUST BE GREATER C THAN OR EQUAL TO 1. C ISAVE - INPUT VECTOR OF LENGTH K+1 CONTAINING, IN ASCENDING C ORDER, THE SUBSCRIPT NUMBERS OF K INDEPENDENT C VARIABLES TO BE INCLUDED IN THE FORTHCOMING REGRES- C SION. C UPON RETURNING TO THE CALLING ROUTINE, THIS VECTOR C CONTAINS, IN ADDITION, THE SUBSCRIPT NUMBER OF C THE DEPENDENT VARIABLE IN K+1 POSITION. C RX - OUTPUT MATRIX (K X K) CONTAINING INTERCORRELATIONS C AMONG INDEPENDENT VARIABLES TO BE USED IN FORTH- C COMING REGRESSION. C RY - OUTPUT VECTOR OF LENGTH K CONTAINING INTERCORRELA- C TIONS OF INDEPENDENT VARIABLES WITH DEPENDENT C VARIABLES. C C REMARKS C NONE C C SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED C NONE C C METHOD C FROM THE SUBSCRIPT NUMBERS OF THE VARIABLES TO BE INCLUDED C IN THE FORTHCOMING REGRESSION, THE SUBROUTINE CONSTRUCTS THE C MATRIX RX AND THE VECTOR RY. C C .................................................................. C SUBROUTINE ORDER (M,R,NDEP,K,ISAVE,RX,RY) DIMENSION R(1),ISAVE(1),RX(1),RY(1) C C ............................................................... C C IF A DOUBLE PRECISION VERSION OF THIS ROUTINE IS DESIRED, THE C C IN COLUMN 1 SHOULD BE REMOVED FROM THE DOUBLE PRECISION C STATEMENT WHICH FOLLOWS. C C DOUBLE PRECISION R,RX,RY C C THE C MUST ALSO BE REMOVED FROM DOUBLE PRECISION STATEMENTS C APPEARING IN OTHER ROUTINES USED IN CONJUNCTION WITH THIS C ROUTINE. C C ............................................................... C C COPY INTERCORRELATIONS OF INDEPENDENT VARIABLES C WITH DEPENDENT VARIABLE C MM=0 DO 130 J=1,K L2=ISAVE(J) IF(NDEP-L2) 122, 123, 123 122 L=NDEP+(L2*L2-L2)/2 GO TO 125 123 L=L2+(NDEP*NDEP-NDEP)/2 125 RY(J)=R(L) C C COPY A SUBSET MATRIX OF INTERCORRELATIONS AMONG C INDEPENDENT VARIABLES C DO 130 I=1,K L1=ISAVE(I) IF(L1-L2) 127, 128, 128 127 L=L1+(L2*L2-L2)/2 GO TO 129 128 L=L2+(L1*L1-L1)/2 129 MM=MM+1 130 RX(MM)=R(L) C C PLACE THE SUBSCRIPT NUMBER OF THE DEPENDENT C VARIABLE IN ISAVE(K+1) C ISAVE(K+1)=NDEP RETURN END C C C C-----------------------------------------------------------------------C SUBROUTINE SPLIN1(N) C-----------------------------------------------------------------------C C C SPLIN SUBROUTINE FROM J.GEORGOPOLIS & G. KNAPP C COMMON X(512),Y(512),COEFF(512),DD(514),BX(4110),BY(4110) C,A(512),B(512),C(512),YY(514) COEFF(1)=0. COEFF(N)=0. IF(N.EQ.2) RETURN NM2=N-2 DO 10 I=1,NM2 H0=X(I) H1=X(I+1) H2=X(I+2) S1=(Y(I+2)-Y(I+1))/(H2-H1) S0=(Y(I+1)-Y(I))/(H1-H0) YY(I)=6.*(S1-S0)/(H2-H0) A(I)=2. B(I)=(H1-H0)/(H2-H0) C(I)=(H2-H1)/(H2-H0) 10 CONTINUE YY(1)=YY(1)/A(1) COEFF(2)=YY(1) IF(N.EQ.3) RETURN C(1)=C(1)/A(1) DO 20 I=2,NM2 A(I)=A(I)-B(I)*C(I-1) C(I)=C(I)/A(I) YY(I)=(YY(I)-B(I)*YY(I-1))/A(I) 20 CONTINUE COEFF(NM2+1)=YY(NM2) DO 30 II=2,NM2 I=NM2-II+1 COEFF(I+1)=YY(I)-C(I)*COEFF(I+2) 30 CONTINUE RETURN END C C C C-----------------------------------------------------------------------C FUNCTION SPLIN2(N,UNK) C-----------------------------------------------------------------------C C C ALSO FROM GEORGOPOLIS & KNAPP C COMMON X(512),Y(512),COEFF(512) I=1 J=2 IF(UNK.LE.X(1)) GO TO 102 I=N J=N-1 IF(UNK.GE.X(N)) GO TO 102 DO 100 I=1,N J=I+1 IF(UNK.LE.X(J)) GO TO 101 100 CONTINUE 101 RI=(UNK-X(I))/(X(J)-X(I)) RJ=(X(J)-UNK)/(X(J)-X(I)) SPLIN2=RJ*Y(I)+RI*Y(J)-(UNK-X(I))*(X(J)-UNK)* 1(COEFF(I)*(RJ+1)+COEFF(J)*(RI+1))/6. RETURN 102 SPLIN2=((Y(J)-Y(I))/(X(J)-X(I))-(X(J)-X(I))* 1COEFF(J)/6.)*(UNK-X(I))+Y(I) RETURN END C-----------------------------------------------------------------------C SUBROUTINE LDIGFIL(N) C-----------------------------------------------------------------------C C C 8407200001-NJZ C 8702290004-NJZ MAJOR CLEAN UP DELETE OLD OPTIONS ADD LOG MODE C C C C COMMON BLOCKS FOR NELS C C COMMON AX(4110),AY(4110),BX(4110),BY(4110),CCY(4110),DDY(4110) COMMON /PARAM0/ NMODE,IRPLT,XM,YM,YLAST,FLAG,NNPT,NPT COMMON /PARAM1/ OFFSET,EVCH,A,R,YY,YN,EX(100),COE(7) COMMON /PARAM2/ IBFLG,IZERO,IBGND,IROI,NEL COMMON /PARAM3/ SYM(20),AEDID(20),ES(20),PK(20),BGD(20) COMMON /CRTLST/ XMIN,DX,YMIN,DY,OFFX,OFFY,XL,YL,XZ,YZ C C C THIS ROUTINE PERFORMS A DIGITAL FILTER OF THE CURRENTLY C DISPLAYED MEMORY THE CURRENT VERSION USES THE FOLLOWING C ALOGRITHM TO FILTER THE SPECTRUM: C C C SKIP FIRST & LAST DATA POINTS AS APPROPRIATE C C C C N=1 REGULAR 1ST ORDER TOPHAT FILTER C N=2 2ND ORDER TOPHAT FILTER C ANEG=1 !RECTANGULAR TOPHAT FILTER N=1 DEFAULT APOS=1 APOS2=0 C C READ IN NEG AND POS LOBE WIDTHS C C WRITE (7,1001) 1001 FORMAT ($, ' TOPHAT OR LOGTOPHAT (T,L) [T] ') READ (5,1002) ANS 1002 FORMAT (1A1) IF(ANS.EQ.' ') ANS = 'T' IF(N.LE.1.OR.ANS.EQ.'L') WRITE(7,20) IF(N.EQ.2.AND.ANS.EQ.'T') WRITE(7,21) 20 FORMAT($,' ENTER WIDTH OF POS. & NEG. LOBES IN CHANNELS C (POS. MUST BE ODD !):') 21 FORMAT($,' ENTER WIDTH OF CENTRAL, NEG, & POS LOBES: ') IF(N.LE.1) READ (5,30) APOS,ANEG IF(N.EQ.2) READ (5,30) APOS,ANEG,APOS2 30 FORMAT(4F10.0) NNEG=INT(ANEG) NPOS=INT(APOS) NPOS2=INT(APOS2) C C DEFINE START AND END POINT TO DO FILTERING C I.E. MAKE SURE WE ARE IN LIMITS OF DATA ARRAYS C ISTRT=NNEG+(NPOS+1)/2 IEND=NPT-NNEG-NPOS/2 IF(N.EQ.2) ISTRT=ISTRT+NPOS2 IF(N.EQ.2) IEND=IEND-NPOS2 C C C LET'S DO IT C C 10 DO 110 I=ISTRT,IEND C C TEMPORARILY STORE SMOOTHED DATA IN AX THEN RECALCULATE C C C C IF USER SPECIFIES WIDTH GREATER THAN ONE CHANNEL THEN C COMPUTE AN AVERAGE VALUE BASED UPON A WINDOW SUM C TFX0=1.E-9 TFX1=1.E-9 TFX2=1.E-9 TFX3=1.E-9 TFX4=1.E-9 C C BEGIN BY COMPUTING THE AVERAGE VALUE IN THE POSITIVE LOBE C DO 111 IJ=-INT(APOS/2),INT(APOS/2) TFX1=TFX1 + BY(I-IJ) 111 CONTINUE C C NOW COMPUTE AVERAGE OF NEGATIVE LOBES C DO 112 IJ=1,INT(ANEG) TFX0=TFX0 + BY(I-INT(APOS/2)-IJ) TFX2=TFX2 + BY(I+INT(APOS/2)+IJ) 112 CONTINUE C C IF 2ND ORDER FILTER ALSO NEED ONE MORE POSITIVE LOBE C IF (APOS2.EQ.0) GO TO 114 DO 113 IJ=1,INT(APOS2) TFX3=TFX3 + BY(I-INT(APOS/2+ANEG)-IJ) TFX4=TFX4 + BY(I+INT(APOS/2+ANEG)+IJ) 113 CONTINUE 114 CONTINUE IF (ANS.NE.'L') GO TO 120 C C USER HAS CHOZEN A LOG DIGITAL FILTER USE CHANNEL NUMBER C AND ADD A SMALL NUMBER SO WE DON'T RUN INTO LOG OF NEGATIVE C NUMBER OR LOG OF ZERO PROBLEMS C X=I X0 =ALOG(X-1) FX0=ALOG(TFX0) X1 =ALOG(X) FX1=ALOG(TFX1) X2 =ALOG(X+1) FX2=ALOG(TFX2) GO TO 130 C C USER HAS CHOZEN NORMAL (LINEAR) FILTER C 120 FX0=TFX0 FX1=TFX1 FX2=TFX2 FX3=TFX3 FX4=TFX4 130 CONTINUE C C FOR LINEAR FILTERS DEFINE THE NORMAL ZERO-AREA FACTOR C FACT=(ANEG-APOS2)/APOS !USED ONLY FOR NON-LOG FILTERS IF(ANS.EQ.'T') AX(I)=(FX3 -FX2 + 2*FACT*FX1 - FX0 + FX4) IF(ANS.EQ.'L') AX(I)=(-FX2*(X1-X0)+FX1*(X2-X0)-FX0*(X2-X1)) C/((X2-X1)*(X2-X1)*(X1-X0)) C C DONE WITH THIS DATA POINT NOW CONTINUE THE LOOP C 110 CONTINUE C C C FINISHED THE FILTERING NOW C STORE ORIGINAL SPECTRUM IN MEMORY #2 'CCY' C STORE DIGITAL FILTERED SPECTRUM IN MEMORY#1 'BY' C COPY DIGITAL FILTERED SPECTRUM INTO DISPLAY MEMORY 'AY' C COPY X AXIS CALIBRATION INTO 'AX' C DO 201 I=1,NNPT CCY(I)=BY(I) !MOVE ORGINAL SPECTRUMN TO CCY BY(I)=AX(I) !MOVE DIG FIL TO BY AY(I)=AX(I) !ALSO MOVE IT TO AY 201 AX(I)=OFFSET + EVCH*I !RECALCULATE X AXIS DO 200 I=NNPT+1,NNPT+10 200 CCY(I)=BY(I) !COPY CALIBRATION CONSTANTS TOO C C LET EVERYBODY KNOW WHAT'S HAPPENING TO THE DATA C C READ(5,13) RES 13 FORMAT(1A2) WRITE(7,1000) 1000 FORMAT(' DIG. FIL. IN MEM#1,ORIG. DATA IN MEM#2: TO CONTINUE') RETURN END C----------------------------------------------------------------------- SUBROUTINE RPDL(XS,NPT) C----------------------------------------------------------------------- C SUBROUTINE TO READ AN ASCII DATA FILE IN EMMPDL DATA FORMAT C INTO AN ARRAY XS(4110) C C FIRST TWO LINES OF OF FILE CONTAIN TEXT ARE FOR IDENTIFICATION C AND CAN BE READ USING A 80A1 FORMAT C THE NEXT LINE CONTAINS THE EXPERIMENTAL DATA STORED AS 10 ASCII C NUMBERS (FLOATING POINT FORMAT) EACH DATA POINT SEPERATED BY A C COMMA AS A DELIMINATOR, C THE FIRST NUMBER OF THE FIRST LINE IS THE TOTAL NUMBER OF DATA C POINTS IN THE SPECTRAL FILE, AND IS LIMITED TO A MAXIMUM VALUE C OF 4095, IF IT'S VALUE IS "NPT" THEN THE MAXIMUM NUMBER OF C LINES OF DATA IN THIS FILE IS GIVEN BY NLINES="NPT/10" +1 C FOLLOWING NLINES OF DATA THERE IS A FINAL LINE OF 10 VALUES C WHICH ARE THE CALIBRATION CONSTANTS OF THE SPECTRUM THEY C ARE DEFINED AS FOLLOWS: C C XS(1) = NPT = NUMBER OF POINTS (CHANNELS) IN THE SPECTRUM C XS(2)-XS(NPT) =DATA C C XS(NPT+1) TO XS(NPT+10) =CALIBRATION CONSTANTS OF SPECTRAL DATA C NPT+1= OFFSET ENERGY OF FIRST CHANNEL C NPT+2= EV/CHANNEL C NPT+3= ACCELERATING VOLTAGE (KV) C NPT+4= INCIDENT BEAM DIVERGENCE (MR) FOR EELS DATA C = X-AXIS (ROD AXIS) TILT (DEGREES) FOR XEDS DATA C NPT+5= SCATTERING ANGLE (MR) C = Y-AXIS (CUP AXIS) TILT (DEGREES) FOR XEDS DATA C NPT+6= LIVE TIME (DWELL*NUMBER OF CHANNELS) FOR EELS OR c SIMPLE LIVE TIME FOR XEDS DATA C NPT+7= DEAD TIME C NPT+8= INCIDENT BEAM CURRENT IN nA C NPT+9= PROBE DIAMETER IN nm C NPT+10= SPECIMEN THICKNESS ALONG BEAM DIRECTION IN nm. C C C C----------------------------------------------------------------------------- C DIMENSION INAME(8) DIMENSION XS(4110),IXL(80) BYTE DAY(9),TIM(8) C C C C READ IN DATA FILE NAME C 204 WRITE(7,1) 1 FORMAT(' ENTER EMMPDL DATA FILE NAME (DEV:NAME.EXT):',/) READ(5,205) (INAME(J),J=1,7) 205 FORMAT(7A2) INAME(8)=0 OPEN (UNIT=11,NAME=INAME,TYPE='OLD',ACCESS='SEQUENTIAL') C C C WRITE IT OUT C C CALL DATE(DAY) CALL TIME(TIM) WRITE(8,206) DAY,TIM,(INAME(J),J=1,7) 206 FORMAT(1X,9A1,' at ',8A1,/,' Input Data File = ',7A2) C C READ IN TWO LINES OF TEXT C WRITE(7,21) WRITE(8,21) 21 FORMAT(' SPECTRUM IDENTIFICATION',/,23('-')) C C FILL XS(1021-4110) WITH ZERO'S C READ(11,2) IXL 2 FORMAT(80A1) WRITE(7,3) IXL WRITE(8,3) IXL 3 FORMAT(1X,80A1) READ(11,2) IXL WRITE(7,3) IXL WRITE(8,3) IXL C C READ FIRST DATA LINE AND DETERMINE THE NUMBER OF CHANNELS C READ (11,1000) (XS(I), I=1,10) 1000 FORMAT (10F12.0) NPT=XS(1) WRITE (7,1001) NPT 1001 FORMAT (' This spectral file contains ',I4,' channels of data') C C CALCULATE THE REMAINING NUMBER OF LINES (YOU HAVE ALREADY READ ONE!) C IF((XS(1)/10. - NPT/10 ).GT.0) NLINES= NPT/10 +1 IF((XS(1)/10. - NPT/10 ).EQ.0) NLINES= NPT/10 C C NOW LOOP THROUGH THE DATA READING EACH LINE C DO 4 I=2,NLINES IF=10*I-9 IL=10*I IF (IL.GE.NPT) IL=NPT C READ(11,1000,END=4,ERR=4) (XS(J), J=IF,IL) 4 CONTINUE C C READ IN CALIBRATION CONSTANTS C IF THEY EXIST C READ(11,1000,ERR=52,END =52) (XS(J), J =NPT+1,NPT+10) 52 CONTINUE IF((XS(NPT+1).EQ.0).AND.(XS(NPT+2).EQ.0)) WRITE(7,53) IF((XS(NPT+1).EQ.0).AND.(XS(NPT+2).EQ.0)) WRITE(8,53) 53 FORMAT(' UNCALIBRATED SPECTRUM:') WRITE(7,7) (XS(J), J=NPT+1,NPT+5) 7 FORMAT(//, C ' OFFSET = ',1F9.3,' eV',/, C ' EV/CHANNEL = ',1F9.3,/, C ' ACCELERATING VOLTAGE = ',1F9.1,' kV',/, C ' INCID. BEAM DIVERGENCE = ',1F9.3,' mR',/, C ' SCATTERING ANGLE = ',1F9.3,' mR',/, C ' LIVE TIME = ',1F9.3,' SEC',/, C ' DEAD TIME = ',1F9.3,' SEC',/, C ' BEAM CURRENT = ',1F9.3,' nA',/, C ' PROBE DIAMETER = ',1F9.3,' nm',/, C ' SPECIMEN THICKNESS = ',1F9.3,' nm',/) 6 WRITE(7,600) 600 FORMAT (' HIT TO CONTINUE') READ(5,2) X CALL CLOSE(11) CALL ERASE !ERASE SCREEN RETURN END C C C C C----------------------------------------------------------------------- SUBROUTINE WPDL(XS,NPT) C----------------------------------------------------------------------- C SUBROUTINE TO WRITE AN ASCII DATA FILE IN EMMPDL DATA FORMAT C FROM AN ARRAY XS(NPT) C C FIRST TWO LINES OF OF FILE CONTAIN TEXT ARE FOR IDENTIFICATION C AND CAN BE READ USING A 80A1 FORMAT C THE NEXT LINE CONTAINS THE EXPERIMENTAL DATA STORED AS 10 ASCII C NUMBERS (FLOATING POINT FORMAT) EACH DATA POINT SEPERATED BY A C COMMA AS A DELIMINATOR, C THE FIRST NUMBER OF THE FIRST LINE IS THE TOTAL NUMBER OF DATA C POINTS IN THE SPECTRAL FILE, AND IS LIMITED TO A MAXIMUM VALUE C OF 4095, IF IT'S VALUE IS "NPT" THEN THE MAXIMUM NUMBER OF C LINES OF DATA IN THIS FILE IS GIVEN BY NLINES="NPT/10" +1 C FOLLOWING NLINES OF DATA THERE IS A FINAL LINE OF 10 VALUES C WHICH ARE THE CALIBRATION CONSTANTS OF THE SPECTRUM THEY C ARE DEFINED AS FOLLOWS: C C XS(1) = NPT = NUMBER OF POINTS (CHANNELS) IN THE SPECTRUM C XS(2)-XS(NPT) =DATA C C XS(NPT+1) TO XS(NPT+10) =CALIBRATION CONSTANTS OF SPECTRAL DATA C NPT+1= OFFSET ENERGY OF FIRST CHANNEL C NPT+2= EV/CHANNEL C NPT+3= ACCELERATING VOLTAGE (KV) C NPT+4= INCIDENT BEAM DIVERGENCE (MR) FOR EELS DATA C = X-AXIS (ROD AXIS) TILT (DEGREES) FOR XEDS DATA C NPT+5= SCATTERING ANGLE (MR) C = Y-AXIS (CUP AXIS) TILT (DEGREES) FOR XEDS DATA C NPT+6= LIVE TIME (DWELL*NUMBER OF CHANNELS) FOR EELS OR c SIMPLE LIVE TIME FOR XEDS DATA C NPT+7= DEAD TIME C NPT+8= INCIDENT BEAM CURRENT IN nA C NPT+9= PROBE DIAMETER IN nm C NPT+10= SPECIMEN THICKNESS ALONG BEAM DIRECTION IN nm. C C C C----------------------------------------------------------------------------- DIMENSION INAME(8) C C C C DIMENSION XS(4110),IXL(80) BYTE DAY(9),TIM(8) C C CALL ERASE !ERASE SCREEN C C READ IN DATA FILE NAME C 204 WRITE(7,1) 1 FORMAT(' ENTER EMMPDL DATA FILE NAME (DEV:NAME.EXT):',/) READ(5,205) (INAME(J),J=1,7) 205 FORMAT(7A2) INAME(8)=0 OPEN (UNIT=11,NAME=INAME,TYPE='NEW',ACCESS='SEQUENTIAL') C C C WRITE IT OUT C C CALL DATE(DAY) CALL TIME(TIM) WRITE(8,206) DAY,TIM,(INAME(J),J=1,7) 206 FORMAT(1X,9A1,' at ',8A1,/,' Output Data File = ',7A2) C C READ IN TWO LINES OF TEXT C WRITE(7,21) 21 FORMAT(' SPECTRUM IDENTIFICATION',/,23('-'),/,' Enter 1st line') READ(5,2) IXL 2 FORMAT(80A1) WRITE(11,3) IXL 3 FORMAT(1X,80A1) WRITE (7,31) 31 FORMAT(' Enter 2nd line') READ(5,2) IXL WRITE(11,3) IXL C C READ FIRST DATA LINE AND DETERMINE THE NUMBER OF CHANNELS C XS(1)=NPT WRITE (7,101) NPT 101 FORMAT (' This spectral file contains ',I4,' channels of data') C C CALCULATE THE REMAINING NUMBER OF LINES (YOU HAVE ALREADY READ ONE!) C C C NOW LOOP THROUGH THE DATA WRITING EACH LINE C C WRITE (11,1000) (XS(J), J=1,NPT) 1000 FORMAT(1X,9(1F9.0,','),1F9.0) C C WRITE CALIBRATION CONSTANTS C 4 CONTINUE WRITE (11,1001) (XS(J), J =NPT+1,NPT+10) 1001 FORMAT(1X,9(1F9.3,','),1F9.3) 52 CONTINUE CALL CLOSE(11) RETURN END