A mining research contract report 

FEBRUARY 1984 

ADAPTIVE-NOISE-CAilCELLATION 

TECHNIOUES 

FOR TH ROUGH-TH Em?1RTH 

ELECTROMAGNETICS 

VOLUME Ill - 4.3 

Green Mountain Radio Research Company 

BUREAU OF MINES 

UNITED STATES DEPARTMENT bF THE INTERIOR 
DISCLAIMER NOTICE 

The views and conclusions contained in this document are those of the author 

and should not be interpreted as necessarily representing the official pol icies or recommendations of the Interior Department's Bureau of Mines or of 

the U.S. Government. 
1 

( 2. 3. Rcripient's Accession No. 1 

SHEET I I I 

4. Tirie and Subtitle IS. Reporr Dare 1 

I Green Mounta l n Rad lo Research Company 

50 Vermont Avenue 

Adaptive-Noise-Cancellation Techniques for Through-tho-Earth Electromagnetlcs, Volume I1 I 

7. Author(s) Frederick H. Raab, Ph.D. 

9. Performing Organization Name and Address 

Fort Ethan Allen 

Winooski, Vermont 05404 

February 1984 

6. Performing Organization Code 

3U750 

8. Performing Organization Repr. 

No- TR84- 1 

10. Project/Task/Work Unit No. 

11. Contracr/Grant No. 1 

United States Department of the Interior, Bureau of Mines 

Washington, DC 20241 

1 

I F l NAL I 

12 Sponsoring Agency Name and Address 

1~~1~ 1983 - February 1984 1 

14. Sponsoring Agency Code 

13. Type of Report & Period 

Covered 

I 

15. Supplementary Notes 

Thls report is Volume Ill of a three-volume set and covers work done under Phase-I1 I of the contract. 

16. Abstracts 

An electromagnetic (EM) system for the location of trapped mine workers at coal-mine (300 m) depths 

has been developed. The use of such a system in deep mines (I-km depth) requlres transmission at extremel y low frequencies (ELF) and a reduction In the effective noise level beyond that whlch can be obtained by acceptable integration times. Automatic-noise-cancellation (ANC) techniques reduce the effective noise level by combining several received signals that contain correlated nolse. Such signals can 

be obtained elther from a remote magnetic antenna or a local electric antenna. 

During Phase Ill of this program, real ELF nolse data were acquired and used to test the ANC algorithms developed during Phase 11. The ELF noise data were obtained from analog recordings made by 

Develco under another Bureau-of-Mines contract. Digltization was accomplished by an A/D converter that 

was added to a GMRR computer. Phase Ill tasks included preparation of software for dlgltization, cow 

puter-computer comnunicatlon, data processing, and statistical analysis. 

Useable data were obtained from eleven half-hour tapes containing nolse from one axis of the primary antenna and three axes of the remote antenna. The overall average improvement In SNR was 4.1 dB. 

However, the average improvements for the lndlviduai X, Y, and Z axes were 4.2, 11.9, and -1.9 dB, which 

roughly follow the observed noise power. Nonlinear processing was not effective with this noise. 

17. Key Words and Document Analysis. 170. Descriptors 

Electromagnetic direction finding 

Trapped-miner location 

Trapped-miner comnunications 

Through-the-earth comnunicat ions 

Near-field locatlon systems 

Adaptive noise cancellation 

Extremely low-frequency (ELF) noise 

Adaptive clipper 

Nonlinear processing 

17b. Identificrs/Open-Ended Terms 

17c. COSATl Fie ld/Group 

21. No. of Pages 

97 

22. Price 

18. Distribution Statement 

Approved for pub1 lc release, distribution uni irnlted. 

FORM NT15-35( 10-70) US COMM-0 C 40329-F 

19. Security Class (This 

Report ) 

mi s 

u 

l'ase) 

UNCLASSIFIEL) 
FORE WARD 

This report was prepared by Green Mountain Radio Research Company, 

Winooski, Vermont under contract 50318070. This contract was initiated under the Coal Mine Health and Safety Program. It was administered under the 

technical direction of the U.S. Bureau of Mines with Steve Shope acting as 

Technical Project Officer. Karen Hatters was the contract administrator for 

the Bureau of Mines. This report is a summary of the work recently completed 

as a part of this contract during the period July 1983 - February 1984. This 

report was submitted by the author in February 1984. 

The manuscript and artwork were prepared by Patricia L. Scott. 
CONTENTS 

Chapter 1. Introduction ............................................ 1.1 Background .............................................. 1.2 ANC Concept ............................................. 1.3 Phase-I Results ......................................... 1.4 Phase-I1 Results ........................................ 1.5 Phase-I11 Objectives and Organization ................... 1.6 References .............................................. Chapter 2. Digitization ............................................ 2.1 Approach ................................................ 2.2 AID-Interface Subroutines ............................... 2.3 SemiDisk-Interface Subroutines .......................... 2.4 Operational Programs .................................... 2.5 Operational Procedures .................................. 2.6 References .............................................. Chapter 3. Acquisition ............................................. 3.1 Noise Data .............................................. 3.2 Storage and Transfer .................................... Chapter 4. Processing .............................................. 4.1 Single-Channel Simulator ................................ 4.2 Quantization ............................................ 4.3 Digitized-Noise Processor ............................... 4.4 Fixed Cancellation ...................................... 4.5 Reference ............................................... Chapter 5. Results ................................................. 5.1 Useable Data ............................................ 5.2 Convergence ............................................. 5.3 Processing Gain ......................................... 5.4 Weighting Vector ........................................ 5.5 Noise Correlation ....................................... 5.6 Analysis ................................................ Chapter 6. Conclusions and Recommendations ......................... 6.1 Accomplishments ......................................... 6.2 Results and Conclusions ................................. 6.3 Recommendations ......................................... Appendix A. Digitization Programs .................................. Appendix B. Communication Program .................................. Appendix C. Processing Program ..................................... Appendix D. Gain Statistics ........................................ 

PAGE 

5 . 12 

5- 6 

7- 9 

9 

9 . 11 

11 

11 . 12 

13 . 22 

13 . 16 

16 . 18 

18 

18 . 21 

2 1 

2 2 

23 . 28 

23 . 26 

27 . 28 

29 . 36 

2 9 

29 . 34 

34 

34 . 36 

34 

37 . 50 

3 7 

3 7 

37 . 41 

41 . 47 

48 

48 . 50 

51 . 53 

5 1 

51 . 52 

52 . 53 

54 . 66 

67 . 72 

73 . 90 

91 . 97 
ILLUSTRATIONS 

PAGE 

Elements of electromagnetic location system ..................... Example adaptive noise-reduction technique ...................... Equipment ....................................................... Filter schematic ................................................ Data packing .................................................... Timing .......................................................... Antenna sites ................................................... Effect of channel upon NLP improvement .......................... Effect of channel upon ANC improvement .......................... Effect of quantization upon NLP improvement. X axis. moderatelevel noise ................................................... Effect of quantization upon ANC improvement. X axis. moderatelevel noise ................................................... Effect of quantization upon NLP improvement. X axis. highlevel noise ................................................... Effect of quantization upon ANC improvement. X axis. highlevel noise ................................................... Effect of cancellation methods. moderate-level noise. X axis .... Effect of cancellation methods. moderate-level noise. Y axis .... Effect of cancellation methods. moderate-level noise. Z axis .... Convergence plot for data set #8-1A (X axis) .................... Convergence plot for data set #12-2A (Y axis) ................... Convergence plot for data set #6-2A (2 axis) .................... Convergence plot for data set #8-1B (X axis) .................... Convergence plot for "disturbed" data set #1-3 (Y axis) ......... Average ANC gains. all axes ..................................... Average ANC gains. X axis ....................................... Average ANC gains. Y axis ....................................... Average ANC gains. Z axis ....................................... Comparison of average total ANC/NLP gain by axis ................ Histogram of average processing gain. all axes .................. Histogram of average processing gain. X axis .................... Histogram of average processing gain. Y axis .................... Histogram of average processing gain. Z axis .................... Comparison of average total ANCINLP and FNC/NLP gains ........... Average weighting kr ............................................. Average noise covariance R .................................... Y Y 

Average noise cross-correlation R ............................. XY 

TABLES 

2-1 Data-identification block ....................................... 20 

3-1 List of tapes. files. and conditions ............................ 24 

3-2 General recordingldigitizing parameters ......................... 26 

4-1 Simulation results ........................................... 31 

5-1 Processing-gain statistics ................................... 42 

5-2 Average covariance. cross correlation. and weighting ............ 47 
5 

CHAPTER 1. INTRODUCTION 

1.1 BACKGROUND 

Electromagnetic (EM) techniques and apparatus for the subsurface location of trapped mine workers have been developed during the past decade. The 

"first-generation" system [l - 4) transmits signals in the 600 - 3000 Hz 

range to maximize the signal-to-noise ratio. Physical movement of the receiver and detection by human ear is required, and the accuracy can be degraded by irregularities in the conducting ground. 

The elements of the "second-generation" EM location system [5] are shown 

in Figure 1-1. The subsurface beacon transmitter is similar to that of the 

first-generation system, but uses a lower transmitting frequency. Uplink 

signals are received by three three-axis magnetic-field antennas (called 

"search coils"). The amplitudes of the nine field components are forwarded 

to a computer in the electronics van. The computer processes the signals and 

determines the position of the subsurface transmitter. The requirements for 

physical movement and detection by human ear are eliminated by this system. 

The deep-mine (1-km depth) application of electromagnetic (EM) trappedminer location techniques [5] and the desire to reduce the effects of the 

conducting ground [6 - 9) mandate the use of lower uplink frequencies. At 

such frequencies (1 - 30 Hz), the ambient-noise levels [5,7,10] are considerably greater than those in the 600-to-3000-Hz band used by the first-generation EM location system. However, intrinsic-safety considerations prevent 

significant increases in the output of the beacon transmitter, and the link 

analysis [5] clearly shows that simple signal integration alone does not produce adequate signal-to-ratio in an acceptable length of time. 

The primary receiving antenna is a three-axis magnetic-field sensor, and 

receives both the uplink signal and noise. Since the uplink signal is a quasi-static magnetic field, it decays very rapidy and has an insignificant associated electric field. In contrast, propagating atmospheric-noise fields 

decay very slowly and include both electric- and magnetic-field components. 

An electric-field antenna or a remote magnetic-field sensor therefore receives the atmospheric noise, but not the beacon signal. 

Since the noises received by the primary and reference antennas are in 

general highly correlated, significant improvement in the signal -to-noi se ratio (SNR) can be achieved by combining their outputs. Since the correlation 

characteristics of the antenna outputs change with time and are therefore not 

known a priori, an adaptive-noise-cancel lation (ANC) algorithm is required to 

combine the outputs to obtain a signal with minimum noise. 

A preliminary investigation of ANC techniques for through-the-earth 

(TTE) EM systems was conducted during Phase I of this program [ll]. A thorough simulation of ELF noise, nonlinear processing, and ANC algorithms was 

conducted under Phase I1 of this program [12]. This report presents the results of Phase 111, during which real ELF noise data were acquired and used 

to evaluate the performance of the nonlinear processing and ANC algorithms 

that were developed during Phase 11. 
H-F I ELD 

NOISE-REFERENCE 

\ /P ELEcTRoN1cs VAN I PR I MARY 

ELF 

MAGNETIC-DIPOLE 

FIELD 

Figure 1-1. 

SUBSURFACE 

BEACON 

LOOP TRANSMITTER 

ANTENNA 

Elements of electromagnetic location system. 
1.2 ANC CONCEPT 

A simplified example of an automatic-noise-cancellation technique is 

shown in Figure 1-2. The output of the primary receiving antenna is represented by x and contains the desired signal corrupted by a significant amount 

of noise. The signal-to-noise ratio (SNR) in this example is about -10 dB. 

Waveform y represents the output of a noise-reference antenna; in the 

present application, this antenna could be either a remote magnetic-field 

sensor or a local electric-field sensor. This received signal contains little or nothing of the desired signal, but contains noise that is highly correlated with that in the output of the primary antenna (e.g., a - 1 and b 

<< 1). 

The output signal z is a linear combination of the primary and noisereference signal s. If the constant c is set properly, much of the noise is 

cancelled, producing the +lo dB SNR shown in Figure 1-2b. 

In general, the noise statistics (hence constant c) are not known a 

priori and are time varying. It is therefore necessary to use an adaptive 

technique to set the value of c for maximum noise cancellation. A variety of 

such techniques have been developed for other applications [I]. 

The transmitters in throuah-the-earth electromaanetic svstems use elec9 loop altpnnas aid cm therefore ha re&w@d. ds m 

8 se of extremely low excitation frequencies to re 

fects of the conducting ground also guarantees that the electric field produced by the subsurface transmitter will be negligible within its region of 

operation. Consequently, the field produced is essentially a pure quasistatic magnetic dipole field having a free-space structure with an inversedistance-cubed amplitude variation. 

Atmospheric noise at extremely 1 ow frequencies is produced primari ly by 

thunderstorms. Nearby thunderstorms produce distinct noise impulses, while 

the aggregate of distant thunderstorms produces a Gaussian background noise. 

Since much of the received noise power is contained in the impulses, it is 

necessary to include processing that edits (blanks) or cancels individual impul ses. 

Because the attenuation rates for ELF signals are very low (e.g., 0.2 

dB/Mm), most of the received noise power (especially after editing) is due to 

distant sources. Since the wavelengths of the signals of interest are 6 Mm 

(3600 mi) or more, noise waveforms tend to be similar at points separated by 

significant distance (15 km). The primary causes of differences in the noise 

waveforms are local geological features (e.g., changes in ground conductivity). Nonetheless, noise waveforms with correlations from 0.9 to 0.9999 (10 

to 40 dB) can be found at sites separated by as much as 15 krn. (This correlation is the basis of remote-reference magnetotell urics). Since the amp1 i - 

tude of the desired signal decreases by 60 dB as distance is increased from 1 

km to 10 km, it is apparent that a remote magnetic-field sensor can provide 

an essentially signal-free source of correlated noise. 
INPUTS OUTPUTS 

(a) Block diagram 

X=6+n 

1 

- - * z=x+cy 

GAIN > 

y = an + bn c 

1 2 

1 

(b) Waveforms 

CONTROL 

CIRCUIT 

Figure 1-2. Example adaptive noise-reduction technique. 

i 

A ' 

t 

NOISE 

POWER 

DETECTOR 
A propagating electromagnetic signal contains both electric and magnetic-field components, which are related by the free-space impedance (377 Q). 

Propagation over conducting ground produces the well-known phenomenon called 

"E-field tilt,'' which causes the electric-field to tilt in the direction of 

propagation. This implies that the horizontal electric-field components are 

strongly correlated with the opposite horizontal magnetic-field components. 

Since the beacon transmitter produces little, if any, electric field, it is 

apparent that a local electric-field antenna can also be used as a source of 

correlated noise. The magnetotelluric technique for geophysical prospecting 

is based upon this phenomenon. A noise-cancellation technique based upon 

electric- and magnetic-field correlation has been proposed for use with trolley phones [13]. 

1.3 PHASE-I RESULTS 

The objective of Phase I was to conduct a preliminary investigation and 

feasibility assessment of ANC techniques for TTE EM systems [I]. The two 

principal areas of investigation were: 

ELF-noise characteristics, and 

ANC algorithms. 

Phase I work also included analysis of correlation estimators. 

The TTE EM location system involves reception of all three magneticfield components and (if ANC is used) possibly all three electric-field components. The incident field includes the vertical electric and horizontal 

magnetic components. The Field-Lewinstein-Modestino model represents singlechannel ELF noise as the sum of Gaussian and power-Ray1 ei gh random-vari abl es. 

Direction-of-arrival characteristics are accommodated by using multiple noise 

sources positioned to correspond to the major thunderstorm areas of the world. 

Horizontal -electric and vertical -magnetic fie1 d components are generated 

by a "geology filter" that accounts for the electromagnetic characteristics 

of the ground. 

The limited amount of published empirical data for the frequencies of 

interest show correlations of both the electric- and magnetic-field components ranging from 0.9 to 0.9999. These correlation levels suggest the possibility of noise reductions from 10 to 40 dB through ANC. 

Candidate ANC algorithm including direct-matrix-inversion (DMI), gradient (least-mean-square - LMS), and recursive (Kalman-fi 1 ter) techniques are 

identified. The DM1 and LMS techniques have the most potential in the TTE EM 

application. 

1.4 PHASE-I1 RESULTS 

The objective of Phase I1 was thorough comparison and evaluation by 

simulation of ANC techniques for the TTE EM application. The specific areas 
of investigation included: 

The Field-Lewinstein-Modestino (FLM) ELF-noise model, 

Di rection-of-arri val (DOA) characteristics, 

The geology fi 1 ter, 

Q Effects of reception bandwidth, 

Nonlinear-processing techniques, 

e ANC algorithms, and 

The combination of nonlinear processing and ANC. 

Software was developed for random-noise generation, including subroutines for uniform, Gaussian, and power-Rayleigh random-variables. The 

Gaussian and power-Rayleigh random variables are obtained by non-1 inear 

transformation of the uniform random variable. Proper operation is verified 

by moment and histogram tests. 

Software for the multicomponent ELF-noise simulation model was implemented. In this model, the output from individual noise sources is obtained by 

summing Gaussian and power-Rayleigh random variables as specified by the FLM 

model. The vertical-electric-field component is obtained by summing the outputs from several such sources. The horizontal-magnetic-field components are 

obtained by transforming the outputs of individual sources into X and Y components and then summing the results. The unit-pulse response of a homogeneous ground is derived and used to implement a geology filter. Proper operation is verified by moment, amplitude-histogram, and DOA-histogram tests. 

The FLM model characterizes ELF noise by power, impulsivity, and spikiness. However, values of the impulsivity and spikiness parameters are published for only a limited number of bandwidths. The effect of bandwidth upon these parameters was therefore analyzed, and a two-moment method for conversion was devel oped. 

A considerable portion of the power in ELF noise is contained in discrete impulses. Nonlinear processing can therefore be used to reduce the effective received noise power. Nonlinear processors based upon editing, clipping, adaptive editing, and adaptive clipping were compared by simulation. 

The adaptive clipper was found to be most effective. 

Narrowband DM1 and LMS ANC algorithms were implemented and compared by 

simulation with Gaussian noise. In this application, the DM1 algorithm must 

both estimate the noi se-covariance matrix and the weighting matrix. Since 

errors in one estimate produce errors in the other, it is effective to apply 

the DM1 algorithm recursively. (Two iterations are sufficient). The simulations show a gain of 0.1 produces the greatest noise reduction by the LMS 

algorithm. The simulations also show the DM1 algorithm to produce significantly greater noise reduction than does the LMS algorithm. 

The DM1 algorithm was modified to operate on complex inputs that represent the narrowband outputs of the nonlinear processor (NLP). The DM1 ANC 

and NLP algorithms were then combined; proper operation of both subsystems 

was verified by simulation. The combined algorithm was then tested under a 

variety of conditions. The benefits obtained with one algorithm alone are 

not additive, but the inclusion of both algorithms is nonetheless desirable. 

Under typical atmospheric conditions, the total noise reduction ranges from 5 
to 30 dB, dependi ng upon the 1 ocal -noi se 1 eve1 . 

1.5 PHASE-I11 OBJECTIVES AND ORGANIZATION 

The objective of Phase I11 is to evaluate the nonlinear-processing and 

adaptive-noise-cancellation algorithms with actual ELF noise. The specific 

areas of work include: 

e Development of software for data digitization, 

Development of software for data transfer, 

Data acquisition, 

Modification of NLP/ANC software algorithms, and 

Data processing and statistical evaluation of the results. 

Chapter 2 describes the data-digitization hardware and software. Digitization is accomplished by an I/O Technology 12-bit A/D converter board 

added to one of GMRR's North Star Horizon computers. The digitization software uses the computer's real-time clock to digitize four analog-input channels at a 100-Hz rate. A semiconductor disk emulator is used for temporary 

data storage to avoid problems with the track-to-track access time of the 

computer's f 1 oppy-di sk drives. 

The ELF-noise data and the data-acquisition procedure are described in 

Chapter 3. Data acquisition was accomplished at Develco's facility in Sunnyvale, CA. This procedure eliminated the possibility of degradation of the 

analog tapes by rerecording or shipment. A total of fourteen half-hour tapes 

were digitized, resulting in some 16 megabytes of data. Digitized data were 

retained by GMRR for subsequent processing and transferred to the Develco 

HP-1000 computer for their possible future use. 

Chapter 4 discusses the data-processing software, which is based upon 

the simulation software developed in Phase 11. Modifications to the simulation software a1 low operation with a single-axis primary antenna, reading 

noise data from a disk file, and extraction of a calibration constant. A 

simulation of fixed (nonadaptive) cancel lation is a1 so added. 

The results of processing the data are presented in Chapter 5. Data 

from each analog tape were partitioned into six segments of 270-s duration. 

Each 270-s data set was used for one simulation run and processed with various clipping thresholds. The best results were obtained with nonlinear processing disabled. An overall average ANC gain of 4.1 dB was obtained for data 

from the eleven "undisturbed" tapes. The average ANC gains for the individual X, Y, and Z axes were 4.3, 11.9, and -1.9 dB, respectively. 

1.6 REFERENCES 

1. H. K. Sacks, "Trapped-miner location and communication systems," hderground Mine Communications, Part 4, Section-to-Place Communications, 

U.S. Bureau of Mines Information Circular 8745, 1977. 
2. R. G. Olsen and A. J. Farstad, "Electromagnetic direction finding experiments for location of trapped miners," IEEE Transactions on Geoscience 

EZectronics, vol. GE-11, no. 4, pp. 178 - 185, October 1973. 

3. J. Durkin, "Assessment of present techniques for the location of trapped 

miners , " Proceedings of the Fifth WVU Conference on Coal Mine nectrotechnology, Morgantown, WV, pp. 1-1 - 1-16, July 30 - August 1, 1980. 

4. R. Gabillard, J. P. Dubus, and F. Cherperell, "Electromagnetic survey 

method applicable to underground quarries," Proceedings of the Thru-theEarth EZectromagnetic Workshop, Golden, CO, (NTIS PB 231,154), pp. 121 - 

129, August 15 - 17, 1973. 

5. L. H. Rorden, T. C. Moore, E. C. Fraser, and L. R. Bulduc, "Phase I development of a system concept for location of trapped miners in deep 

mines by an electromagnetic method," Phase I Report, contract J0199009, 

Develco, Inc., Sunnyvale, CAY October 1979. 

6. T. W. H. Caffey, "Locating a buried magnetic dipole," Proceedings of the 

Zectromagnetic Guided Waves Workshop, Boulder, COY March 28 - 30, 1978. 

7. F. H. Raab and P. K. Hansen, "Electromagnetic retransmission system for 

locating trapped mine workers," Final Report, contract H0188071, (PB81- 

215329) Polhemus Navigation Sciences, Inc., Essex Junction, VT, February 

1980. 

8. F. H. Raab and P. K. Hansen, "Electromagnetic technique for subsurface 

position determination," Proceedings of the IEEE Industry Applications 

Society Annual Meeting, Philadelphia, PA, pp. 111 - 118, October 5 - 9, 

1981. 

9. F. H. Raab, "Quasi-static magnetic-field technique for determining position and orientation," IEEE Transactions on Geoscience and Remote Sensing, vol. GE-19, no. 4, pp. 235 - 243, October 1981. 

10. E. L. Maxwell, "Atmospheric noise from 20 Hz to 30 kHz," Radio Science, 

vol. 2, (new series), no. 6, pp. 637 - 644, June 1967. 

11. F. H. Raab, "Adaptive-noise-cancellation techniques for through-theearth electromagnetics, Volume I," Final Report GMRR-TR82-1, Green 

Mountain Radio Research Company, Burlington, VT, January 1982. 

12. F. H. Raab, "Adaptive-noise-cancellation techniques for through-theearth electromagnetics, Volume 11," Final Report GMRR-TR83-1, Green 

Mountain Radio Research Company, Winooski, VT, September 1983. 

13. R. L. Lagace, D. A. Curtis, J. D. Foulkes, and J. L. Rothery, "Technical 

service for mine communications research; Transmit antennas for portable 

VLF to MF wireless mine communications," Contract No. H0346045 - Task 

Order No. 1, Arthur D. Little, Inc., Cambridge, MA, May 1974 - May 1977. 
CHAPTER 2. DIGITIZATION 

ELF noise data suitable for use in evaluating the performance of ANC algorithms were collected by Develco, Inc., as part of their work on the deepmines EM location system under Bureau-of-Mines contract 50199009. The desired 

noise data reside in a set of analog recordings made by a TEAC R-61 portable 

cassette data recorder. Four channels were recorded simultaneously, and the 

input bandwidth of each channel is about 30 Hz. 

To evaluate performance of the various ANC algorithms on GMRR's North 

Star computers, it was therefore necessary to digitize the analog recordings 

and to store the resultant data on North Star floppy disks. This chapter 

describes the software and procedures developed to accomplish the data digitization and transfer. 

2.1 APPROACH 

The approach (Figure 2-1) selected for acquisition of the noise data required equiping a GMRR computer with an A/D converter, taking it to Develco, 

and digitizing the data there. This approach allowed processing of the original analog recordings and also allowed transfer of the digitized data to the 

Develco computer. In contrast to other possible approaches, it avoided 

equipment-rental fees, noise added in analog re-recording, and the possibility of damaging the original recordings during shipment. 

A North Star Horizon computer is used to control data digitization and 

storage. This computer has a 280 microprocessor, S-100 bus, 56 kilobytes of 

RAM, and two 338-ki 1 obyte-capaci ty fl oppy-di sk dri ves. This computer i s a1 so 

equipped with a 512-kilobyte SemiDiskm (semiconductor disk emulator), and 

uses the CP/M-80@ operating system. 

Fi 1 ter 

The four-channel 1 ow-pass f i 1 ter assembly suppresses broadband noi se (0 

- 1000 Hz) in the recorder output. This noise cannot be suppressed in software because of the 100-Hz sampling frequency. 

The filter assembly was built by Develco engineer David Schleicher, and 

each channel employs a two-pole active circuit (Figure 2-2) similar to that 

used in Develco's deep-mine EM receiver. Each channel has unity gain at low 

frequency, both poles at 29-Hz, and a theoretical response of 

Measurements of the filter response show 22 dB attenuation at 100 Hz and not 

more than 2" - 3" phase shift in the 0 - 30-Hz passband. 
NORTH STAR HORIZON 

~igufe 2-1. Equipment. 

Digitization 

Data digitization is accomplished by an 1/0 Technology S-100 A/D/A converter board installed in the North Star computer. This board provides eight 

single-ended channels and has a built-in auto-zero capability for elimination 

of dc-offset error. 

The board requires 12 us per sample when sampling one channel with a 

constant gain and 24 us per sample when changing channels and/or gains; its 

speed is therefore more than adequate for the desired 100-Hz sampling frequency. Its 12-bit accuracy corresponds to a 72-dB dynamic range, which 

should also be more than adequate, since the specified SNR of the TEAC R-61 

recorder is only 50 dB. The converter has gain settings of 0.5, 2.0, 8.0, 

and 32.0 and a maximum measurement voltage (after application of gain) of 

25 V; these levels are compatable with the '1.5-V maximum output of the TEAC 

R-61. The 10-MQ input impedance of the A/D/A board should also be compatable 

with the 6004 output of the TEAC R-61. 

A 100-Hz sampling frequency was selected to ensure that no information 

in the 30-Hz-bandwidth analog signals is lost in the sampling process. Since 

the 12-bit output from the AID converter can be packed into 1.5 bytes, simultaneous digitization of four channels with a 100-Hz sampling frequency produces data at the rate of 600 bytes/s. A forty-minute analog tape therefore 

produces a total 1,440,000 bytes. 

Storage 

The North Star floppy-disk drive is capable of storing data at the rate 

of 695 bytes/s on the average. However, the specified track-to-track access 

time of 40 ms is four times the 10-ms sampling interval. It is therefore 

impractical to store A/D output directly on floppy disk. 

Data can be written to the SemiDisk by the operating system at an 

average rate of 3360 bytes/s. However, data are transferred to a disk only 

in 128-byte records, and 12 ms are consumed in each such transfer. 

It is therefore necessary to transfer data to the SemiDisk on a byteby-byte basis using 1/0 calls. Tests show that the average time required to 

transfer one byte via this mechanism is 0.09 ms. The worst-case time required to transfer one byte is 0.2 ms, and occurs when a new sector and/or 

track must be specified through additional 110 calls. Transfer of the six 

bytes generated by one set of four samples therefore consumes at the most 

only 1.2 ms. Since the time required is only a small fraction of the 10-ms 

sampling period, this mode of access is satisfactory. 

2.2 A/D-INTERFACE SUBROUTINES 

The set of subroutines described in this section is used to control the 

1/0 Technology A/D/A converter board [l] and to interpret its output. The 

subroutines are listed in Appendix A. 
Access 

The converter board can be jumpered to allow control by either memory or 

110 references. The required number of addresses (128) is rather large for 

the 8-bit I/O address space of 8080-type microprocessors. Since the North 

Star Horizon disk control ler uses memory addresses in the range of E899H to 

EB99H, memory from ECggH to FFFFH is not generally utilized. The memory 

block from ECPgH to EFFFH can therefore be disabled, allowing the use of 

these addresses for control of the AID converter. The use of memory control 

results in slightly faster access than is possible with I/O control. 

Control Bytes 

Subroutine ADSET assembles the two bytes required to order an A/D conversion, and computes the addresses for each. The "upper byte" contains the 

bias used for offset-voltage correction. The "lower byte" includes the channel number, gain, and whether single-ended or differential measurements are 

to be made. 

Conversion 

Subroutine ADCNV orders an A/D conversion and returns the result. To 

allow adequate settling time after changing channel or gain settings, two 

conversions are ordered and only the output from the second is actually 

returned. Provision is included for waiting for the conversion-complete 

signal. However, in practice the various computer instructions that must 

be carried between ordering the conversion and reading out the result take up 

more time than the conversion. (The North Star Horizon running a program 

written in Microsoft FORTRAN-80 requires 26 PS to do a PEEK and 9 PS to do an 

IF.) 

Interpretation - of Output 

Subroutine ADVOLT converts the two-byte result of a conversion into a 

floating-point voltage. Eight bits from the lower byte and three bits from 

the upper byte are assembled into an integer. The integer is then converted 

into a floating-point variable and scaled. The sign is determined from bit 

four of the upper byte and the voltage is modified if necessary. Finally, 

the voltage is scaled to include the converter gain. 

-- Auto Zero 

Subroutine ADBIAS automatically finds the bias byte that brings the 

output voltage closest to zero. A binary-search procedure with a total of 

ten iterations is used. The analog-input lines must be grounded for the 

auto-zero procedure to be effective. 
Gain Selection 

Subroutine ADBIAS queries the user for the maximum anticipated input 

voltage and sets the gain to the maximum value that can be used without saturating the converter. 

2.3 SEMIDISK-INTERFACE SUBROUTINES 

The SemiDiskm is a 512-kilobyte random-access memory that responds to 

1/0 calls rather than memory references, thus allowing its use in an 8080- 

type envi ronrnent with 16-bi t memory addressing [2]. The SerniDi sk is accessed 

through four 1/0 addresses, which are used to specify data byte, byte number 

(0 - 127), sector number (0 - 15), and track number (0 - 255). The set of 

five subroutines listed in Appendix A and explained below provide convenient, 

rapid access to the SemiDisk. 

Initialization 

The SemiDisk is initialized by subroutine SDINIT. SDINIT sets up the 

1/0 addresses (8gH - 83H in this system), sets the byte, sector, and track 

counters to zero, and (optionally) zeros all SemiDisk memory. 

--- Read and Write 

Subroutine SDREAD reads one byte of data from the SemiDisk at the location currently specified by the byte, sector, and track counters. Subroutine 

SDWRIT similarly writes one data byte to.the SemiDisk. 

-- Set and Increment 

Subroutine SDSET is used to write desired byte, sector, and track numbers to the SemiDisk. Subroutine SDINCR increments these counters in an appropriate manner after a byte is read or written. 

2.4 OPERATIONAL PROGRAMS 

Operational use and checking of the A/D system is accomplished by a set 

of three programs: 

"AD" digitizes four analog inputs and stores the results in the 

SemiDi sk. 

"SF" transfers data from the SemiDisk to floppy disk at the end of 

a digitizing session. 

0 "FD" decodes data from a floppy disk and displays the corresponding 

set of voltages for each sample set. 
Program listings are given in Appendix A. Disk backup is accomplished by the 

standard CP/M utilty "COPY". 

Diaitization 

Upon start-up, the "AD" program requires the user to respond to five 

questions regarding 

o Physical -channel s e l e c t i o n , 

Maximum input voltage, 

Bias-correction method, 

Duration of processing, and 

e Sampling interval. 

Any set of four of the eight AID-converter input channels can be selected for digitization. The gain is set to the maximum value consistant with 

the maximum anticipated input voltage. Bias correction (auto-zero) can be 

accomplished by assumption, measurement, or reading correction bytes from a 

disk file. The values of the correction bytes are displayed by the program, 

and must be placed in the first line of the disk file in FORTRAN 415 format 

if that option is to be used in subsequent runs. 

The duration of the processing is specified by the number of 100-sample 

data sets (four-channels) to be processed. For the desired 100-Hz sampling 

frequency, this number is the number of seconds of data to be processed. The 

sampling frequency is set by specifiying the number of 3.3-ms intervals in 

one sampling interval. Setting this number to 3 produces a 10-ms sampling 

interval, hence the desired 100-Hz sampling frequency. 

Subroutine PACK compresses the eight-byte converter output for a fourchannel data set into six bytes for storage. As shown in Figure 2-3, the 

four lower bytes (in which all eight bits are utilized) are stored without 

change. However, the low-order nibbles from four upper bytes are packed into 

two bytes. 

Timing is controlled by the real-time clock in the North Star Horizon, 

which is jumpered to generate 3.3280-ms intervals; three such intervals and 

the re1 ated instructions produce the desired 10-ms sampl i ng interval (Figure 

2-4). The total time required for sampling, packing, and storing the data is 

2.0 ms (1.4 ms for sampling, 0.6 ms for packing and storing). To allow for 

any possible variation in time requirements, sampling is undertaken during 

the first 3.3-ms interval and packing/storing is undertaken during the second 

3.3-111s interval. One or more 3.3-ms intervals are added to obtain a sampling 

interval of 10-ms or longer duration. Tests show the resultant 100-Hz sampling frequency to be accurate to within '0.001 Hz. 

Each time it is used, the program requests from the operator an ASCII 

data-identification label, whose length can be up to 64 characters. The 

first 128 bytes of data stored on the SemiDisk (hence subsequently on floppy 

disk) include this label and other parameters such as gain, number of samples, and duration of the sampling interval. The format of this data block 

is shown in Table 2-1. 
CLOCK FLAG J-A-l-l 

TIME 

BYTE NUMBERS 

Figure 2-4. Timing. 

USE 

Data-identification label - ASCII text 

Gain code (0: 0.5, 1: 2.0, 2: 8.0, 3: 32.0) 

Number of 100-sample sets 

Number of 3.3-111s intervals per sample 

Unused, set to zero 

Table 2-1. Data-identification block. 
rranstr>r to -Flo_pey -- - 

Disk 

I'r-c~cjrdi*, "5F " transfers data from the Sem1l)isk to a floppy disk. The 

proyrairi rellllires the operator to specify the nurnber of 100-sample sets and 

t-he fileriami. of the floppy-disk file. Data are transferred to the disk in 

packed for111 to rlaxirni zc its storage capabi 1 i ty. The operating system irnposes 

d read-after-write check to ensure data integrity. 

- Data Decoding 

The data stored on a floppy disk by prograln "SF" can be decoded and displayed by prograq "FD". The data-identification block is decoded first, and 

the label, gain, number of 100-sample sets, number of 3.34s samples, and 

sarnpl iny frequency are displayed on the console. The program then begins 

unpacking and decoding the data. Each set of four voltages is displayed 

along with the sample-set and sample numbers. 

2.5 OPERATIONAL PROCEDURES 

The mdxirnum length of a digitizing session is limited by the 338-kilobyte (346,112-byte) capacity of the floppy disk. For a sampling frequency of 

100 Hz, the maximum session duration is 576 s (9.61 minutes). Prior to beginning to digitize a new tape, the operator must determine the duration of 

the recording and break it into a suitable set of two or three intervals. 

Convenient session lengths are 540 s (9.0 minutes) and 570 s (9.5 minutes). 

To digitize the data, the operator loads the program, answers its questions, and readies the recorder. The operator then starts the recorder and 

immediately thereafter presses "RETURN" on the console to start digitizing 

the data. As the end of the session approaches, the operator must be ready 

to stop the recorder upon the bell signal from the console. 

Data must be transferred to the floppy disk before any additional data 

are digitized. The operator must be sure to supply "FD" with the correct 

number of data sets. Transfer of data to the floppy disk is only slightly 

faster thar; d~gitization (e.g., 555 s vs. 575 s). If a disk-write error 

occurs, the transfer must be repeated. 

The operator should then make a dup:icate (backup) of the disk usi~y the 

CP/M "COPY" program. A disk copy takes about two minute:. 

Aft,er dl1 data are digitized and stored in duplicate on floppy disks, 

they can be transferred to the Develco computer. This transfer is accorns1 i shed via RS-232 ports and the communication program described in Appendix 

3. 
2.6 REFERENCES 

1. A/D/A Analog to Digital and Digital to Analog Converter Board User 

Manual, P/N 52748-900, Revision 5, 3/83, I/O Technology, Canyon Country, 

CA, 1982. 

2. The 4931 512K S-100 SemiDisk and The 4931-1 lMbyte S-100 SerniDisk 

Manual, SemiDisk Systems, Beaverton, OR, 1982. 
CHAPTER 3. ACQUISITION 

ELF noise data were acquired during a visit by the author to the 

Develco, Inc. plant in Sunnyvale, California. The principal objectives were: 

o To digitize the analog noise recordings, 

o To store the data on North Star diskettes for use by GMRR, and 

e To transfer the data to an HP-1000 for use by Develco. 

The digitization hardware, software, and procedure are described in Chapter 2. 

This chapter describes the noise data and the software for their transfer to 

the HP-1000. 

3.1 NOISE DATA 

Fourteen analog tapes (Table 3-1) containing two-site noise data were 

digitized. The first three channels of these tapes correspond to the X, Y, 

and Z axes of the remote-antenna. The fourth channel is either the X, Y, or 

Z axis (Table 3-1) of the primary antenna. 

Each tape contains approximately thirty minutes of noise data, and is 

digitized in three 575-s (9 min, 35 s) segments and one 120-s (2 min) segment. All tapes begin with a sequence of 2 1-V calibration signals. Digitization was accomplished with the converter gain set at 2.0 to allow for 

the maximum ' 1.5-V output of the recorder. 

The primary and remote-reference antennas were located at Stone Canyon 

and Bear Valley, respectively. These sites are near Hollister, CA, and are 

separated by approximately eight kilometers (Figure 3-1). General recording 

and digitization parameters are 1 isted in Table 3-2 

The fifteenth data set shown in Table 3-1 contains three examples of 

system noise. Set # 15-1 is a digitization of a recording of the "0-V" 

calibration signal of the recorder. Set # 15-2 is a digitization of a 

recording made with no calibration signal and the recorder inputs shorted. 

Set # 15-3 is a digitization of the intrinsic converter noise with the converter inputs shorted. All three system-noise data sets contain 575-s of 

data and all four channels were subjected to the same conditions. 

The voltage output from the recorder can be converted into received 

magnetic-field intensity in Amperes per meter through a two-step process. 

First, the voltage must be divided by a factor of approximately 0.29, which 

corresponds to the measured calibration-signal voltage. The result is then 

multiplied by 12.10~ to obtain the magnetic field at the antenna. 

Develco engineers anticipated typical noise 1 eve1 s of -120 dB relative 

to 1 A/m in a 1-Hz bandwidth. However, the observed noise levels were in the 

considerably less; e.g., -141 dB, -137 dB, and -147 dB for the X, Y, and 

axes, respectively. Tests using two nearby antennas to cancel ambient noise 

show a system-noise level of -153 dB. 
TAPE 10 PRIMARY PRIMARY BACKUP SECONDS FILE- TAPE CALI BRAT I ON 

AX IS Dl SK Dl SK RECD NAME COUNTER S l GNAL 

1 Y Jan. 14, 1983, 16:45. Good weather, some disturbance. 

52 7 479 575 ND1-1 Y 

2 Z Jan. 14, 1983, 18:OO. Good weather. 

53 1 483 575 ND2-1 000 -219 Y 

53 2 484 575 ND2-2 219 - 376 N 

533 485 575 ND2-3 376 - 507 N 

534 486 120 ND2-4 507 - end N 

Jan. 5, 1983, 16:OO. 

53 5 487 575 ND3- 1 000 - 196 Y 

536 488 575 ND3-2 196 - 351 N 

53 7 489 575 ND3-3 351 - 482 N 

538 490 120 ND3-4 482 - end N 

Jan. 11, 1983, 05:30. Good weather. 

539 49 1 575 ND4-1 000 - 220 Y 

540 492 575 ND4-2 220 - 378 M 

54 1 493 575 ND4-3 378 - 509 N 

54 2 4 94 120 ND4-4 509 - end N 

Jan. 11, 1983, 06:20. Good weather. 

543 495 575 ND5-1 000-216 Y 

544 496 575 ND5-2 216 - 373 N 

545 497 575 ND5-3 373 - 503 N 

546 4 98 120 ND5-4 503 - end N 

Jan. 5, 1983, 18:OO. 

547 4 99 575 ND6- 1 000 - 219 Y 

54 8 500 575 ND6-2 219 - 377 N 

54 9 50 1 575 ND6-3 377 - 508 N 

550 50 2 120 ND6-4 508 - 535 N 

Jan. 11, 1983, 07:30. Good weather. 

55 1 503 575 ND7- 1 000 - 218 Y 

552 504 575 ND7-2 218 - 375 N 

553 505 575 ND7-3 375 - 505 N 

554 506 120 ND7-4 505 - 533 N 

Jan. 7, 1983, 09:40. Good weather. 

555 50 7 575 ND8- 1 000 - 21 9 

556 508 575 ND8-2 219 - 378 

557 509 575 ND8-3 378 - 510 

558 5 10 120 ND8-4 510 - 536 

Jan. 5, 1983, 11:lO. 

559 51 1 575 ND9-1 000 - 220 

560 512 575 ND9-2 220 - 378 

56 1 513 575 ND9-3 378 - 510 

562 514 120 ND9-4 510 - 536 

Jan. 7, 1983, 11:20. Good weather. 

563 51 5 575 ND10-1 000 - 219 Y 

564 516 575 ND10-2 219 - 377 N 

56 5 517 575 ND10-3 377 - 507 N 

566 518 120 ND 10-4 507 - end N 

Table 3-1. List of tapes, fi IPS, and conditions 
TAPE ID PRIMARY PRIMARY BACKUP SECONDS FILE- TAPE CAL l BRAT I ON 

AXIS Dl SK Dl SK RECD NAME COUNTER S l GNAL 

11 X Jan. 6, 1983, 19:15. 

56 7 519 575 NDll-1 000 - 220 Y 

568 520 575 NDll-2 220 - 378 N 

569 52 1 575 ND 1 1-3 378 - 509 N 

570 522 120 ND 1 1 -4 509 - end N 

12 Y Jan. 12, 1983, 19:30. Good weather. 

57 1 523 575 ND12-1 000 - 219 Y 

572 524 575 ND12-2 219 - 376 N 

573 52 5 575 ND12-3 376 - 506 N 

574 526 120 ND 12-4 506 - 533 N 

13 X Jan. 14, 1983, 14:OO. Good weather, disturbed data. 

575 58 7 575 ND13-1 000 - 219 Y 

576 588 575 ND13-2 219 - 376 N 

577 589 575 ND13-3 376 - 507 N 

578 590 120 ND13-4 507 - end N 

14 Y Jan. 7, 1983, 10:30. Good weather. 

579 59 1 575 ND14-1 000 - 219 Y 

580 592 575 NDt4-2 219 - 376 N 

58 1 593 575 ND14-3 376 - 507 N 

582 594 120 ND 14-4 507 - end N 

15 Noise diagnostic measurements - same conditions all four channels 

583 

N/A 584 

595 575 ND15-1 Ca I i b ftO-Vw Y 

N/A 585 

596 575 ND15-2 Rcdr, shorted Y 

N/A 597 575 ND15-3 UD, shorted N 

Table 3-1. List of tapes, files, and conditions (continued). 
Figure 3-1. Antenna sites. 

AID CHANNELS 

# 1: Remote antenna, X axis 

# 2: Remote antenna, Y axis 

# 3: Remote antenna, Z axis 

# 4: Primary antenna, specified axis 

ANTENNA SITES 

GAIN 

Gain = 2.0 

Remote: Bear Valley, CA 

Primary: Stone Canyon, CA 

Table 3-2. General recording/digitizing parameters. 
3.2 STORAGE AND TRANSFER 

The digitized data from each analog tape were stored on a set of four 

floppy disks. The initial portion of each data subset was validated by inspection with the file-display program (FD). 4 backup copy of each disk wds 

made. A total of 16.4 megabytes of digitized noise data was produced. 

The data were transferred to Develco's HP-1000 corqputer via an 2s-232 

interface. While this proved to be rather slow, the total time consumed was 

less than that required to implement a high-speed parallel interface. 

Communication Protocol -- 

A binary-bl ock protocol was used to maxirrlize the data-transfer rate. 

This protocol is based upon the transmission ot a block of 128 8-bit data 

bytes, followed by a single parity-check byte. The check byte is obtained by 

exclusive-0Ring the 128 data bytes, and detects all sinyls errors. The receiving computer compares the received check byte with one computed from the 

received data bytes and issues a positive or negative acknowledgec~ent. The 

transmission is repeated up to two additional times when a negative acknow1 edgement i s recei ved. 

To accommodate the requirements of the HP-1000 (receiving computer), the 

protocol, as seen by the transmitting (North Star Horizon) comput~r, ccnsi;ts 

of the following steps: 

0 Wait for D, from HP-1000. 

o Send one character: 

3 to start data block 

sB to end file 

E, to end program. 

o Send 128 data bytes. 

e Send check byte. 

o Await E~, send +. 

B Await E~, send %. 

0 Receive 

A for ~uccessful receptiorl. K 

for unsuccessful reception. 
Await 5, send +. 

0 Await 5, send \. 

Await 5, send %. 

Program Operation 

The transmitter program is listed in Appendix B. 

The transmitter program must be started first. Before data can be 

transferred, the operator must supply I/O-port parameters in response to the 

questions asked by the program. The I/O-port information is given in the 

CP/M BIOS listing and the computer manual. 

During operation, the program displays the number of the disk record 

that is currently being transferred. 

The program can be permanently configured for a particular computer by 

using the CP/M-80 "SAVE" command. The required steps are: 

o Find the program size (in 128-byte disk records) using "STAT". 

o Divide the number of disk records by two to obtain the number of 

256-byte memory pages. 

Run the program and supply the desired 1/0 parameters. 

e Interrupt the program by typing "5;" (control-C, ASCII #3) in 

response to "FILENAME = ". 

o Use "SAVE" to save the memory image in a disk file. 

A program that has been "SAVEd" will never again request 1/0 information from 

the operator. 
CHAPTER 4. PROCESSING 

Software for testing the wideband NLP/ANC (nonlinear-processing / adaptive-noise-cancellation) algorithm with digitized ELF noise was developed 

through modification of the simulation program [l, Chapter 71 developed in 

Phase-I1 of this program. The specific steps in this effort were: 

e Modificaticn of the NLP/ANC simulation program for a single-channel 

primary antenna, 

* Inclusion of an optional simulation of quantization by the A/D 

converter, 

e Further modification of the single-channel simulation program to 

out-put simulated digitized noise to a disk file, 

o Further modification of the single-channel simulation prograrn to 

obtain its noise inputs from a disk file, and 

o Modification of the disk-file processing prograln to add a fixed 

(nonadapti ve) noi se-cancel1 ation processor. 

In addition to the development of the processing program, cursory 

examinations of the differences between channels, effects of quantization, 

and effects of cancellation method were made. 

4.1 SINGLE-CHANNEL SIMULATOR 

The digitized ELF noise (Chapter 3) includes all three axes of the remote-reference antenna, but only one axis is the primary antenna. The wideband NLP/ANC simulation program [l, Chapter 71 was therefore modified for 

operation with a selectable (X, Y, or 2) single-axis primary-antenna input. 

A limited number of simulations were conducted using parameters analogous to those of Figure 7-24 of [I]. However, the signal and noise amplitudes were scaled to resemble those of the digitized ELF noise. The results 

are shown in Figures 4-1 and 4-2 and Table 4-1. It is apparent that the behaviors of estimators for the X and Z axes are similar. However, any correlation between the behavior of the Y-axis estimate and those of th~ other two 

axes is rather vague. 

4.2 QUANT I ZAT I ON 

The single-channel NLP/ANC-simulation program includes optif~r~al s~mulation of qudntization and limiting by the A/U converter (Chant~r 2). The results are shown in Table 4-1 and Fiqures 4-3 - 4-6. 

The convergence of anquantized and qudntized simulations are, in general, similar. Differences in the processing gains are generally no rlore than 

a v. - 

13 w 

- .- 

U 

wow 

>wo 

a "1 

2 r C c. -- WC 

Yl rO V) 

.- 3 C 

0 Lr .- zcm 3 m 

-4.4~ mmm -om u.+w u>-ID .++a .................. 

000 om0 wr-u uom m-.+ 000 

NN N I 4 mm NO mm 

NNO cwo m~r- mmr- wmm oh- .................. 

be0 had Urn- -ON ah4 000 4- 4- NN 4- NN 

000 h~rn ~mw mmm =curd -0-4 

.................. 

wwc N-I- mow om- mmo ooo 4- a + 4 I .-,A -4, 

N+N em- m-N .+win NCUU 0p.m .................. 

NNO 4-40 NmO mh?- WmO hWO -+ 1 00 - m~ I NN I mm I 

mw~ am+ m-N hem mwr- mhw .................. 

moo hm- NU- ~ww N-o mmo I 44 NhI 44 I NNI 

wuo wwo ohm UNN m~e -mcv .................. 

mmo mew- omo N++ ma0 mwo .-,A, - I -A, -- .+-I 

mm mu mm urn mm uu 00 OC 00 00 00 00 8, I I I, t, , , 4, g WW WW WW WLLJ WW 

CCN -N wm m NN 

rn 0-3 Um U3N hN COW m m mm ~m mw -m mm 4m NN r-+ 4- be . . . . . . . . . . . . mm om mm mm mm no 00 00 00 00 00 00 

, I I I I I i # I I I I 

WW WW WW WW ww ww Oh Om 030 Ud kn.3 Um m N- -N mr- a-w .om we m m m- ww mw m ?- mu m mu ww . . . . . . , . . . - . . NN NN NN NN NN NN 

CO 00 00 00 00 00 

f I 8 I I, I I t I , I 

WW WW WW WW WW WW 

mn m wo m- m- -m mm m mi. NO mm mw ua N- N- -N r-m mm em -4 om NN NCL 

. . . . . . . . . . . . 

30 30 30 3W 30 20 

'TW IS EZ EX ZZ 22 

XX XX >t XX NN XX 

-- - - ---- -.NC me0 mwo hmc moo ,+NL 4 4 4 

-. 

W u 

4 

> 

a 

3, 

&. 

u 

z> 

a 

u 

h 

0 

LA 

- 

LL 

-x 

u 

w 

4% 

Zb 

La 

W 

U 

1 

13 

2 Z 

V1 - 

X 

< 

u. 

m~ 

wz 

C 


a deci be1 , especial ly for Gaussian noise and low-level atmospheric noise. 

For moderate- and high-level atmospheric noise, somewhat larger differences occur, especially in the ANC-processing gain. The principal cause of 

the increased differences appears to be clipping of large-amplitude noise 

spikes by the A/D converter. It is apparent in Table 4-1 that the estimation 

error for linear processing is nearly identical for Gaussian and low-level 

atmospheric noise, but differs for moderate- and high-level atmospheric 

noise. 

Analog noise was digitized with an A/D-converter gain of 2.0, which 

results in a converter clipping level of 2.5 V (Chapter 2). Since the maximum voltage recordable on the TEAC R-61 is 1.5 V, the A/D converter introduced no additional clipping of the ELF noise. 

4.3 DIGITIZED-NOISE PROCESSOR 

The NLP/ANC-processor program for digitized ELF noise was obtained by 

splitting the simulation program into two programs (Appendix C). 

The first of these two programs generates simulated ELF noise as in the 

single-channel program. The noise is then quantized, packed, and written to 

disk in the format used to store the digitized ELF noise. 

The second program reads digitized, packed noise data from a disk file. 

The data are then unpacked, converted to floating-point numbers, and processed as in the simulation program. 

The program for processing digitized noise includes a subroutine for 

processing calibration signals when they are present. The operator must 

first use the file-display program (Chapter 2) to determine the duration (in 

seconds) and approximate voltage level of the 1-V calibration signal. The 

subroutine reads and decodes the data file for the specified number of 

seconds. Voltages fa1 ling within ten (10) percent of the specified value 

during the specified time interval are averaged to estimate the true level 

of the 1-V calibration signal in each channel. All subsequent voltages are 

divided by the appropriate averaged calibration voltage to scale them to 

true voltages. 

Proper operation of these prograrns was verified by comparing their 

results with those obtained by the single-channel sirnulation program. 

4.4 FIXED CANCELLATION 

A fixed-noise-cancellation (FNC) processor was added to enable comparison of adaptive- and nonadaptive-cancellation techniques. The linear FNC 

processing occurs after independent nonlinear processing is applied to both 

the primary and reference inputs. The output from a selected axis of the 

reference antenna is then subtracted from the output of the primary antenna. 
The results from one set of simulations with moderate-level atmospheric 

noise are shown in Figures 4-7 - 4-9. Convergence patterns are quite similar 

for both ANC and FNC. ANC usual ly (but not always) produces slightly more 

noise reduction than FNC. However, this simulation does not show the true 

advantages of ANC over FNC, since it employs identically oriented antennas 

and identical ground filters at each site. 

4.5 REFERENCE 

1. F. H. Raab, "Adaptive-noise-cancellation techniques for through-theearth electromagnetics, Volume 11," Final Report GMRR TR83-1, Green 

Mountain Radio Research Company, Winooski, VT, September 1983. 
a. 

9 omo 

. c., 0 ." 

0 I v C 

. * CU 

1 & YI 

9.- u 

. d 2-w 

. . .m ?4 OYW 

h - urn", U" E 

Y - m - _LN. 

9 ". ,: ,Ux: 0 N Y1100 - "LO w " . -4 -- .mu - -00-I h > .u .. . " I OW 0- C EV 

- - !in: 2: a r ., *IC. VI 
CHAPTER 5. RESULTS 

Data from eleven of the fourteen ELF-noise tapes were used to evaluate 

the performance of the nonl inear-processing (NLP) and adapti ve-noise-cancellation (ANC) algorithms. The characteristics of the data are discussed in 

Chapter 3, and the processing algorithms are described in Chapter 4. This 

chapter presents the results of the evaluation. 

5.1 USEABLE DATA 

Three of Develco's tapes (#I, #3, and #13) were "disturbed" by anomolies 

such as communication outages and cows rubbing against one of the receiving 

antennas. Since specific times of the disturbances are not known, data from 

these three tapes were not used in assessing average performance of the signal -processing a1 gorithms. 

The data-transfer procedure divided each analog tape into three useable 

disk files, each containing 540 s of ELF-noise data. Since obtaining reliable average-performance figures requires a large number of processing runs, 

each of the three useable disk files is subdivided into two 270-s runs. This 

division produces a total of 66 processing runs (24 X-axis, 18 Y-axis, and 24 

2-axi s). 

Because the fourth disk file for each tape contains no more than 120 s 

of data and the end of the tape creates a "disturbance," the fourth disk 

files are not processed. Since data set #9-1 contains a very long cal ibration signal, it is divided into two 260-s runs. 

5.2 CONVERGENCE 

The convergence of the NLP/ANC algorithm with real ELF noise is similar 

to that with simulated ELF noise. Examples of "normal" convergence trajectories are shown in Figures 5-1 - 5-3. 

Convergence trajectories are, however, highly varied. For example, the 

trajectories of Figures 5-1 and 5-4 are derived from adjacent data segments 

of data set #8-1. An example of a "disturbance" and the algorithm's recovery 

can be seen in Figure 5-5. 

5.3 PROCESSING GAIN 

The averages for each processing run are assembled into a file, which is 

then processed to obtain averages, standard deviations, minima, and maxima. 

These results are computed as applicable for each noise tape, each axis, and 

overall. The results for individual data sets are presented in Appendix D. 

A number of tests (not shown) were also conducted to verify that the signal 

amplitude has an insignificant effect upon the signal-processing gains. 
Figure 5-1. Convergence plot for data set #8-1A (X axis). 

G 

n O Gf 

Figure 5-2. Convergen.ce plot for data set #12-2A (J' axis). 
Figure 5-3. Convergence plot for data set #6-2A (2 axis). 

. G 0 

a 

Figure 5-4. Convergence plot for data set #8-1B (X axis). 

The overall results are shown in Table 5-1 and Figures 5-6 - 5-11. Figure 5-6 shows the signal-processing gain achieved by nonlinear processing and 

ANC as a function of the adaptive-clipped threshold T. NLP gain Gn is, by 

itself, insignificant. ANC gain Ga reaches a maximum value of 4.1 dB with T 

= to disable nonlinear processing. Setting -r to finite values (e.g., the 

0.707 expected to optimize processing for moderate-level atmospheric noise) 

resuits in considerable degradation of Ga and total signal-processing gain 

Average signal-processing gains for each axis are shown in Figures 5-7 - 

5-9 and compared in Figure 5-10. For the X axis, performance is similar to 

the all-axis performance, with Ga reaching a maximum valile of 4.5 dB. For 

the Y axis, the maximum value of Ga is 11.9 dB; however, for the Z axis, the 

maximum value of Ga is -1.9 dB (a 1.9-dB inerease in estimation error). Some 

nonlinear processing produces a very slight increase in total signal-processing gain for the X axis; however, in general nonlinear processing only reduces the total gain. 

Histograms of the average ANC gains are presented in Figures 5-11 - 

5-14. It is apparent from the histograms (and the standard deviations in 

Table 5-1) that there is a good deal of variability in the average gains. 

This is especially true for the all-axis averages, which are based upon inherently different individual-axis averages. For example, the standard deviation of the 4.1-dB all-axis average gain is 7.1 dB, whereas the standard deviation of the 11.9-dB Y-axis average gain is 6.1 dB. The averages are therefore significant but their values must not be considered exact. 

The gains achieved by adaptive noise cancel lation (ANC) and fixed noise 

cancellation (FNC) are compared in Figure 5-15. 0verall;ANC is superior to 

FNC by 0.3 dB. However, for the Y axis (which has more noise power to be 

cancelled), Ga exceeds G by 1.6 dB. A gain of 1.6 dB is equivalent to a 30- f 

percent reduction in signal-integration time. Under adverse conditions, Ga 

values as low as -6.7 dB are observed; in contrast, G values as low as f 

-15.0 dB are observed. Under good conditions, Ga can be as high as 26.2 dB, 

while the maximum observed value of G is 22.6 dB. f 

5-4 WEIGHTING VECTOR 

The average of the ANC weighting vectors is depicted in Figure 5-16 and 

Table 5-2. It is apparent that cancellation of the noise in the X and Y axes 

relies heavily upon the X and Y axes of the remote-reference antenna. However, the use of the other axes is not insignificant. The weighting vector 

for cancellation of Z-axis noise varies greatly and on the average is distributed evenly among all three remote-reference inputs. 
Statistical analysis of disk-ANC numerical data, Program # FHR-354- 

NUMERICAL DATA - AVERAGES 

Single-channel NLP/ANC algorithm with disk-file input 

Program # FHR-353E 

Noise-file parameters: A/D gain = 2.0 Sampling frequency = 100.0 

Total seconds in noise file = 575 

Processing parameters: Signal = .I414 +j 0.000 

Nsec = 270 Kout = 2 Nrun = 2 

Nset = 2 F3dB = 1.00 ............................................................................. ............................................................................. 

Threshold = infinite 

Set Axis Date Time Comnent NSC NSS Gn Ga Gta G f Gtf 

AVERAGES 

X axis ( 24 data sets) 0.0 4.3 4.3 3.6 3.6 

Y axis ( 18 data sets) 0.0 11.9 11.9 10.3 10.3 

Z axis ( 24 data sets) 0.0 -1.9 -1.9 -.8 -.8 

Overal l ( 66 data sets) 0 .O 4.1 4.1 3.8 3.8 

STANDARD DEVIATIONS 

Xaxis (24 data sets) 0.0 5.2 5.2 7.8 7.8 

Y axis ( 18 data sets) 0.0 6.1 6.1 5.5 5.5 

Z axis ( 24 data sets) 0.0 2.4 2.4 2.4 2.4 

Overal l ( 66 data sets) 0.0 7.1 7.1 7.1 7.1 

MINIMA 

X axis ( 24 data sets) 

Y axis ( 18 data sets) 

Z axis ( 24 data sets) 

Overall ( 66 data sets) 

MAX I MA 

X axis ( 24 data sets) 0.0 12.6 12.6 13.9 13.9 

Y axis ( 18 data sets) 0.0 26.2 26.2 22.6 22.6 

Z axis ( 24 data sets) 0.0 2.9 2.9 5.3 5.3 

Overall ( 66 data sets) 0.0 26.2 26.2 22.6 22.6 ............................................................................. 

............................................................................. 

Threshold = 2.828 

AVERAGES 

X axis ( 24 data sets) .1 4.2 4.2 3.6 3.7 

Y axis ( 18 data sets) .3 10.7 11 .O 5.1 5.4 

Z axis ( 24 data sets) -1.2 -1.5 -2.7 -1. -2.2 

Overal l ( 66 data sets) -.3 3.9 3.6 2.3 2 .O 

STANDARD DEVIATIONS 

X axis ( 24 data sets) .6 5.6 5.8 7.9 8.1 

Y axis ( 18 data sets) 1.6 7.2 6.6 6.9 6.3 

Z axis ( 24 data sets) 2.7 2.3 2.7 2 .0 2.8 

Overal l ( 66 data sets) 1.9 7.1 7.5 6.5 6.9 

MINIMA 

X axis ( 24 data sets) 

Y axis ( 18 data sets) 

Z axis ( 24 data sets) 

Overall ( 66 data sets) 

MAX I MA 

X axis ( 24 data sets) 1.6 12.9 13.1 14.7 15.1 

Y axis ( 18 data sets) 6.3 25.3 25.2 12.8 12.7 

Z axis ( 24 data sets) 4.4 3.3 3.2 3.2 3.2 

Overal l ( 66 data sets) 6.3 25.3 25.2 14.7 15.1 ............................................................................. ............................................................................. 

Note: Data sets #1, #3, and #13 excluded due to distur 

Table 5-1. Processing-gain statistics. 
Threshold = 1.414 

Set Ax is Date Time Comnent NSC NSS Gn Ga Gta G f Gtf 

X axis ( 24 data sets) .4 4.2 4.5 2.8 3.1 

Y axis ( 18 data sets) .9 9.0 9.9 7.4 8.3 

Z axis ( 24 data sets) -1.1 -1.6 -2.6 -1.3 -2.4 

Overal l ( 66 data sets) -.O 3.4 3.4 2.5 2.5 

STANDARD DEVIATIONS 

X axis ( 24 data sets) 

Y axis ( 18 data sets) 

Z axis ( 24 data sets) 

Overall ( 66 data sets) 

MINIMA 

X axis f 24 data sets) 

Y axis ( 18 data sets) 

Z axis ( 24 data sets) 

Overall ( 66 data sets) 

MAX l MA 

Xaxis (24 data sets) 7.8 19.5 20.9 21.5 22.9 

Y axis ( 18 data sets) 9.0 19.1 19.8 17.2 19.4 

Z axis ( 24 data sets) 8.1 3.8 7.0 3.8 4.9 

Overall ( 66 data sets) 9.0 19.5 20.9 21.5 22.9 ............................................................................. ............................................................................. 

Threshold = 0.707 

Set Axis Date Time Comnent NSC NSS Gn Ga Gta G f Gtf 

X axis ( 24 data sets) 

Y axis ( 18 data sets) 

Z axis ( 24 data sets) 

Overall ( 66 data sets) 

STANDARD DEV I ATIONS 

X axis ( 24 data sets) 

Y axis ( 18 data sets) 

Z axis ( 24 data sets) 

Overal l ( 66 data sets) 

MINIMA 

X axis 

Y axis 

Z axis 

Overal I 

MAX I MA 

X axis 

Y axis 

Z axis 

Overa l I 

( 24 data sets) 

( 18 data sets) 

( 24 data sets) 

( 66 data sets) 

( 24 data sets) 

( 18 data sets) 

( 24 data sets) 

( 66 data sets) 

Note: Data sets ,#I, 83, and #13 excluded due to disturbances. 

Table 5-1. Processing-gdin statistics. (Continued) 


u 

Ul 

r. .- YI 

1 4 

N 

9 - N E - II m 

VI 

rn 

01 u, - c P .- 

9 z vr al b- U 

0 P 

a c. 

& aJ n U1 E m a L 

c TY > - rn 

CI 4 'C 

0 0 + a L 

07 

CI m 

0 - I 7 ." 

;r: :: "u I 

u L" .- 

xi! - a 4 c- > .- L I. 

QI +- -. c 
Averaging of correlation data, Program 8 FHR-355- 

Total number of data sets = 33 

Numbers of samples for local X, Y, and Z axes = 12 9 12 

Kyy AVC = 

Rxy AVG - 

W AVG = 

i ,ihl11 3-2, i!vpr;li_jp covariavc e, ,:rosl; correla t iori, dnc; wei ghtinq. 
5.5 NOISE CORRELATION 

The output file produced from the noise-data processing program (Appendix C) includes the estimated remote-reference covariance matrix R crossYY ' 

correlation vector and weighting vector w for each run. These data 

XY , were assembled into a file and then averaged. 

The resultant average covariance and cross-correlation matrices are 

given in Table 5-2 and represented graphically in Figures 5-17 and 5-18. On 

the average, the Y-axis noise power exceeds the X-axis noise power by 1.3 dB 

and the Z-axis noise power by 10.7 dB. While the X-axis noise power exceeds 

the Y-axis noise power in some data sets, the Z-axis noise power is consistently about 10-dB lower than the X-axis and Y-axis noise powers. 

Axis-to-axis correlation is readily observable in the Figures 5-17 and 

5-18. To compare correlations, it is necessary to normalize the off-diagonal 

terms of k by the corresponding noise powers; e.g., Y Y 

The results thus obtained from the values in Table 5-2 are: 

The cross-correlation matrix (Figure 5-18) is generally similar to the 

covariance matrix (Figure 5-17). The lack of 2-axis-to-Z-axis correlation is 

apparent. Phase shift in the X-X and Y-Y elements is also apparent. The 

magnitude of the quadrature component of these elements is a measure of the 

effect of the local geology. It is interesting to note that the magnitudes 

of the imaginary part of the X-X and Y-Y elements is comparable to the magnitude of the 2-Z element. This is not suprising, since the Z-axis magnetic 

field is produced by the local geology. 

5.6 ANALYSIS 

The ANC algorithm works as it should, reducing the effective noise level 

by 4.1 dB on the average. The improvement is greater when the noise level is 

higher (e.g., 11.9 dB for the Y axis) and lower when the noise level is lower 

(e.g., -1.9 dB for the Z axis). 

Overall, the ANC gain exceeded the FNC improvement by only 0.3 dB. However, as atmospheric noise power (hence correlation) increases, the difference between ANC and FNC gains also increases. The most significant advantage of ANC is its ability to prevent disastrous SNR degradation under adverse conditions. ANC is therefore definitely desirable. 

Nonlinear processing was on the average harmful for this set of noise 

data. The gain Gn was relatively small and more than overcome by the associated decrease in Ga. 

The character of the noise should be determined from statistical 

analysis. However, some conclusions about its character can be drawn from 

t,he behavior of the ANC and NLP algorithms. The small or negative gains 

achieved by nonlinear processing suggest that the Z-axis noise is Gaussian 

and that the X- and Y-axis noises are only slightly impulsive (not even 

"low-level ," however). The observed correlation matrices and ANC gains are 

consistant with local, Gaussian noise whose amplitude is 10-dB below that of 

the atmospheric noise. These observations on the character of the noise are 

consistant with those made by Develco engineers (Chapter 3). 

Because of the relatively high values of the standard deviations, the 

averages derived here must be considered approximate indications of what is 

to be expected. A considerably greater number of noise-data sets will be 

required to improve the accuracy of the averages. 
CHAPTER 6. CONCLUSIONS AND RECOMMENDATIONS 

The extension of through-the-earth electromagnetic location techniques 

to deep-mine depths (1-km) requires more noise reduction than can be achieved 

(in an acceptable amount of time) by simple signal integration. Intrinsicsafety requirements preclude significant increases in the power of the underground transmitter. The development of a signal-processing t~chnique such as 

adaptive noise cancellation is therefore clearly necessary. 

6.1 ACCOMPLISHMENTS 

Phase I of this prograrn was a preliminary theoretical investigaticn of 

ELF-noise characteristics and ANC techniques. The nature of the multicomponent El-F-noise model was defined. The available noise-correlation data were 

reviewed and the results suggested that from 10 to 40 dB of noise reduction 

could be achieved by the use of ANC. DM1 and LMS techniques were iaentifled 

as candidate ANC algorithms. 

Phase I1 of this program was a comparison and evaluation of techniques 

through simulation. Simulation software was implemented for both the multicomponent ELF-noise model and the NLPIANC signal processor. Nonlinear processing techniques were compared and the adaptive clipper was found to offer 

the greatest processing gain. Adaptive-noise-cancellation techniques were 

compared and DM1 was found to be superior to LMS. Simulation results showed 

that it is desirable to include both NLP and ANC in the signal processor, but 

the gains of the two processors are not additive. Noise reductions ranging 

from 5 to 25 dB were predicted for typical conditions. 

The objective of Phase I11 is to evaluate the signal -processing algorithm with real ELF noise. The principal accomplishments during Phase I11 of 

this prograrn are: 

o Implementation of hardware and software to digitize ELF noise, 

o Digitization of fourteen half-hour analog tapes of ELF noise, 

Modification of the Phase-I I signal -processing a1 gori thm to work 

with a single-channel primary antenna and to accept disk-file input, 

Processing of the digitized data with various clipping thresholds, 

and 

s Compilation and analysis of the results. 

6.2 RESULTS AND CONCLUSIONS 

The all-axis average processing gain is 4.1 dB. When processing gains 

are averaged on an axis-by-axis basis, ANC gains of 4.3, 11.9, and -1.9 dB 

are obtained for the X, :', and Z axes, respectively. The ANC gain roughly 
fol lows the received noise power in each axis. 

Nonlinear processing by itself produces a small gain in some cases. 

However, the resultant degradation of the ANC gain is far larger than the NLP 

gain. The greatest total gain was achieved with nonlinear processing disabled. 

The number of data sets available is somewhat small, all data sets were 

recorded during the same season, and the noise levels were in general 15 to 

20 dB lower than nominal. The results obtained should therefore be regarded 

as pessimistic approximate indications of what can be achieved. However, 

the results clearly show that ANC works. Under high-noise-power conditions, 

ANC can make a significant improvement in the SNR of through-the-earth electromagnetic systems operating at extremely low frequencies. 

6.3 RECOMMENDATIONS 

An effective signal-processing technique offering SNR enhancement in 

TTE EM systems has been developed by this program. However, there remain 

several unanswered questions and conjectures, as well as a need for more 

ELF noise data. 

Determination of the characteristics of the digitized ELF noise was not 

part of the Phase-111 efforts. From the performance of the signal-processing 

algorithm, the noise appears to be only slightly impulsive (somewhere between 

Gaussian and low-level impulsivity, in contrast to the anticipated moderatelevel impulsivity). A statistical evaluation of the digitized noise would 

allow validation of the multicomponent model and determination of its parameters for use in future simulations. This evaluation should include determination of histograms, moments, best fit to the Field-Lewinstein-Modestino 

model , and di rection-of-arri val characteristics. 

The digitized noise is believed to represent quiet atmospheric conditions. Additional measurements made at various times throughout the year 

are required to ascertain what to anticipate for normal and worst-case operational environments. 

Additional data sets would also enable a more accurate assessment of the 

capabilities of the signal-processing system. It is desirable (although not 

absolutely necessary) to have all three axes of the primary antenna recorded 

simultaneously. 

The results presented in Chapter 5 show the potential adverse impact 

of an incorrect setting of the adaptive-clipper threshold. Disabling the 

nonlinear processing is not advisable, since operation may be required in 

the presence of thunderstorms or machine-generated impulsive noise. In 

addition, ELF noise is expected in general to be somewhat more impulsive 

than that used in this evaluation. A possible solution is to estimate both 

the second and fourth moments of the received noise and set the thresold to 

correspond to both noise power and impulsivity. Research prn this two-deqreeof-freedom adaptive clipper is needed. 
The signal -processi ng a1 gori thm developed by thi s program employs a 

narrowband ANC algorithm. Nonlinear processing is applied independently to 

each wideband input, which is then averaged. The ANC processing then uses 

the estimated noise covariance and cross-correlation to combine the averages 

to minimize the expected estimation error. 

An alternative approach employs a wideband ANC processor. In this approach, wideband inputs are combined by the ANC processor, whose output is 

then applied to the nonlinear processor. The wideband ANC system requires 

a recursive algorithm, since the matrix inversions required by the DM1 algorithm make it unsuitable for real-time use. It is not anticipated that the 

wideband algorithm would make the ANC and NLP gains additive, since the ANC 

processor should reduce greatly the input to the nonlinear processor. Additional study and simulation are required to determine whether or not the 

performance of the wideband algorithm can be expected to be better than that 

of the narrowband a1 gori thm devel oped by thi s program. 
APPENDIX A. DIGITIZATION.PROGRAMS 

AID-Interface Subroutines 

C Program # FHR-342C, 18/28/83 

C Interface subroutines for 1/0 technology A/D/A board 

C 

SUBROUTINE ADVERS(LUN1 

C Print version number 

IMPLICIT LOGICAL (L) 

WRITE(LUN,10) 

10 FORMAT(' A/D/A Interface, Program # FtR-342C1 ) 

RETURN 

END 

SUBROUTINE ADSET(LCHAN,LGAIN,LSD,LBIAS,LCBU,LCBL) 

Form control bytes and set addresses 

Configured for access from memory, base address = ABASE 

LSD = 

0 for sl ng l e-ended operat ion 

32 for differential operation 

LSR = (applicable only if D/A Is used) 

0 for single 

64 for refresh 

IMPLICIT LOGICAL (L) 

IMPLICIT INTEGOi (A-K,M-Z) 

LOGICAL PEEK 

DIMENSION LBIT(8) 

COWON /ADA/ LBIT, LSR, ABASE, ACBU, ACBL 

DATA ABASE /X1 ECOB1/ 

DATA LBIT / X1Oll, X102', X104', X108', 

> Xl 10' , X1 20' , X1401 , X1 80' / 

Addresses for upper and lower control bytes 

ACBU = ABASE + 1 

ACBL = ABASE 

Form control bytes 

LCBU = LBl AS 

LCBL = LCHAN .OR. (LGAIN*8) .OR. LSD .OR. LSR 

Add start-conversion bl t 

LCBL = LCBL .OR. LBIT(8) 

RETURN 

END 

SUBROUTINE ADCNV(LCBU,LCBL,LR2,LRl) 

C Read A/D converter 

C Configured for access from memory 

C Must ca l l ADSET prior to f l rst use of ADCNV 

C Test result, North Star Horlzon: < 0.4 ms per reading 

IMPLICIT LOGICAL (L) 

IMPLICIT INTEGER (A-K,M-Z) 

LOG1 CAL PEEK 

DIMENSION LBIT(8) 

COl+lON /ADA/ LBIT, LSR, ABASE, ACBU, ACBL 

C 

C Start first conversion 

C 

CALL POKE(ACBU,LCBU) 

CALL POKE(ACBL,LCBL) 

C 

C Read upper byte and test for settl lng 
C 

100 LTEST = PEEK(ACBU1 .AND. LBIT(8) 

IF(LTEST.EQ,B) GO TO 100 

C 

C Start second conversion 

C 

CALL POKE( ACBL, LC8L) 

C 

C Read upper byte and test for completion 

C 

200 LR2 = PEEK(ACBU1 

LTEST = LR2 .AND. LBIT(8) 

IF(LTEST.EQ.B) GO TO 200 

C 

C Conversion completed, read lower byte 

C 

LRl = PEEK(ACBL) 

L 

RETURN 

END 

SUBROUTINE ADVOLT(LR2,LRl,LGAIN,LCHAN,V) 

C Interpret output of A/D converter 

IMPLICIT LOGICAL (L) 

IMPLICIT INTEGER (A-K,M-Z) 

REAL V 

C 

C Get channel no 

ti 

LCHAN = (LR2 .AND. (X170')) / 16 

C 

C Get lower byte - 

ti 

I1 = LR1 .AND. (X'7F1 

LTEST = LRl .AND. (X180' 

IF(LTEST,NE.B) I1 = I1 .OR. (X10080' 

C 

C Get upper byte 

C 

LTEST = LR2 .AND. (X807' 

I2 = LTEST 

12 = 12 * 256 

ti 

C Assemble Into Integer 

C 

IV = I1 .OR* 12 

c Get sign bit 

C 

LSlGN = LR2 .AND. (X108') 

C 

C Convert to real 

C 

V = 2847-5 - IV 

IF(LSIGN.EQ.0) V = V - 2048.0 

C 

C Include full-scale and gain factors n 

b 

V = V * 4.8828327E-03 / 4'"LGAIN 

C 

RETURN 

END 

SUBROUTINE ADSIM(LGAIN,LCHANDV,LR2,LR1) 

C Simulated A/D output bytes for specified gain, channel, and voltage 

IWLICIT LOGICAL (L) 

IMPLICIT INTEGER (A-K,M-Z) 

REAL V, VS 

C 
C Check for saturation, then scale voltage 

C 

VS = 0.5 * V * 4**LGAIN 

IF(VS .GT. 4.9988) VS = 4.9988 

I F(VS .LT. -4.9988) VS = -4.9988 

VS = VS * 409.5983 

L 

C Extract sign, convert to integer 

C 

LSIGN = 1 

IF(VS.LT.0.0) LSIGN = 0 

IF(VS.LT.0.0) VS = VS + 2047.5 

IV = IFIX(2047.5 - VS) 

C 

C Form bytes 

C 

IV = IV .AND, (Xl07FFl) 

LRl = IV .AND. (XtOOFFt ) 

LR2 = (IV / 256) .OR. (LCHAN * 16) 

IF(LSIGN.EQ.1) LR2 = LR2 .OR. (X1081) 

C 

RETURN 

END 

SUBROUTINE ADBIAS(LCHAN,LGAIN,LSD,LBIAS2,VBIAS2) 

C Find bias constant required to zero A/D output 

C Ground or short inputs during zeroing 

IMPLICIT LOGICAL (L) 

IWLlClT INTEGER (A-K,M-Z) 

REAL VBIASl, VBIAS2, VBIAS3 

LOG I CAL PEEK 

DIMENSION LBiT(8) 

Tr 

COMMON /ADA/ LBIT, LSR, ABASE, ACBU, ACBL 

CI 

C Set up end points and find correspondlng readings 

C 

NIT = 0 

NITMAX = 10 

IBIASI = 0 

LBIASI = IBlASl 

CALL ADSET(LCHAN,LGAIN,LSD,LBIASl,LCBU,LCBL) 

CALL ADCNV(LCBU,LCBL,LR2,LRl) 

CALL ADVOLT(LR2,LRl,LGAlN,LCHD,VBIASl) 

IBIAS3 = 255 

LBIAS3 = IBIAS3 

CALL ADSET(LCHAN,LGAIN,LSD,LBIAS3,LCBU,LCBL) 

CALL ADCNV(LCBU,LCBL,LR2,LRl) 

CALL ADVOLT(LR2,LRl,LGAIN,LCHD,VBIAS3) 

c C Begin iteration 

C 

100 NIT = NIT + 1 

IBIAS2 = (IBIASI + IBIAS3) / 2 

LBIAS2 = IBIAS2 

CALL ADSET(LCHAN,LGAIN,LSD,LBIAS2,LCBU,LCBL) 

CALL ADCNV(LCBU,LCBL,LR2,LRl) 

CALL ADVOLT(LR2,LRl,LGAIN,LCHD,VBIASZ) 

C 

C Select new endpoint or quit 

C 

IF(NIT.EQ.NITMAX) GO TO 200 

IF(VBIAS2.GT.0) GO TO 130 

C 

IBlASl = IBIAS2 

LBlASl = LBIAS2 

VBlASl = VBIAS2 

GO TO 100 

C 

130 IBIAS3 = IBIAS2 

LBIAS3 = LBIASZ 
VBIAS3 = VBI AS2 

GO TO 100 

L 

200 RETURN 

END 

SUBROUTINE AM;AIN(LGAIN) 

C Set gain accord lng to maximum intended Input voltage 

IMPLICIT LOGICAL (L) 

DIMENSION G(4) 

DATAGhf.5, 2.0, 8.0, 32.0/ 

C 

10 WRITE(3,20) 

20 FORMAT(' Maximum input voltage = l) 

CALL ATR(VMAX,IERR) 

GMAX = 4.9988 / VMAX 

IF(GMAX.LT.G( 1)) GO TO 10 

C 

N=2 

40 IF(GMAX.LT.G(N)) GO TO 80 

N=N+l 

IF(N.EQ.5) GO TO 80 

GO TO 40 

C 

80N=N-1 

LGAlN =N - 1 

WRITE(3,90) G(N) 

90 FORMAT(' Gain = l ,F5.1) 

I, 

RETURN 

END 

SemiDisk-Interface Subroutines 

C Program # FHR-330-, 07/28/83 

C Subroutines for direct access to SemlDisk 

C WARNING: Cannot use these and SemiDisk as disk at same time! 

C 

C Counters: 

C Byte: 0 - 127 

C Sector: 0 - 15 

C Track : 0 - 255 

SUBROUTINE SDINIT(LCLR) 

Initialization for direct read-write to SemlDisk 

Set LCLR = 0 to set parameters 

= 1 to set parameters and clear SemlDisk 

IMPLICIT LOGICAL (L) 

LOGICAL INP 

COMMON /SDPARM/ LPD ATA 

Set for SemiDisk base 1/0 port address = 80H. 

LPDATA = 128 

LPBYTE = LPDATA + 1 

LPTRAK = LPDATA + 2 

LPSECT = LPDATA + 3 

Set counters to zero and set SemiDisk 

, LCTRAK, LCSECT 

I; 

LCBYTE=O 

LCTRAK=0 

LCSECTd 

CALL SDSET 
C 

C Clear Semidlsk - Write zeros untll track counter returns to zero 

C 

IF(LCLR.NE.1) RETURN 

DO 30 IT=l, 256 

DO 20 IS=l, 16 

DO 10 IB=1, 128 

CALL SDWRIT(0) 

10 CONTINUE 

20 CONTI NUE 

30 CONTINUE 

C 

C Set counters to zero and set SemiDisk 

C 

LCBYTEd 

LCTRAK4 

LCSECTd 

CALL SDSET 

C 

RETURN 

END 

SUBROUT I NE SDSET 

C Set SemlDisk to desired byte, sector, and track 

IWLlClT LOGICAL (L) 

LOGICAL INP 

COWON /SDPARM/ LFQATA,LPBYTE,LPTRAK,LPSECT,LCBME,LCTRAK,LCSECT 

CALL OUT(LPBYTE.LCBYTE) 

CALL OUT( LPSECT;LCSECT) 

CALL OUT(LPTRAK,LCTRAK) 

RETURN 

END 

SUBROUTINE SDWRIT(LCHR) 

C Write one byte (LCHR) to SemlDisk and lncrement location 

IWLICIT LOGICAL (L) 

LOGICAL INP 

CWON /SDPARM/ LPDATA,LPBYTE,LPTRAK,LPSECT,LCBYTE,LCTRAK,LCSECT 

C 

C Write byte into current location 

C 

CALL OUT(LPDATA,LCHR) 

C 

C lncrement location If necessary 

C 

CALL SDI NCR 

C 

RETURN 

END 

SUBROUTINE SDREAD(LCHR) 

C Read one byte (LCHR) to SemiDisk and lncrement location 

IWLICIT LOGICAL (L) 

LOGICAL INP 

COMMON /SDPARM/ LPDATA,LPBYTE,LPTRAK,LPSECT,LCBYTE,LCTRAK,LCSECT 

C 

C Read byte from current location 

C 

LCHR=INP(LPDATA) 

C 

C lncrement location if necessary 

I; 

CALL SDI NCR 

C 

RETURN 

END 
SUBROUTINE SDINCR 

C Increment SemlDlsk location pointers and set SemiDlsk if necessary 

IMPLICIT LOGICAL (L) 

LOG1 CAL INP 

COMMON /SDPARM/ LPDATA,LPBYTE,LPTRAK,LPSECT,LCBYTE,LCTRAK,LCSECT 

C 

C Next byte 

C 

LCBYTE=LCBYTE + 1 

iF(LCBYTE.GT.0) RETURN 

C 

C Next sector 

C 

LCBYTE 4 

LCSECT=LCSECT + 1 

IF(LCSECT.LE.15) GO TO 10 

I; 

C Next track n 

C, 

LCSECT4 

LCTRAK=LCTRAK+l 

C 

C Set SemlDisk for new byte, sector, and track 

C 

10 CALL SDSET 

C 

RETURN 

END 

- A/D Operation 

C Program # FHR-345-, 09/02/83 

C Four-channel A/D converter with output to SemiDisk 

C 

IMPLICIT LOGICAL (L) 

IMPLICIT INTEGER (A-K,M-Z) 

LOGICAL PEEK 

REAL V(4), VMAX, GAIN, GMAX, FREQ 

DIMENSION LBITCB), LCBU(41, LCBL(41, LBIAS(4), LCHAN(41, 

> LID(128), LR1(6), LRZ(4) 

COMMON /ADA/ LBIT, LSR, ABASE, ACBU, ACBL 

COMMON /SDPARM/ LPDATA,LPBYTE,LPTRAK,LPSECT,LCBYTE,LCTRAK,LCSECT 

C 

STOP 4 

WRITE(3,10) 

10 FORMAT(' Four-channel A/D converter with SemiDisk output1, 

> I, Program # FHR-345-I) 

CALL ADVERS(3) 

C 

C Get A/D parameters 

C 

20 WRITE(3,30) 

30 FORMAT(/,' Physical channel numbers (0 - 7) : I) 

DO 40 N = 1, 4 

CALL ATI(1,IERR) 

IF((IERR.NE.0) .OR. (I.LT.0) .OR. (I.GT.7)) GO TO 20 

LCHANCN) = I 

40 CONTINUE 

C 

CALL ADGAIN(LGA1N) 

C 

LSD = 0 

LSR = 0 

I; 

C Offset bias correction 

C 

CALL OBC(LCHAN,LGAIN,LSD,LBIAS) 

C 

C Form control bytes and addresses, set A/D converter 
I; 

DO 110 N=l, 4 

CALL ADSET(LCHAN(N),LGAiN,LSD,LBIAS(N),LCBU(N),LCBL(N)) 

110 CONTINUE 

C 

C Get run parameters 

C 

WRITE(3,120) 

120 FORMAT(/,' Number of 100-sample sets = I) 

CALL AT1 (NSR, IERR) 

WRITE(3,130) 

130 FORMAT(' Number of 3.3-ms intervals per sample = ) 

CALL ATi(NINT,IERR) 

NEXTRA = NI NT - 2 

FREQ = 300.0 / NlNT 

WRiTE(3,135) FREQ 

135 FORMAT(' Sampling frequency = I,F10.3) 

WRITE(3,140) 

140 FORMAT(/,* Recording ID = I) 

READ(3,150) (LID(N), N=1,64) 

150 FORMAT(64Al) 

C 

C Initialize semidisk 

C 

CALL SDiNiT(1) 

C 

C Prepare identification data and write to SemlDisk 

C 

C LID(1) - LID(64) are ASCII message from operator 

LiD(65) = LGAlN 

CALL ILC(NSR,LID(66),LiD(67)) 

CALL ILC(NINT,LiD(68),LID(69)) 

DO 230 N = 70, 128 

LID(N) = 0 

230 CONTINUE 

DO 240 N = 1, 128 

CALL SDWRIT(LID(N)) 

240 CONTINUE 

C 

C Advise user to beg in 

C 

CALL OUTPUT(3,7) 

WRITE(3,300) 

300 FORMAT(//,I NOW READY - Press RETURN to start1 ) 

READ(3,310) 

310 FORMATO 

NSAMP = 0 

C 

C Begln taking measurements - reset clock 

C 

CALL OUTPUT(3,7) 

CALL CMSET 

DO 600 NS = 1, NSR 

DO 500 N1 = 1, 100 

CALL INTRPT(LINT1 

IF(LiNT.EQ.27) GO TO 900 

C 

C Digitize, then wait for end of first 3.3-ms Interval c. 

C, 

DO 450 N=l, 4 

CALL ADCNV(LCBU(N),LCBL(N),LR2(N),LRl(N)) 

450 CONTINUE 

CALL TIMER 

L 

C Pack, then wait for end of second 3.3-ms interval 

C 

CALL PACK(LR2,LRl) 

CALL TIMER 

C 

C Transfer to SemiDisk, then wait for end of additional 3.3-n 

DO480 N = 1. 6 
CALL SDWRIT(LRl(N)) 

480 CONTINUE 

CALL WAIT(NEXTRA1 

500 CONTINUE 

600 CONTINUE 

CALL OUTPUT(3,7) 

C 

900 CONTINUE 

STOP 

END 

SUBROUTINE OBC(LCHAN,LGAIN,LSD,LBIAS) 

C Offset bias correction 

IMPLICIT LOGICAL (L) 

IWLICIT INTEGER (A-K,M-Z) 

REAL V(4) 

DIMENSION LBIT(8), LC8U(4), LCBL(41, LBIAS(4), LCHAN(4) 

COmON /ADA/ LBIT, LSR, ABASE, ACBU, ACBL 

C 

100 WRITE(3,110) 

110 FORMAT(/,' Method for offset-bias correction I, 

> '(1-assume, 2-measure, 3-disk) : '1 

CALL ATI(METHOD,IERR) 

lF((METHOD.LT.l).OR.(METHOD.GT.3)) GO TO 100 

IF(METHOD.EQ.2) GO TO 130 

IF(METHOD.EQ.3) GO TO 160 

C 

C Assumed bias corrections 

DO 120 N=l, 4 

LBIASCN) = 127 

120 CONTINUE 

GO TO 180 

C 

C Measured bias corrections 

130 CONTINUE 

DO 140 N=l, 4 

CALL ADBIAS(LCHAN(N),LGAIN,LSD,LBIAS(N),V(N)) 

140 CONTINUE 

WRITE(3,150) (V(N), N=l, 4) 

150 FORMAT(' Blas voltages = ',4F9.4) 

GO TO 180 

C 

C Blas corrections from disk file 

160 CONTINUE 

CALL FOPENl(6,IERR) 

READ(6,170) (LBIAS(N), N=l, 4) 

170 FORMAT(415) 

ENDFILE 6 

C 

C Display bias corrections 

180 WRITE(3,198) (LBIASCN), N=l, 4) 

190 FORMAT(' Bias bytes = ',419) 

C 

RETURN 

END 

SUBROUTINE ILC(I,Ll,L2) 

C Converts one integer to two bytes 

C L1 is high-order byte 

C Memory assignment for Microsoft FORTRAN-88 

IWLICIT LOGICAL (L) 

LOGICAL PEEK 

MI=MLOC$(I) 

L24'EEK(Ml) 

Ml=Ml+l 

Ll-PEEK(M1) 

RETURN 

END 
SUBROUTINE PACK(LR2,LRl) 

C Pack 8 bytes into 6 bytes 

iMPLlClT LOGICAL (L) 

IMPLICIT INTEGER (A-K,M-Z) 

DIMENSION LBIT(81, LRl(1). LR2(1) 

COMMON /ADA/ LBIT, LSR, ABASE, ACBU, ACBL 

C 

LRl(5) = LR2( 1) .AND. (X1OF1 ) 

LRl(6) = LR2(3) .AND. (X1OF1 

C 

DO101 =1,4 

J = 1 +4 

LTEST = LR2(2) .AND. LBITCI) 

IF(LTEST.NE.0) LRl(5) = LRl(5) .OR. LBIT(J) 

LTEST = LR2(4) .AND. LBIT(1 

I F(LTEST.NE.0) LRl(6) = LRl(6) .OR. LBIT( J) 

10 CONTINUE 

C 

RETURN 

END 

SUBROUTINE CLKSET 

C Reset real-time clock In North-Star Horizon 

CALL OUT(6,X150') 

RETURN 

EM) 

SUBROUTINE WAITCN) 

C Walt N time intervals 

DO 10 1=1, N 

CALL TIMER 

10 CONTINUE 

RETURN 

END 

SUBROUTINE TIMER 

C 3.3-millisecond timer for North-Star Horizon 

C Note: Jumper real-time clock for 3.3 ms 

C Resets clock at end of interval 

C 

IMPLICIT LOGICAL (L) 

LOG1 CAL I NP 

C 

C Wait for clock flag 

10 L2=1NP(6) .AND. 4 

IF(L2.EQ.0) GO TO 10 

C 

C Reset clock to start new interval 

CALL OUT(6,X' 50' ) 

C 

RETURN 

EM) 

SemiDi sk - to Floppy-Di sk Transfer 

C Program # FHR-346A, 09/02/83 

C SemlDlsk to file transfer 

C 

IMPLICIT LOGICAL (L) 

IMPLiClT INTEGER (A-K,M-Z) 

LOG1 CAL PEEK 

REAL V(41, GAIN 

DIMENSION LBIT(8), LCBU(4), LCBL(4), LBIAS(4), LCHAN(41, 

> LID(128), LR1(6), LR2(4) 

COWON /ADA/ LBIT, LSR, ABASE, ACBU, ACBL 
COMMON /SDPARM/ LPDATA,LPBYTE,LPTRAK,LPSECT,LCBYTE,LCTRAK,LCSECl 

C 

WRiTE(3,lB) 

10 FORMAT(' SemlDisk to f i le transfer,', 

> 1 Program # FHR-346A1 ) 

CALL ADVERS(3) 

WRITE(3,20) 

20 FORMAT( 

CALL FOPENO(7,IERR) 

CALL SDINIT(~) 

WRITE(3,30) 

30 FORMAT(' Number of 10&samp le sets = ) 

CALL ATI(NSR.IERR) 

C 

C Identification data 

C 

WEC = 1 

NLOC = 1 

LEND = 0 

W 110 N = 1, 128 

CALL SDREAD(LCHR) 

CALL PUTCHR(NREC,NLOC,LCHR,LEND) 

118 CONTINUE 

C 

C Measurements 

C 

W 600 NS = 1, NSR 

DO 500 N1 = 1, 100 

C 

C Transfer six bytes 

C 

DO 450 N=l, 6 

CALL SDREAD(LCHR1 

CALL PUTCHR(NREC,NLOC,LCHR,LEND) 

450 CONTINUE 

C 

500 CONTINUE 

600 CONTINUE 

C 

C Zero rest of file and close 

C 

LEND = 1 

CALL PUTCHRC NREC, NLOC, LCHR, LEND) 

ENDFILE 7 

STOP 

END 

SUBROUTINE PUTCHR(NREC,NLOC,LCHR,LEND) 

C Unidirectional single-character write to disk file 

IMPLICIT LOGICAL (L) 

DIMENSION LBUF(128) 

C 

IF(LEND.NE.0) GO TO 10 

C 

LBUF(NL0C) = LCHR 

NLOC = NLOC + 1 

I F( NLOC. LE. 128 RETURN 

C 

WRITE(7, REC=NREC) LBUF 

NLOC = 1 

NREC = NREC + 1 

C; 

10 CONTINUE 

W 20 N = NLOC, 128 

LBUF(N) = 0 

20 CONTINUE 

WRI TE(7, REC-NREC) LBUF 

C 

RETURN 

END 
Data Decoder 

C Program # FHR-3470, 09/02/83 

C Decoder for A/D converter data 

C 

IMPLICIT LOGICAL (L) 

IMPLICIT INTEGER (A-K,M-Z) 

LOGICAL PEEK 

REAL V(4), GAIN, FREQ 

DIMENSION LBIT(8), LCBU(4), LCBL(41, LBIAS(41, LCHAN(4); 

> LID(128), LRl(61, LR2(4) 

COWON /ADA/ LBIT, LSR, ABASE, ACBU, ACBL 

C 

WRITE(3,10) 

10 FORMAT(' Decoder for A/D converter data', 

> 1, Program # FHR-3470') 

CALL ADV ERS ( 3 

WRITE(3,20) 

20 FORMAT( 

CALL FOPENI(6,IERR) 

C 

C Get identification data 

C 

NIEC = 1 

NLOC = 1 

DO 110 N = 1, 128 

CALL GETCHR(NREC,NLOC,LID(N),LEM)) 

IF(LEND.NE.0) GO TO 900 

110 CONTINUE 

C 

WRITE(3,120) (LIDCN), N=l, 64) 

120 FORMAT(/,' Data identification : ',64A1) 

C 

LGAlN = LID(65) 

GAIN = 0.0 

I F(LGAIN.EQ.0) GAIN = 0.5 

IF(LGAIN.EQ.1) GAIN = 2.0 

IF(LGAIN.EQ.2) GAlN = 8.0 

IF(LGAIN.EQ.3) GAlN = 32.0 

WRITE(3,130) GAIN, LGAlN 

130 FORMAT(' Gain, gain code = ',F5.2,15) 

C 

CALL LIC(LID(66),LID(67),NSR) 

WRITE(3,140) NSR 

140 FORMAT(' Number of 100-sample sets = I, 18) 

C 

CALL LIC(LID(68),LID(69),NINT) 

WRITE(3,150) NlNT 

150 FORMAT(' Number of 3.3-ms intervals / sample = ',I81 

IF(NINT.EQ.0) GO TO 170 

FREQ = 300.0 / NlNT 

WRITE(3,160) FREQ 

160 FORMAT(' Sampling frequency = ',F10.3) 

170 CONTINUE 

C 

C Measurements 

L 

WRITE(3,300) 

300 FORMAT(/,' Nsetl ,4XJ1N1' ,6XJ1V1' ,7XJ1V2' ,7XJ1V3', 

> 7X,'V4',/) 

DO 600 NS = 1, NSR 

DO 500 N1 = 1, 100 

C 

C Get six bytes from disk 

C 

DO 450 N=l, 6 

CALL GETCHR(NREC,NLOC,LRl(N),LEND) 

IF(LEND.NE.0) GO TO 900 
450 CONTINUE 

C 

C Unpack and print 

C 

CALL UNPACK (LR2, LR 1 

DO 460 N=l, 4 

CALL ADVOLT(LR2(N),LRl(N),LGAIN,LCHD,V(N)) 

460 CONTINUE 

WRITE(3,470) NS, N1, (V(N), N=1, 4) 

470 FORMAT(' ',216,4F9.3) 

C 

500 CONTINUE 

600 CONTINUE 

L 

900 WRITE(3,910) 

910 FORMAT(/,' EOF encountered1 ,/) 

ENDFILE 6 

STOP 

END 

SUEROUTINE UNPACK(LR2,LRl) 

C Unpack 6 bytes into 8 bytes 

IMPLICIT LOGICAL (L) 

IMPLICIT INTEGER (A-K,M-Z) 

DIMENSION LBIT(8), LR1(1), LR2(1) 

COMMON /ADA/ LBIT, LSR, ABASE, ACBU, ACBL 

C 

LR2(1) = LRl(5) .AND. (X1OFt) 

LRZ(3) = LRl(6) .AND. (X1OF1) 

C 

LRZ(2) = 0 

LRZ(4) = 0 

DO 10 l = 1, 4 

J = 1 +4 

LTEST = LRl(5) .AND. LBIT(J) 

IF(LTEST.NE.0) LR2(2) = LR2(2) .OR. LBIT(I) 

LTEST = LRl(6) .AND. LBIT(J) 

IF(LTEST.NE.0) LR2(4) = LRZ(4) .OR. LBIT(I) 

10 CONTINUE 

C 

RETURN 

END 

SUBROUTINE LIC(Ll,L2,1) 

C Converts two bytes to one integer 

C L1 is high-order byte 

C Memory assignment for Microsoft FORTRAN-80 

IMPLICIT LOGICAL (L) 

Mt=MLOCS( I) 

CALL POKE(Ml,L2) 

MI =Ml+l 

CALL POKE(M1,Ll) 

RETURN 

END 

SUBROUTINE GETCHR(NREC,NLOC,LCHR,LEND) 

C Unldirectional single-character read from disk file 

IMPLICIT LOGICAL (L) 

DIMENSION LBUF(128) 

C 

IF(NLOC.EQ. 1) READ(6, REC=NREC, END=l0) LBUF 

C 

LCHR = LBUF(NLOC) 

LEND = 0 

NLOC = NLOC + 1 

IF(NLOC.LE.128) RETURN 

C 

NLOC = 1 
NREC = NREC + 1 

RETURN 

C 

10 LCHR = 0 

LEND = 1 

RETURN 

END 
APPENDIX B. COMMUNICATION PROGRAM 

Program # FHR-350B, 09/27/83 

Transmlssion of binary files to HP-1000 via RS-232 

HP binary-block protocol : 

Send : 

STX (start block), SUB (EOF), or EOT (all done) 

128 8-bit binary data bytes 

Check byte 

Answer: 

ACK or NAK 

Repeat up to two times if NAK received 

Intersperse HP protocol for HP terminals 

I WL I C I T l NTEGER ( A-K ,M-Z) 

IMPLICIT LOGICAL (L) 

DIMENSION LBUF(128) 

COWON /lo/ LUNC,LUNI,LUNO,LPAR,LISP,LISM,LIDP,LIDM, 

> LOSP,LOSM, LODP,LODM, LBUF 

DATA LREAD /-I/ 

DATA LUN I , LUNO /0,0/ 

DATA LISP,LISH,LIDP,LIDM /0,0,0,0/ 

DATA LOSP,LOSM,LODP,LODM /0,0,0,0/ 

DATA NEMAX /3/ 

CALL GETLUN(LUNC,LUND) 

WRITE(LUNC,10) 

10 FORMAT(' HP binary-block transmitter, Program # FHR-350B',/, 

> Copyright (C) 1983 by Green Mountain Radio Research' ,/) 

Get I /O parameters 

I; 

IF(LREAD.NE.-1) GO TO 290 

LREAD=l 

C lnput channel 

100 WRITE(LUNC,lIB) 

110 FORMAT(' lnput channel :I 

CALL GETPRT(LUNC,LISP,LISM,LIDP,LIDM) 

C Output channel 

200 WRITE(LUNC,210) 

210 FORMAT(' Output channel : I) 

CALL GETPRT(LUNC,LOSP,LOSM,LODP,LODM) 

C 

C Clear 1/0 port 

C 

290 WRITE(LUNC,310) 

31 0 FORMAT(/,' Press RETURN when ready1 ) 

CALL OUTPUT(LUNC,7) 

READ(LUNC,320) 

320 FORMAT ( ) 

CALL CLRlO 

C 

C Begin transfer - open disk file 

C 

300 CALL OUTPUT(LUNC,7) 

CALL FOPENI (LUND,iERR) 

IREC=l 

C 

C Get disk record and transmit 

C 

400 CONTINUE 

C 

C Wait for HP to send DCl 

C 

CALL GETDCl 

I; 

C Read record 

C 

READ(LUND, REC=IREC, END=l000) LBUF 

NERROR = 0 
410 WRITE(3,420) IREC 

420 FORMAT(' # ',I61 

CALL SNDREC 

CALL HP1 

CALL GETACK(LERR0R) 

CALL HP2 

IF(LERROR.EQ.0) GO TO 490 

CALL OUTPUT(LUNC,7) 

NERROR = NERROR + 1 

IF(NERROR.LT.NEMAX) GO TO 410 

GO TO 1000 

C 

490 IREC=IREC+l 

GO TO 400 

C 

C Done - close file and see whether to continue 

C 

1000 ENDFILE LUND 

WRITE(LUNC, 1010) 

1010 FORMAT(/,' Continue ? '1 

READ(3.1020) LANS 

1020 FORMAT ( A 1 ) 

LANS=LUPPER(LANS) 

I F(LANS.NE.89) GO TO 1100 

C 

C End of file but more files - send SUB, then set of null 

CALL XCZ 

GO TO 300 

C 

C No more f i les - send EOT and set of nu1 Is 

1100 CALL XEOT 

CALL HALT 

STOP 

END 

SUBROUTINE GETLUN(LUNC,LUND) 

C Get console and dl sk LUN 

IMPLICIT LOGICAL (L) 

LUNC=3 

LUND =ti 

RETURN 

END 

SUBROUTINE GETPRT(LUNC,LISP,LISM,LIDP,LIDM) 

C Get 1/0 port parameters 

IWLlClT LOGICAL (L) 

10 FORMAT(2Al) 

C 

20 WRITE(LUNC,3B) 

30 FORMAT(' Status port (2-digit hex) = 

READ(LUNC,l0) LH2,LHl 

LH4 =' 0' 

LH3 =I 0' 

CALL HIC(LH4,LH3,LH2,LHl,I,LERR) 

IF(LERR.NE.0) GO TO 20 

LISP=\ 

L 

40 WRITE(LUNC,50) 

50 FORMAT(' Status mask (2-digit hex) = '1 

READ(LUNC,lg) LH2,LHl 

LH4 =I 0' 

LH3 =' 0' 

CALL HIC(LH4,LH3,LH2,LHl,I,LERR) 

IF(LERR.NE.0) GO TO 40 

LISM=I 

C 

60 WRITE(LUNC,70) 

70 FORMAT(' Data port (2-digit hex) = ') 

READ(LUNC,l@) LH2,LHl 
LH4 =I 0' 

LH3=' 0' 

CALL HIC(LH4,LH3,LH2,LHl,IJLERR) 

IF(LERR.NE.0) GO TO 60 

LIDP=I 

C 

80 WRITE( LUNC,90) 

90 FORMAT(' Data mask (2-dig it hex) = 

READ(LUNC,lB) LH2,LHl 

LH4 =' 0' 

LH3=' 0' 

CALL HIC(LH4,LH3,LH2,LHl, LLERR) 

IF(LERR.NE.0) GO TO 80 

LIDM=I 

C 

RETURN 

END 

SUBROUTINE SNDREC 

IMPLICIT INTEGER (A-K,M-Z) 

IMPLICIT LOGICAL (L) 

DIMENSION LBUF(128) 

COWION /lo/ LUNC,LUNI,LUNO,LPAR,LISP,LISM,LIDP,LIDM, 

> LOSP,LOSM,LODP,LODM,LBUF 

DATA LSTX / 2 / L 

C Send start-transmlssion character 

C 

CALL OUTPUP(LSTX) 

C 

C Send 128 data bytes 

C 

DC) 10$ 1=1, 128 

CALL OUTPUP(LBUF(1)) 

100 CONTINUE 

C 

C Compute and send parity check 

C 

LPAR = 0 

DO 200 1 = 1, 128 

LPAR = LPAR .XOR. LBUF( I) 

208 CONTINUE 

CALL OUTPUP( LPAR) 

C 

RETURN 

END 

SUBROUTINE XCZ 

C Transmit control-Z to signal end of transmission 

IMPLICIT INTEGER (A-K,M-Z) 

IMPLICIT LOGICAL (L) 

DIMENSION LBUF(128) 

COMMON /lo/ LUNC,LUNI ,LUNO,LPAR,LISP,LISM,LIDP,LIDM, 

> LOSP,LOSM,LODP,LODM,LBUF 

C 

CALL OUTPUP ( 26) 

DO 10 I = 1, 129 

CALL OUTPUP(0) 

10 CONTINUE 

C 

CALL HPl 

C 

RETURN 

END 

SUBROUTINE XEOT 

C Transmit control-D to signal all fin,ished 

IMPLICIT INTEGER (A-K,M-Z) 
IMPLICIT LOGICAL (L) 

DIMENSION LBUF(128) 

COWON /lo/ LUNC,LUNI,LUNO,LPAR,LISP,LISM,LIDP,LIDM, 

> LOSP,LOSM,LODP,LODM,LBUF 

C 

CALL OUTPUP(4) 

DO 10 1 = 1, 129 

CALL OUTPUP (0 ) 

10 CONTINUE 

CALL HPl 

C 

RETURN 

END 

SUBROUTINE GETACK(LERR0R) 

C Get ac know l edgement f rom rece i ver 

IWLICIT INTEGER (A-K,M-Z) 

IMPLICIT LOGICAL (L) 

DIMENSION LBUF(128) 

COMMON /lo/ LUNC,LUNI,LUNO,LPAR,LISP,LISM,LIDP,LIDM, 

> LOSP,LOSM,LODP,LODM,LBUF 

DATA LACK /6/ 

C 

CALL INPUP(LRCV) 

LERROR = 0 

IF(LRCV.EQ.LACK) RETURN 

LERROR = 1 

C 

RETURN 

END 

SUBROUTINE INPUP(LCHR) 

C Direct input from computer 1/0 port 

I WL I C IT l NTEGER ( A-K ,M-Z) 

IWLICIT LOGICAL (L) 

DIMENSION LBUF(128) 

COWON /lo/ LUNC,LUNI,LUNO,LPAR,LISP,LISM,LIDP,LIDM, 

> LOSP,LOSM,LODP,LODM,LBUF 

C 

C Check status, wait for something to happen 

100 LCHR=INP(LISP) 

LCHR=LCHR. AND .L I SM 

IF(LCHR.EQ.0) GO TO 180 

C 

C Something came in; get it 

LCHR=INP(LIDP) 

C 

RETURN 

END 

SUBROUTINE OUTPUP(LCHR) 

C Direct output to computer 1/0 port of 8080 computer 

IMPLICIT INTEGER (A-K,M-Z) 

IMPLICIT LOGICAL (L) 

DIMENSION LBUF(1281, LBIT(8) 

COMMON /lo/ LUNC,LUNI,LUNO,LPAR,LISP,LISM,LIDP,LIDM, 

> LOSP,LOSM,LODP,LODM,LBUF 

DATALBlT/l, 2, 4, 8, 16, 32, 64, X180'/ 

C 

C Check status - wait untl l OK 

C 

10 LSTAT=INP(LOSP) .AND. LOSM 

IF(LSTAT.EQ.0) GO TO 10 

C 

C Send character 

L 

CALL OUT(LODP,LCHR) 
RETURN 

END 

SUBROUTINE CLRlO 

Clear/initialize 1/0 ports 

IWLICIT INTEGER (A-K,M-Z) 

IWLlClT LOGICAL (L) 

DIMENSION LBUF(128) 

COIJMON /lo/ LUNC,LUNI,LUNO,LPAR,LISP,LISM,LIDP,LIDM, > LOSP,LOSM,LODP,LODM,LBUF 

LCHR=INP(LIDP) 

RETURN 

END 

SUBROUTINE GElDC1 

IWLICIT INTEGER (A-K,M-Z) 

IMPLICIT LOGICAL (L) 

DIMENSION LBUF(128) 

COMMON /lo/ LUNC,LUNI,LUNO,LPAil,LISP,LISM,LIDP,LIDM, 

> LOSP, LOSM, LODP, LODM, LBUF 

DATA LDCl /17/ 

CALL INPUP(LCHR) 

IF(LCHR.EQ.LDC1) RETURN 

WRITE(LUNC,10) 

FORMAT(' *** FIRST CHAR NOT DC1 ***I 

CALL HALT 

RETURN 

END 

SUBROUTINE HP1 

HP-1000 transmlssion terminator 

IMPLICIT INTEGER (A-K,M-Z) 

IMPLICIT LOGICAL (L) 

DIMENSION LBUF(128) 

COMMON /lo/ LUNC,LUNI,LUNO,LPAR,LISP,LISM,LIDP,LIDM, 

> LOSP,LOSM,LODP,LODM,LBUF 

Receive ENQ and send ACK (twice) 

CALL ENQACK 

CALL ENQACK 

RETURN 

END 

SUBROUTINE HP2 

HP-1000 reception terminator 

IMPLICIT INTEGER (A-K,M-Z) 

IMPLICIT LOGICAL (L) 

DIMENSION LBUF(128) 

COMMON /lo/ LUNC,LUNI,LUNO,LPAR,LISP,LISM,LIDP,LIDM, 

> LOSP,LOSM,LODP,LODM,LBUF 

Get CR LF 

CALL INPUP(LCHR) 

CALL INPUP(LCHR) 

Receive ENQ and send ACK (twice) 

CALL ENQACK 

CALL ENQACK 

CALL ENQACK 
L 

RETURN 

END 

SUBROUTINE ENQACK 

C Recel ve ENQ, send ACK 

IMPLICIT INTEGER (A-K,M-Z) 

IMPLICIT LOGICAL (L) 

DIMENSION LBUF(128) 

CWON /lo/ LUNC,LUNI,LUNO,LP~,LISP,LISM,LIDP,LIDM, 

> LOSP,LOSM,LODP,LODM,LBUF 

DATA LENQ, LACK /5, 6/ 

C 

CALL INPUP(LCHR1 

IF(LCHR.NE.LENQ) GO TO 100 

CALL OUTPUP(LACK) 

RETURN 

C 

C Error 

100 WRITE(LUNC, 110) LCHR 

110 FORMAT(' *** ENQ NOT RECEIVED ***',I51 

CALL HALT 

RETURN 

EM) 
APPENDIX C. PROCESSING PROGRAM 

Digitized-Noise Simulator 

C Program # FHR-352A, 18/28/83 

C ELF-noise generator with disk-fi le output 

C 

IMPLICIT LOGICAL (L) 

DIMENSION P(4),GAMMAL(4) ,A(4) ,CHX(4) ,CHY(4) ,PSI (4),R(4), 

> SIGX(4), HXOLD( 100), WOLD( 100), FW( 100), 

> ZW(3), YW(31, ZR(3), Z1(3), YR(3), Y1(3), 

> LR1(6), LR2(4), LID(128) 

CWON /BLOCKl/ IR~,~R~,P,A,CHX,CHY,PSI,R,SIGX,PSIMIN,PS~D 

COMON /BLOCK2/ NFI LT, FW 

COmON /PCPARM/ NREC, NBYTE 

DATA PSI (1) ,PSI (2),PSI (3) ,PSI (4) /45.8,135.8,-70*6n8.8/ 

DATA LGAIN, LCHAN /1,0/ 

10 FORMAT(' ELF-noise generator with disk-file output, l, 

> 1 Program # FHR-352A1 ,/) 

C 

C Fixed parameters and constants 

C 

DR4.017453293 

TWOP1=6.2831852 

SIGMAG=l.B 

KMCD=l00 

FSAMP=l00.0 

CF I LT=l .0 

FSIG=10.0 

l LS=l 

I LY =3 

I LY2=I LY"2 

N3MS=3 

NREC=l 

NBYTE=l 

C 

C Announce name and get other parameters 

C 

WRITE(3,10) 

CALL ADVERS(3) 

CALL ADSET(LCHAN,LGAIN,LSD,LBIAS,L~U,L~L) 

90 WRITE(3,lBB) 

100 FORMAT(/, 1 INPUT: ~oi se-source parameters1 ,/, 

> 1 #',%,lP, Gamna, A, Pslmin, Pslmax (degrees)l) 

DO 120 1=1, 4 

WRITE(3.110) I 

110 FORMAT(' ',12,%) 

CALL ATR(P(I),IERR) 

IF( IERR.EQ.27) GO TO 90 

CALL ATR(GAMMAL( I), IERR) 

IF(IERR.EQ.27) GO TO 90 

CALL ATR(A(I),IERR) 

IF(IERR.EQ.27) GO TO 90 

120 CONTINUE 

CALL ATR(PSIMIN,IERR) 

IF(IERR.EQ.27) GO TO 90 

CALL ATR(PSIMAX,IERR) 

IF( IERR.EQ.27) GO TO 90 

PSID=PSIMAX-PSIMIN 

C 

C Type of reference 

WRITE(3,130) 

130 FORMAT(/,' Reference (E or M) ? '1 

READ(3,135) LREF 

135 FORMAT(A1) 

LREF=LUPPER(LREF) 

C 

C Local noise level 

WRI TE(3,140) 
140 FORMAT(' Plocal = '1 

CALL ATR(PLOCAL,IERR) 

SIGLCL=SQRT(PLOCAL) 

C 

C Pr imary-antenna channel 

WRITE(3,155) 

155 FORMAT(' Primary-antenna axis (X, Y, or Z) = 1) 

READ(3.156) LPCH 

LPCH = LUPPER(LPCH1 

156 FORMAT(A1) 

IPCH = LPCH - 87 

C 

C Other parameters 

WRITE(3,200) 

200 FORMAT (/, ' Nsec = ' 

CALL ATI(NSEC,IERR) 

WRITE(3.220) 

220 FORMAT(' Nfilter, Yzx, Yzy, Yzz = 

CALL ATI(NFILT,IERR) 

CALL ATR(YZX,IERR) 

CALL ATR(YZY,IERR) 

CALL ATR(YZZ, I ERR) 

WRITE(3,240) 

240 FORMAT(' Seeds = ' 

CALL ATI(ISEED1,IERR) 

CALL ATI(ISEED2,IERR) 

IF((ISEEDl.NE.0).AND.(ISEED2.NE.0)) GO TO 250 

CALL SEED(ISEED1) 

CALL SEED(ISEED2) 

250 lRl=lSEEDl 

IR2=1 SEED2 

WRITE(3,260) 

260 FORMAT(' Output voltages ? ) 

READ(3,156) LOV 

LOV = LUPPER(L0V) 

C 

C Initial lze output fl le and write ID to it 

C 

CALL FWENO(6,IERR) 

CALL IDPREP(LGAIN,NSEC,N3MS,LID) 

DO 270 1=1, 128 

CALL PUTCHR(LID(1)) 

270 CONTINUE 

C 

C Wrlte parameters and headlngs 

C 

WRlTE(3,300) 

300 FORMAT(23X,'Nol se-Source Data1 ,/,' 1 ' ,5X,'P' ,9X,'Gamna',9X, 

> 'A1 ,6X,1Ps1/Psimin1 ,2X,1Psimax1 ,/) 

DO 340 1 =1,4 

IF(I.EQ.4) GO TO 320 

WRlTE(3,310) I, P(I), GAMMALtI), A(I), PSI(I) 

310 FORMAT(' #', 12,3G12.4,F8.1,2X,F8.1) 

GO TO 330 

320 WRITE(3,310) I, P(I), GAMMAL(I), A(I), PSIMIN, PSIMAX 

330 CONTINUE 

PSI (I )PSI (I )*DR 

CHX(I)-SIN(PSI(I)) 

CHY(I).COS(PSi(I)) 

340 CONTINUE 

PSIMINPSIMIN*DR 

PSI MAXPSI MAXVR 

PSI DPS I D"DR 

PT.B.0 

DO 345 1 =1,4 

PTPT+P( I ) 

345 CONTINUE 

PT4 PT*4 .0 

PT9=PT*9.0 

WRITE(3,346) PT 

346 FORMAT(' T 1,612.4) 
C 

C Output headlngs 

C 

IF(LOV.EQ.89) WRITE(3,405) 

405 FORMAT(/,ZX,fNsec' ,W,'N1Omsl ,6X,'Vxrf ,6X,"4yr1, 

> 6X,1Vzr' ,6X,'Vp1 ,/) 

C 

C Initialize parameters 

C 

ISEEDl=IRl 

I SEED2=IR2 

CALL RANDU(IR2,lRl) 

CALL RANDU(IR2,lRl) 

DO 410 1=1,4 

IF(P(I).GT.0.0) CALL ELFPRM(P(I),GAMMAL(I),A(I),SIGX(I),R(I)) 

410 CONTINUE 

DO 440 J=l ,NFILT 

HXOLD(J)=8.8 

HYOLD(J)=B.0 

FW(J)=PR(CFILT,J) 

440 CONTINUE 

C 

C Begln seconds loop 

c 

DO 7000 KN=l, NSEC 

WRITE(3,490) KN 

490 FORMAT(' T = '.I51 

L 

C Begin 10-ms loop P 

b 

DO 1000 KM = 1, KMCU 

C 

CALL iNTRPT(LINT) 

IF(LINT.EQ.27) GO TO 9999 

L 

C Beg In nolse generation 

C 

C Generate incident nolse fields 

CALL ELFINC(ZW(l),ZW(2),EZ) 

C 

C Horizontal electrlc flelds 

YW(l)=FILTER(ZW(2),HYOLD) 

YW(2)=FILTER(ZW(l),HXOLD) 

C 

C Vertlcai magnetlc field 

ZW(3)=YZX*ZW(l) + YZYUZW(2) + YZZ*EZ 

C 

C Vertical electrlc field 

YW(3)=EZ 

C 

C Magnetic reference antenna? 

IF(LREF.EQ.69) GO TO 510 

DO 500 I=1,3 

YW( I )=ZW( I) 

5Q0 CONTINUE 

510 CONTINUE 

C 

C Generate and add local noise to reference-antenna output 

C 

DO 520 1=1, 3 

CALL GAUSS(IR2,1Rl,0.0,SIGLCL,GN) 

YW(I )=YW( I) + GN 

520 CONTINUE 

C 

C Select primary channel 

ZW( 1) = ZW( IPCH) 

C 

C Quantize 

DO 530 1 =I, 3 

CALL ADSIM(LGAIN,LCliAN,YW(I),LR2(1). LRl(I)) 

530 CONTl NUE 
CALL ADSlM(LGAIN,LCHAN,ZW(l),LR2(4), LRl(4)) 

DO 588 1 = 1, 3 

CALL ADVOLT(LR2(1),LRl(I),LGAIN,LCHANO,YW(I)) 

580 CONTl NUE 

CALL ADVOLT(LR2(4),LR1(4),LGAIN,LCHANO,ZW(l)) 

C 

C Pack and output 

IF(LOV.EQ.89) WRITE(3,780) KN, KM, (YW(I), l=1,3), ZW(1) 

780 FORMAT(' ,216,4F9.3) 

CALL PACK(LR2,LRl) 

DO 718 1=1, 6 

CALL PUTCHR(LRl(1)) 

710 CONTINUE 

C 

C End 10-ms loop 

C 

900 CONTINUE 

1888 CONTINUE 

C 

6990 CONTINUE 

7880 CONTINUE 

C 

9808 CONTINUE 

C 

C 

9999 CALL PCEND 

ENDFILE 6 

STOP 

END 

SUBROUTINE ELFINC(HX,HY,EZ) 

C INCIDENT ELF FIELDS 

DIMENSION P(4),A(4),CHX(4),CHY(4),PS1(4),R(4), 

> SIGX(4) 

COMMON /BLOCKl/ IR2,IRl,P,A,CHX,CHY,PSI,R,SIGX,PSIMIN,PSID 

100 Hx4.0 

HY4.8 

EZ4.0 

DO 318 1 =1,4 

IF(P(I).LE.0.0) GO TO 300 

IF( I .LT.4) GO TO 200 

C GET RANDOM AZIMUTH ANGLE FOR SOURCE 14 

CALL RANDU(IR2,lRl) 

I0=IABS(IRl) 

PSl(4)=(l0/32767.0)*PSlD + PSlMlN 

CHX(4)-SIN(PSI (4)) 

CHY(4)=COS(PSI(4)) 

C 

C COMPUTE, DECOWOSE, AND SUM OUTPUTS FROM NOISE SOURCES 

C ELF nolse via FLM model 

280 CALL GAUSS(IR2,1Rl,0.@,SIGX(I),Z) 

CALL PWFWAY(IR2,1Rl,R(I),A(I),Y) 

H=Z+Y 

C Incident flelds 

HX=M+H*CHX ( I 

HY =HY+H*CHY ( I 

EZ=EZ+H 

380 CONTINUE 

310 CONTINUE 

RETURN 

END 

SUBROUTINE ELFPRM(PZ,GAMMAL,A,SIGMAX,R) 

C CALCULATION OF PARAMETERS FOR SUBROUTINE ELFNSE 

C FOR GAUSS l AN COMPONENT 

C GAMMAL=IMPUSLIVITY; GAMMA=GAMMA FUNCTION 

PX=PZ/( 1.0 + GAMMAL**Z) 

Sl GMAX=SwT(PX) 

C FOR POWER-RAYLEIGH COMPONENT 
R=SQRT((PZ-PX)/GPMMA(I.B + 2.0/A)) 

RETURN 

END 

FUNCTION FILTER(H,HOLD) 

DIMENSION HOLD(100), FW(100) 

COWON /BLOCKZ/ NF I LT, FW 

C SAVE H 

DO 10 J=2,NFILT 

JM=J- 1 

HOLD( J ) +OLD( JM) 

10 CONTINUE 

HOLD( 114 

C COMPUTE RESPONSE 

F I LTER4.0 

DO 20 J=l ,NFILT 

FILTER=FILTER + HOLD(J)*FW(J) 

28 CONTINUE 

RETURN 

END 

FUNCTION PR(C,K) 

C PULSE RESPONSE OF FLAT, HOMOGENEOUS CONDUCTING GROUND 

IF(K.GE.2) GO TO 200 

IF(K.EQ.1) GO TO 100 

C K<=0 

PR4.0 

RETURN 

C K=l 

100 PRX 

RETURN 

C K>=2 

200 PR.C*(1.0/(1.0*K)**0.5 - 1.0/(1.0*K-1.0)**0.5) 

RETURN 

END 

SUBROUTINE ILC(I,Ll,LZ) 

C Converts one Integer to two bytes 

C L1 is high-order byte 

C Memory assignment for Microsoft FORTRAN-8Q 

IMPLICIT LOGICAL (L) 

LOGICAL PEEK 

MI .MLOCS( I ) 

LZ-PEEK(M1) 

Ml=Ml+l 

Ll=PEEK(Ml) 

RETURN 

END 

SUBROUTINE PACK(LR2,LRl) 

C Pack 8 bytes into 6 bytes 

IMPLICIT LOGICAL (L) 

IMPLICIT INTEGER (A-K,M-Z) 

DIMENSION LBIT(8), LRl(l), LR2(1) 

COWON /ADA/ LBIT, LSR, ABASE, ACBU, ACBL 

C 

LRl(5) = LR2( 1) .AND. (X1OFl ) 

LRl(6) = LR2(3) .AND. (X1OF1) 

C 

w 101 =1,4 

J =I +4 

LTEST = LRZ(2) .AND. LBIT(1) 

IF(LTEST.NE.0) LRl(5) = LRl(5) .OR. LBIT(J) 

LTEST = LR2(4) .AND. LBIT(I) 

IF(LTEST.NE.0) LRl(6) = LRl(6) .OR. LBIT(J) 

10 CONTINUE 

C 
RETURN 

END 

SUBROUTINE IDPREP(LGAIN,NSR,NINT,LID) 

C Prepare Identification data 

IWLICIT LOGICAL (L) 

DIMENSION LID(128) 

C 

C LID( 1) - LID(64) are ASCl l message from operator 

WRITE(3,140) 

140 FORMAT(/,' Recording ID = '1 

READ(3,150) (LIDCN). N=1,64) 

150 FORMAT(64Al) 

C 

C Gain (LGAIN), number of seconds (NSR), number of samples/second (NINT) 

LID(65) = LGAIN 

CALL ILC(NSR,LID(66),LID(67)) 

CALL ILC(NINT,LID(68),LID(69)) 

C 

DO 230 N = 70, 128 

LID(N) = 0 

230 CONTINUE 

C 

RETURN 

END 

SUBROUTINE PUTCHR(LCHR) 

C Put single character to disk file 

C Must be sequential 

C Set MEC = NBYTE = 1 to start 

IMPLICIT LOGICAL (L) 

DIMENSION LBUF(128) 

COWON /PCPARM/ NREC, NBYTE 

C 

C Put character In buffer 

C 

LBUF(NBYTE1 = LCHR 

NBYTE = NBYTE + 1 

IF(NBYTE.NE.129) RETURN 

C 

C Put buffer on disk 

C 

WRITE(6, REC-NIEC) LBUF 

NBYTE = 1 

NREC = NREC + 1 

C 

RETURN 

END 

SUBROUTINE PCEND 

C Dunp buffer for PUTCHR 

IMPLICIT LOGICAL (L) 

WON /PCPARM/ NREC, NBYTE 

C 

10 IF(NBYTE.EQ.1) GO TO 20 

CALL PUTCHR(0) 

GO TO 10 

C 

20 RETURN 

END 
Digitized-Noise Processor 

C Program # FHR-353E, 11/11/83 

C Single-channel NLP/ANC algorithm with disk-file input 

C 

IMJLlClT LOGICAL (L) 

DIMENSION 

> SR(3), SI (3), ZW(31, YW(3), ZR(3), Z1(3), YR(3), Y1(3), 

> PHATZ(31, PHATY(3), NNLZ(31, NNLY(31, 

> ZBARR(3), ZBARI (31, YBARR(31, WAR1 (3), 

> SIGZR(3). SIGZI(3), SIGYR(3), SIGY1(3), 

> SIGWR(3,3), SIGW1(3,3), SIGZYR(3,3), SIGZY1(3,3), 

> RWR(3,3), RWI (3,3), RZYR(3,3), RZY1(3,3), 

> CONJ(3), DUM1(3,3), DUM2(3,3), DUM3(3,3), 

> SHATR(3), SHAT1(3), RXYR(3,3), RXY1(3,3), 

> DUM4(6,6), DUM5(6,6), DUM6(6,6), RINVR(3,3), RINV1(3,3), 

> WR(3,3), W1(3,3),ZRL(3),YRL(3),ZIL(3),YIL(3), 

> SIGZRL(3), SIGYRL(31, SIGZIL(3), SIGYIL(3), 

> ZBARRL(31, ZBARI L(3), YBARRL(31, YBARIL(3), 

> LRl(61, LR2(4), LID(128), LBIT(B), CAL(41, 

> ZRF(3), ZIF(31, ZBARRF(31, ZBARI F(3), SIGZRF(3), SIGZIF(3) 

COFMON /FLTR/ FGAIN, FTHR, TOTSMP, COSWAV, SINWAV 

COWON /ADA/ LBIT, LSR, IABASE, IACBU, IACBL 

DATA LGAIN, LCHAN /1,0/ 

10 FORMAT(' Slngle-channel NLP/ANC algorithm wlth disk-', 

> 1 f i le input' ,/,I Program # FHR-353Et ,/) 

C 

C Flxed parameters and constants 

C 

DR4.017453293 

TWOPI 4.2831852 

S I GMAG=l,B 

KMCU -1 00 

FSAMP=100.0 

CF I LT=1.0 

FSIG=l0.0 

I LS=l 

I LY =3 

I LY2=I LY"2 

C 

C Announce name, get number of disk files to process 

C 

WRiTE(3,10) 

CALL ADVERS(3) 

WRITE(3.20) 

20 FORMAT ( ) 

CALL ADSET(LCHAN,LGA.IN,LSD,LBIAS,LCBU,LCBL) 

C - 

C Get number of files to process, open files 

C 

CALL FOPENO(I~~,IERR) 

25 WRITE(3,30) 

30 FORMAT(/,' Number of disk f l les to process = '1 

CALL ATI(NDMAX,IERR) 

IF((NDMAX.LT.l).OR.(NDMAX.GT.4)) GO TO 25 

DO 40 ND=l , WMAX 

LUN = ND + 5 

CALL FGPENI(LUN,IERR) 

40 CONTINUE 

C 

C Start disk file 

C 

DO 9500 ND = 1, WMAX 

LUN = ND + 5 

NREC=l 

NLOC=l 

WRITE(3.50) ND 

WRITE(10,50) NO 

50 FORMAT(/,' ',80('='),/,' Disk # ',Il,/,) - 

I; 

C Get data ID, get callbratlon constants 
C 

DO 80 1=1, 128 

CALL GETCHR(LUN,NREC,NLOC,LID(I),LEND) 

IF(LEND.NE.0) GO TO 85 

80 CONTINUE 

CALL GETIDD(LID,LGAIN,GAIN,NSR,NINT,FREQ) 

CALL GETCAL(LUN,NREC,NLOC,ND,LGAIN,CAL,NSC,NSS,LERR) 

IF(LERR.NE.0) GO TO 85 

GO TO 90 

85 WRITE(3,86) 

86 FORMAT( 1 *** BAD FORMAT - l NPUT F l LE ***' ) 

GO TO 9999 

C 

90 CONTINUE 

C 

C Slgnal parameters 

IF(ND.NE.1) GO TO 300 

MltTE(3,1611) 

160 FORMAT(/,' Signal components (I, Q) = 1) 

CALL ATR(SR(l),IERR) 

CALL ATR(SI(l),IERR) 

WRITE(3,170) 

170 FORMAT(' Axis (X, Y, or Z) for flxed cancellatlon : l) 

READ(3,180) LAXlS 

180 FORMAT(A1) 

LAXlS = LUPPER(LAXIS) 

IAXIS = LAXIS - 87 

C 

C Other parameters 

WRITE(3,200) 

200 FORMAT ( Nsec, Kout = ) 

CALL ATI(NSEC,IERR) 

CALL ATI(KOUT,IERR) 

IJRITE(3,210) 

210 FORMAT(' Nset, F3dB, Thresh = l) 

CALL ATI(NSET,IERR) 

CALL ATR(F3DB,IERR) 

CALL ATR(FTW,IERR) 

KMAX=NSEC+NSET 

KAVG=NSET + NSEC/2 

WRITE(3.230) 

230 FORMAT(l Number of runs = CALL AT1 (NRUN, IERR) 

C 

C Wr I te parameters and head lngs to output f l les 

C 

300 WRITE( 10,101 

C 

WRITE( 10,305) (LID( I), 1 =1,64) 

305 FORMAT(' Nolse-file parameters:l,//,' ID: nl,64A1,1n1) 

WRITE(10,310) GAIN, FREQ 

310 FORMAT(' A/D gain = l,F5,1,' Sampling frequency = l,F5.1) 

WRITE(10,320) NSR 

320 FORMAT(' Total seconds in noise file = ',I51 

WRITE( 10,330) (CAL( I), 1~1.4) 

330 FORMAT(' I-V Cal l brat Ion level s = ',4F9.3) 

WRITE( 10,340) NSC, NSS 

340 FORMAT(' Callbratlon, start tlmes (s) = 1,218) 

C 

WRITE(10,370) %(I), SI (1) 

370 FORMAT(/,l Processing pararneters:l,//, 

> 1 Signal = l, 612.4, ' +j l, 612.4) 

WRITE(l0,380) LAXlS 

380 FORMAT(l Axis for flxed cancellation: l,Al) 

WRITE( 10,390) NSEC, KOUT 

390 FORMAT(' Nsec = 1,16,1 Kout = l.16) 

WRITE(10,391) NSET, DDB, FTHR 

391 FORMAT(' Nset = 1,16,' F3dB = l,F5.2, ' Thresh = l,G12.4) 

C 

C Start run 
C Output headings 

C 

WRITE(3,405) 

WRITE( 10,405) 

405 FCRMAT(/,~X,'T~ ,6X,'MSEL1 ,9X,'Gn1 ,6X,'Ga1 ,6X,'Gta1, 

> 5X,'Gf1 ,6X,'Gtf1 ,/) 

C 

C Initialize parameters 

C 

W 430 J=1,3 

ZRL(J)4.0 

ZIL(J)4.0 

YRL(J)4.0 

YIL(J)4.0 

ZR(J)4.0 

Zl(J)4.0 

YR( J )4.0 

Yl (J)4.0 

SlGZR(J)4.0 

SIGZI (J1d.0 

SlGYR(J)4.0 

SIGYI (J)4.0 

SiGZRL(J)4.0 

SIGZIL(J)4.0 

SlGYRL( J14.0 

SIGYIL(J)4.0 

SlGZRF(J)4.0 

SIGZIF(J)4.0 

DO 420 1 =1,3 

SIGWR(I,J)4.0 

SIGYYI (I,J)d.0 

SIGZYR(I,J)4.0 

SIGZYI(l,J)4.0 

420 CONTINUE 

PHATZ( J ) 4.0 

NNLZ( J Id 

PHATY(J)4.0 

NNLY(J)d 

430 CONTINUE 

FGAIN=2.0 * TWOPI * F3DB / FSAMP 

OD=TWOPI*FSIG/FSAMP 

G0AVG4.0 

GlAVG=B.0 

G2AVG4.0 

G3AVGd.B 

G4AVG4.0 

C Total samples, samples after settl ing, number used to compute gain 

NG-0 

TOTSW4.0 

TOTSS4.0 

C 

C Begin seconds loop 

C 

DO 7000 KN=1, KMAX 

C 

C Begln 10-ms loop 

C 

DO 1000 KM = 1, KMOD 

C 

CALL INTRPT(LINT1 

IF(LINT.EQ.27) GO TO 9999 

OMEGAT=OD*KM 

IF(KM.EQ.KM00) OMEGAT4.0 

COSWAV4OS (OMEGAT) 

SINWAV=SiN(OMEGAT) 

TOTSW=TOTSMP + 1.0 
C 

C Get slx bytes from input f l l e, unpack, and convert to voltages 

C 

DO 530 1=1, 6 

CALL GETCWI(LUN,NREC,NLOC,LRl(I),LEND) 

IF(LEND.NE.0) GO TO 85 

530 CONTINUE 

CALL UWACK(LR2,LRl) 

DO 580 1 = 1, 3 

CALL ADVOLT(LR2(I),LRl(I),LGAIN,LCHANO,YW(I)) 

YWCl) = YW(l) / CAL(I) 

580 CONTINUE 

CALL ADVOLT(LR2(4),LR1(4),LGAIN,LCHANO,ZW(I)) 

ZW( 1) = ZW( 1) / CAL(4) 

C 

C Add signal 

C 

ZW(I)=ZW(l) + SR(l)*COSWAV + SI(l)*SINWAV 

C 

C Detect wlthout nonlinear processlng 

Do 600 1=1,3 

CALL LP(YW(I),YRL(I),YIL(I)) 

600 CONTINUE 

DO 605 14, ILS 

CALL LP(ZW(I),ZRL(I),ZIL(I)) 

605 CONTINUE 

C 

C Apply nonllnear processlng and obtaln narrowband waveforms 

C 

DO 610 1=1,3 

CALL NLP(YW(I),PHATY(I),NNLY(I),YR(I),YI(I)) 

610 CONTINUE 

DO 620 1-1, ILS 

CALL NLP(ZW(I),PHATZ(I),NNLZ(I),ZR(I),ZI(I)) 

620 CONTINUE 

C 

C Fixed nolse cancellation 

L 

ZRF(1) = ZR(1) - YR(IAXIS) 

ZIF(1) = ZI(1) - YI(IAXIS) 

C 

C If settled, accumulate measurement statistics 

C 

IF(KN.LE.NSET) GO TO 900 

TOTSS=TOTSS + 1 

C 

C Blas Zbarl 

CALL MADD( ILS,l,SIGZRL,ZRL,SIGZRL) 

CALL MADD(ILS,l,SIGZIL,ZIL,SIGZIL) 

C Blas Ybarl 

CALL MADD(3,1,SIGYRL,YRL,SIGYRL) 

CALL MADD(3,1,SIGYIL,YIL,SIGYIL) 

C Bias Zbar 

CALL MADO(ILS,l,SIGZR,ZR,SIGZR) 

CALL MADD(ILS,l,SIGZI,ZI,SIGZI) 

C Blas Ybar 

CALL MAW(3,1,SIGYR,YR,SIGYR) 

CALL MADD(3,l,SIGYI,YI,SIGYI) 

C Blas Zbarf 

CALL MADD(ILS,l,SIGZRF,ZRF,SIGZRF) 

CALL MADD(ILS,l,SIGZIF,ZIF,SIGZIF) 

C 

C Covariance (z*) * (yT) 

CALL SMULT(ILS,l,ZI,-1.8,CONJ) 

CALL CMVMLT(ILS,I,ILY,ZR,CONJ,YR,YI,DUMl,DUM2,DUM3) 

CALL MADD(ILS,ILY,SlGZYR,DUMl,SIGZYR) 

CALL MADD(ILS,ILY,SIGZYI,DUMZ,SIGZYI) 

C 

C Covariance (y*) * (yT) 

CALL SMULT(ILY,l,YI,-1.8,CONJ) 

CALL CMVMLT(ILY,l,ILY,YR,CONJ,YR,YI,DUMl,DUMZ,DUM3) 
CALL MADD(iLY,ILY,SIGYYR,DUMl,SIGYYR) 

CALL MADD(ILY,ILYJSIGYYI,DUM2,SIGYYI) 

C 

C End 10-ms loop 

C 

900 CONTINUE 

1000 CONTINUE 

C 

C Zero NLP counters at end of settling period 

C 

IF(KN.NE.NSET) GO TO 1100 

DO 1050 1=1,3 

NNLZ( I ) 4 

NNLY( I )4 

1058 CONTINUE 

1100 CONTINUE 

L 

C Time to make estimate? 

C 

IF(KN.LE.NSET) GO TO 6990 

KSET = KN - NSET 

IF(KSET/KOUT*KOUT.NE.KSET) GO TO 6990 

C 

C Estimate bias and covariances 

C 

CK=1.0/TOTSS 

C 

C Bias vectors 

CALL SMULT(ILS,l,SIGZRL,CK,ZBARRL) 

CALL SMULT(ILS,l,SIGZIL,CK,ZBARIL) 

CALL SMULT(3,1,SIGYRL,CKJYBARRL) 

CALL SMULT(3,1,SIGYIL,CK,YBARIL) 

CALL SMULT(ILS,l,SIGZRF,CK,ZBARRF) 

CALL SMULT(ILS,l,SIGZIF,CK,ZBARIF) 

CALL SMULT(ILS,l,SIGZR,CK,ZBARf?) 

CALL SMULT(ILS,l,SIGZI,CK,ZBARI) 

CALL SMULT(3,1,SIGYRJCK,YBARR) 

CALL SMULT(3,1,SIGYI,CK,YBARI) 

C 

C Covar iance Ryy 

CALL SMULT(3,3,SIGYYRJCK,RYYR) 

CALL SMULT(3,3,SIGYYIJCK,RWI) 

C 

C Covariance Rzy 

CALL SMULT(ILSJ3,SIGZYR,CK,RZYR) 

CALL SMULT(ILS,3,SIGZYIJCK,RZYI) 

C 

C Linear estlmate 

C 

CALL MSE(ILS,SR,SI,ZBARRL,ZBARIL,RMSEL) 

C 

C NLP / Ne-ANC estimate 

C 

CALL SMULT(ILS,I,ZBARR,I.0,SHATR) 

CALL SMULT(ILS,l,ZBARI,1.0,SHATI) 

CALL MSE(ILS,SR,SI,SHATR,SHATI,~SEN) 

C 

C Fixed-cancel lation estimate 

C 

CALL MSE(ILS,SR,SI,ZBARRF,ZBARIF,RMSEF) 

C 

C Invert estlmated covariance Ryy of reference 

C 

CALL CMI(ILY,ILYZ,RWR,RWI,RINVR,RINVI,DUM4,DUM5,DUM6,DET) 

iF(DET.NE.0.0) GO TO 5600 

CALL MZER0(3,3,RINVR) 

CALL MZERO(3,3,RINVI) 

5600 CONTINUE 

L 

C DM1 est imate #0 

C 
CALL DMI(ILS,ILY,ILY2,RZYR,RZYI,ZBARR,ZBARI, > YBARR,YBARI,RXYR,RXYI,SHATR,SHATI,WR,WI,DUM3, 

> DUM4,DUM5,DUM6,RINVRJRINVI) 

CALL MSE(ILS,SR,SI,SHATR,SHATI,RMSE~) 

C 

C DM1 estimate #1 

L 

CALL DMI(ILS,ILY,ILY2,RZYR,RZYI,ZBARR,ZBARI, 

> YBARR,YBARI,RXYR,RXYI,SHATR,SHATI,WR,WI,DUM3, > DUM4,DUM5,DUM6,RINVR,RINVI) 

CALL MSE(ILS,SR,SI,SHATR,SHATI,RMSE1) 

C 

C DM1 estimate #2 

L 

CALL DMI(ILS,ILY,ILY2,RZYR,RZYI,ZBARR,ZBARI, 

> YBARR,YBARI,RXYR,RXYi,SHATR,SHATI,WR,WI,DUM3, 

> DUM4,DUM5,DUM6,RINVRJRINVI) 

CALL MSE(ILS,SR,SI,WATR,SHATI,WSE2) 

C 

C Compute improvement factors and output results 

C 

T=TOTSS/FSAMP 

G0=0.0 

GI 4.0 

G24.0 

G34.0 

G44.0 

IF(RMSEL.EQ.0.0) GO TO 6760 

IF(RMSEN.EQ.0.0) GO TO 6740 

C lmprovement due to clipping 

G0=20.0 * ALOGIB(RMSEL/RMSEN) 

6740 IF(RMSE2.EQ.0.0) GO TO 6760 

C Improvements due to ANC 

GI =20.0 * ALOGl0( RMSEN/RMSE2) 

C Total improvement - NLP + ANC 

6760 G2X0 + G1 

C lmprovement due to flxed cancellation 

IF(RMSEF.EQ.0.0) GO TO 6765 

G3=20.0 * ALOGlB(RMSEN/RMSEF) 

6765 G4X0 + G3 

IF(KN.LE.KAVG) GO TO 6770 

L 

C Average lmprovements over last half of test 

C 

NG*l 

G0AVGX0AVG + G0 

GIAVGXIAVG + GI 

GZAVG42AVG + 62 

G3AVGX3AVG + G3 

G4AVGX4AVG + G4 

6770 CONTINUE 

WRITE(3.6800) T, RMSEL, G0, GI, G2, 63, 64 

WRITE(10,6800) T, WSEL, G0, GI, 62, 63, 64 

6800 FORMAT(' ',F6.lJGl2.4,5F8.1) 

C 

6990 CONTINUE 

7000 CONTINUE 

L 

C Output final averages and covariances 

C 

G0AVG X0AVG/NG 

G 1 AVG X 1 AVG/NG 

GZAVG=GzAVG/NG 

G3AVG=G3AVG/NG 

G4AVG=G4AVG/NG 

WRITE(3,7900) GBAVG, GIAVG, G2AVG, G3AVG, G4AVG 

WRITE(10,7900) GBAVG, GIAVG, GZAVG, G3AVG, G4AVG 

7900FORMAT(l AVG1,12X,5F8.1) 

WRITE(10,8000) 

8000 FORMAT(/,' Zbar = 

WRITE(10,8010) ZBARR(l), ZBARI(1) 
8010 FORMAT(' ',G12.4,' +j ',G12.4) 

WRITE( 10,8020) 

8020 FORMAT(/, ' Y bar = ' 

WRITE(10,8010) (YBARR(I), I=1,3), (YBARI(I), 1=1,3) 

WRITE( 10,8100) 

8100 FORMAT(/,' Ryy = ') 

DO 8130 I=l,ILY 

CALL CMW(10,3,3,I,RYYR,RWl) 

8130 CONTINUE 

WRITE(10,8148) DET 

8140 FORMAT(/,' Det(Ryy) = ',G12.4) 

WRITE(18,8288) 

8200 FORMAT (/ , ' Rxy = ' ) 

CALL CMW(10,1,3,1,RXYR,RXYI) 

WRITE( 10,8300) 

8300 FORMAT(/,' W = ') 

CALL CMW( 18,l ,3,1 *WR,WI ) 

C 

C End run 

C 

9000 CONTINUE 

C 

C Enddisk file 

C 

9500 ENDFl LE LUN 

C 

C 

9999 ENDFILE 18 

STOP 

END 

SUBROUTINE DMI(ILS,ILY,ILYZ,RZYR,RZYI,ZBARR,ZBARI,YBARR,YBARI, 

> RXYR,RXYI,SHATR,SHATI,WR,WI,DUM3,DUM4,DUM5,DUM6,RINVR,RINVI) 

C Estimate of signal using DM1 ANC algorithm 

DIMENSION RZYR(ILS,ILY), RZYI(ILS,ILY), 

> ZBARR(ILS), ZBARI(ILS), YBARR(ILY1, 

> YBARi(ILY), SHATR(ILS), SHATI(ILS), WR(ILS,ILY), 

> WI(ILS,ILY), DUM4(ILYZ,ILYZ), DUM5(ILYZ,ILY2), 

> DUM6(ILYZ,ILY2), RINVR(ILY,ILY), RINVI(ILY,ILY) 

C 

C Initial estimate of covariance Rxy (WR and WI are dunnies) 

C 

CALL SMULT(ILS,l,SHATI,-1.0,SHATI) 

CALL CMVMLT(ILS,l,iLY,SHATR,SHATI,YBARR,YBARI,WR,Wi,DUM3) 

CALL SMULT(ILS,ILY,WR,-1.0,WR) 

CALL SMULT(ILS,ILY,WI,-1.0,WI) 

CALL MADD(ILS,ILY,RZYR,WR,RXYR) 

CALL MADD(ILS,ILY,RZYI,WI,RXYI) 

C 

C Compute weighting matrix W 

C 

CALL CMVMLT(ILS,ILY,ILY,RXYR,RXYI,RINVR,RINVI,WR,WI,DUM3) 

C 

C Calculate improvement vector 

C 

CALL CMVMLT(ILS,iLY,l,WR,WI,YBARR,YBARI,SHATR,SHATI,DUM3) 

C 

C Estimate 

C 

CALL SMULT(ILS,l,SHATR,-1.0,SHATR) 

CALL SMULT(ILS,l,SHATI,-1.0,SHATI) 

CALL MADD(ILS,l,SHAR,ZBARR,SHATR) 

CALL MADD(ILS,l,SHATI,ZBARI,SHATI) 

RETURN 

END 

SUBROUTINE MSE(i,SR,SI,ER,EI,WSE) 

C Mean-square error In estimate of complex vector 

DIMENSION SR(I), SI(I), ER(I), EI(I) 
RMSE4.0 

DO 10 J=l,l 

RMSE=RMSE + (SR(J) - ER(J)Iu*2 

RMSE=RMSE + (SI(J) - EI(J))**2 

10 CONTINUE 

WSE=SQRT(RMSE/I) 

RETURN 

END 

SUBROUTINE NLP(VI,PHAT,NNL,SIHAT,SQHAT) 

COWON /FLTR/ GAIN, T, TOTSMP, COSWAV, SINWAV 

C Signal-nulling loop with clipper 

C 

C Null signal uslng prevlous estimates 

VIN=VI - S1HAT"COSWAV - SQHAT*SINWAV 

C 

C Average power, find rms, set threshold 

PHAT=PHAT + VINU*2 

RMS=SQRT (PHAT/TOTSW 

TAW= RMS*T 

TAUN -RMS*T 

C 

C Cl ip Input and Increment nonl inear event counter 

IF(VIN.LT.TAUP) GO TO 10 

V I N =TAW 

NNL=NNL+ 1 

GO TO 100 

10 IF(VIN.GT.TAUN) GO TO 100 

V I N =TAUN 

NNLENNL+ 1 

GO TO 100 

C 

C Estimate I and Q components of signal 

100 SlHATSlHAT + VIN*COSWAV*GAIN 

SQHAT=SQHAT + ViN*SINWAV*GAIN 

RETURN 

END 

SUBROUTINE LP(VI,SIHAT,SQHAT) 

COWON /FLTR/ GAIN, T, TOTSMP, COSWAV, SINWAV 

C Signal-nulling loop without clipper - L 

C Nu1 l signal uslng previous estimates 

VINIVI - SIHAT*COSWAV - SQHAT*SINWAV 

C 

C Estimate I and Q components of slgnal 

100 SIHAT=SIHAT + VIN*COSWAV*GAIN 

SQHATSQHAT + VIN*SINWAV*GAIN 

RETURN 

END 

SUBROUTINE CMVMLT(I,J,K,A,B,C,D,E,F,DUM) 

C Multlpllcatlon of complex matrices 

C E + jF = (A + JB) * (C + JD) 

DIMENSION A(I,J), B(I,J), C(J,K), D(J,K), E(I,K), F(I,K), 

> DUM(I,K) 

C Real part 

CALL MVMLT(I,J,K,A,C,E) 

CALL MVMLT(I,J,K,B,D,DUM) 

CALL SMULT(I,K,DUM,-1.0,DUM) 

CALL MADD(I,K,E,DUM,E) 

C Imaginary part 

CALL MVMLT(I,J,K,A,D,F) 

CALL MVMLT(I,J,K,B,C,DUM) 

CALL MADD(I,K,F,DUM,F) 

RETURN 

END 
SUBROUTINE CMI(N,N2,A,B,E,F,G,DUMl,DUM2,DET) 

C Inversion of complex matrix 

C R = A + jB; Inverse(R) = E + jF 

C Scratch arrays: G, DUMl, DUM2; dimensions N2*N2, N2=N*2 

DIMENSION A(N,N), B(N,N), E(N,N), F(N,N), 

> G(N2,N2), DUMl(N2,N2), DUM2(N2,N2) 

C 

C Load G 

C 

Do 20 1=1, N 

DO 10 J=l, N 

Il=N + I 

Jl=N + J 

G(I,J)=A(I,J) 

G(Il,Jl)=A(I,J) 

G(I,Jl)=B(I,J) 

G(Il,J)=-B(I,J) 

10 CONTINUE 

20 CONTINUE 

C 

C Invert G 

C 

CALL INVERT(G,N2,DET,DUMl,DUM2) 

C 

C Recover E and F 

C 

DO 120 1=l,N 

DO 110 J=l,N 

J 1 =J+N 

E(I,J)=G(I,J) 

F(I,J)=G(I,Jl) 

110 CONTINUE 

120 CONTINUE 

RETURN 

END 

SUBROUTINE CMW(LUN,N,M,I,CR,CI.) 

C Wrlte out 3-element complex vector or row of matrix 

IMPLICIT LOGICAL (L) 

DIMENSION CR(N,M), CI (N,M) 

C 

WRITE(LUN,10) (CR(I,J), J=1,3), (CI(I,J), J=1,3) 

10 FORMAT(' I' ,3G12,4,' I +j 1' ,3G12.4,'3 ) 

C 

RETURN 

END 

SUBROUTINE UNPACK(LR2,LRl) 

C Unpack 6 bytes Into 8 bytes 

IMPLICIT LOGICAL (L) 

IMPLICIT INTEGER (A-K,M-Z) 

DIMENSION LBIT(81, LRlCl), LRZ(1) 

COMMON /ADA/ LBIT, LSR, ABASE, ACBU, ACBL 

C 

LR2(1) = LRl(5) .AND. (X1OF1 

LR2(3) = LRl(6) .AND, (X1OF1) 

C 

LR2(2) = 0 

LR2(4) = 0 

DO 10 1 = 1, 4 

J=1+4 

LTEST = LRl(5) .AND. LBIT(J) 

IF(LTEST.NE.0) LR2(2) = LR2(2) .OR. LBIT(1) 

LTEST = LRl(6) .AND. LBIT(J) 

IF(LTEST.NE.0) LR2(4) = LR2(4) .OR. LBITCI) 

10 CONTINUE 

C 

RETURN 
END 

SUBROUTINE LIC(Ll,L2,1) 

C Converts two bytes to one integer 

C Ll is high-order byte 

C Memory assignment for Microsoft FORTRAN-80 

IMPLICIT LOGICAL (L) 

Ml MLOCS ( I 

CALL POKE(M1,LZ) 

MI M1+1 

CALL POKE(M1,Ll) 

RETURN 

END 

SUBROUTINE GETCHR(LUN,NIEC,NLOC,LCHR,LEND) 

C Unidirectional single-character read from disk file 

IWLlClT LOGICAL (L) 

DIMENSION LBUF(128) 

C 

I F(NLOC.EQ. 1) READ(LUN, REC=WEC, END=10) LBUF 

C 

LCHR = LBUF(NL0C) 

LEND = 0 

NLOC = NLOC + 1 

IF(NLOC.LE.128) RETURN 

C 

NLOC = 1 

FREC = NREC + 1 

RETURN 

C 

10 LCHR = 0 

LEND = 1 

RETURN 

END 

SUBROUTINE GETIDD(LiD,LGAiN,GAIN,NSR,NINT,FREQ) 

C Decode and display identification data in array LID 

IWLICIT LOGICAL (L) 

DIMENSION LID( 128) 

C 

WRITE(3,120) (LID(N), N=1, 64) 

120 FORMAT(/,' Data identification : ',64A1) 

C 

LGAIN = LID(65) 

GAIN = 0.0 

IF(LGAIN.EQ.0) GAlN = 0.5 

IF(LGAIN.EQ.1) GAIN = 2.0 

IF(LGAIN.EQ.2) GAlN = 8.0 

IF(LGAiN.EQ.3) GAlN = 32.0 

WRITE(3,130) GAIN, LGAIN 

130 FORMAT(' Gain, gain code = ',F5.2,15) 

C 

CALL LIC(LID(66),LID(67),NSR) 

WRITE(3,140) NSR 

140 FORMAT(' Number of 100-sample sets = ', 18) 

C 

CALL LIC(LID(68),LID(69),NINT) 

WRITE(3,150) NlNT 

158 FORMAT(' Number of 3.3-ms intervals / sample = ',I81 

IF(NINT.EQ.0) GO TO 170 

FREQ = 300.0 / NlNT 

WRITE(3,160) FREQ 

160 FORMAT(' Sampling frequency = ',F10.3) 

170 CONTINUE 

C 

RETURN 

END 
SUBROUTINE GETCAL(LUN,WEC,NLOC,ND,LGAiN,CAL,NSC,NSS,LERR) 

C Extract calibration constants from noise-data file 

IWLICIT LOGICAL (L) 

DIMENSION CAL(4), NCAL(4), LR1(6), LR2(4), V(4) 

LERR = 0 

C 

C First disk of set? 

C 

IF(ND.EQ.1) GO TO 5 

NSC = 0 

GO TO 620 

C 

C Do data contain calibration signal? 

C 

5 WRiTE(3,10) 

10 FORMAT(/,' Calibration signal included In data ? I) 

READ(3,20) LCAL 

20 FORMAT( A1 

LCAL=LUPPER(LCAL) 

iF(LCAL.EQ.89) GO TO 200 

C 

C Calibration signal not present, get values from operator 

C 

NSCd 

WRITE(3,110) 

110 FORMAT(' Level of 1-V cal lbration signals for',/, 

> Channels 1, 2, 3, 4 = 

DO 130 N = 1, 4 

CALL ATR(CAL(N),IERR) 

130 CONTINUE 

GO TO 660 

C 

C Calibration signal present, find values 

C 

C Get duration and approximate val ue 

200 CONTINUE 

WRITE(3,210) 

210 FORMAT(' Duration of calibration slgnal in seconds = I) 

CALL ATI(NSC,IERR) 

WRITE(3,220) 

220 FORMAT(' Approximate calibration-signal voltage = I) 

CALL ATR (CALAPX , I ERR 

C 

C Set window and zero averages and counters 

W1 = 0.9 * CALAPX 

W2 = 1.1 * CALAF'X 

DO 230N=l,4 

NCAL(N) = 0 

CAL(N) = 0.0 

230 CONTINUE 

C 

C Process calibration signal 

C 

DO 400 NS = 1, NSC 

Do 390 N1 = 1, 100 

C 

C Get and unpack data set 

C 

C Get six bytes from disk 

Do 300 N=1, 6 

CALL GETCWI(LUN,NREC,NLOC,LRl(N),LEND) 

iF(LEND.NE.0) GO TO 9Q0 

300 CONTINUE 

C 

C Unpack and print 

CALL UNPACK(LR2,LRl) 

DO 310 N=l, 4 

CALL ADVOLT(LR2(N),LRl(N),LGAlN,LCHD,V(N)) 

310 CONTINUE 

C 
C Average if within window 

DO 370 N = 1, 4 

V(N) = ABS(V(N)) 

IF((V(N).LT.Wl) .OR. (V(N).GT.WZ)) GO TO 360 

CAL(N) = CAL(N) + V(N) 

NCALCN) = NCAL(N) + 1 

360 CONTINUE 

370 CONTINUE 

C 

390 CONTINUE 

400 CONTINUE 

C 

C Scale averages 

C 

DO 420 N = 1, 4 

IF(NCAL(N).LE.a) GO TO 800 

CAL(N) = CAL(N) / NCAL(N) 

420 CONTINUE 

WRITE(3,430) (CAL(N), N=1,4) 

430 FORMAT(' 1-V Cai ibration levels = l.4F9.3) 

C 

C Get start time 

C 

600 CONTINUE 

WRITE(3,610) 

610 FORMAT(' Start time = ') 

CALL ATI(NSS,IERR) 

C 

C Skip to start 

C 

620 NSKlP = NSS - NSC 

IF(NSKlP.EQ.0) RETURN 

DO 650 NS = 1, NSKIP 

DO 640 N1 = 1, 100 

DO 630 N=l, 6 

CALL GETCHR(LUN,N?EC,NLOC,LRl(N),LEND) 

IF(LEND.NE.0) GO TO 900 

630 CONTINUE 

640 CONTINUE 

650 CONTINUE 

C 

C Done - normal exit 

C 

RETURN 

C 

C Error 

L 

800 WRITE(3,810) N 

810 FORMAT(/,' *** NO CALIBRATION SIGNAL ON CHANNEL ',15) 

900 LERR=l 

RETURN 

END 
9 1 

APPENDIX D. GAIN STATISTICS 

NUMERl CAL DATA - AVERAGES 

Single-channel NLP/ANC algorithm with disk-file Input 

Program # FHR-353E 

Noise-file parameters: 

A/D gain = 2.0 Sampling frequency = 100.0 

Total seconds In noise file = 575 

Processing parameters: 

Signal = .I414 +J 0.000 

Nsec = 270 Kout = 2 Nrun = 2 

Nset = 2 F3dB- = 1.00 Thresh = inf i n lte 

Set Axis Date Time Comnent NSC NSS Gn Ga Gta Gtf 

1-1 Y 01/14 16:45 GWX, SD 10 17 0.0 9.5 

1-1 0.0 10.4 

1-2 0 17 0.0 7.8 

1-2 0.0 10.7 

1-3 0 17 0.0 16.3 

1-3 8.8 1.3 

2-1 Z 01/14 18:00 GWX 10 12 0.0 2.9 

2- 1 0.0 .9 

2-2 0 12 0.0 -a8 

2-2 0.0 -.8 

2-3 0 12 8.8 -3.9 

2-3 0.0 -3.6 

01/11 05:30 GWX 11 

0 

0 

01/11 06:20 GWX 16 

0 

0 

01/11 07:30 GWX 9 

0 

0 
8-1 X 01/07 09:40 GWX 11 18 0.0 6.4 

8- 1 8.8 12.6 

8-2 0 18 0.0 2.5 

8-2 8.8 -.9 

8-3 0 18 0.0 9.1 

8-3 0.0 10.8 

01/05 ll:10 GWX 528s 40 45 0.0 

0.0 

0 10 0.0 

0.0 

0 10 0.0 

0.0 

01/07 11:20 GWX 9 11 0.0 

0.0 

0 10 0.0 

0.0 

0 10 0.0 

0.0 

X 01/06 19:15 8 10 0.0 

0.0 

0 10 0.0 

0.0 

0 10 0.0 

0.0 

Y 01/12 19:30 GWX 3 4 0.0 

0.0 

0 4 0.0 

0.0 

0 4 0.0 

0.0 

X 01/14 14:00 Disturbed 3 6 0.0 

0.0 

0 6 0.0 

0.0 

0 6 0.0 

0.0 

Y 01/07 10:30 GWX 10 14 8.8 

0.0 

0 14 0.0 

0.0 

0 14 0.0 

0.0 

NUMERl CAL DATA - AVERAGES 

Single-channel NLP/ANC algorithm with disk-f 1 le input 

Program 1 FHR-353E 

Noise-f I le parameters: 

A/D gain = 2.0 Sampling frequency = 100.0 

Total seconds in noi se f i le = 575 

Processing parameters: 

Signal = .I414 +j 0.000 

Nsec = 270 Kout = 2 Nrun = 2 

Nset = 2 F3dB = 1 .OO Thresh = 2.828 

Set Axis Date Time Comnent NSC NSS Gn Ga Gta Gtf 

1-1 Y 01/14 16:45 GWX, SD 10 17 .4 9.5 

1 - 1 -.4 3.2 

1-2 0 17 -2.5 3.8 
GWX 18 

0 

0 

GWX 11 

0 

0 

GWX 16 

0 

0 

GWX 9 

0 

0 

GWX 11 

0 

0 

GWX 528s 40 

0 

0 

GWX 9 

0 

0 
01/12 19:30 GWX 3 4 .4 15.1 

-.3 11.7 

0 4 -.4 22.8 

.I 12.8 

0 4 -.I 19.5 

-. 1 9.9 

13-1 X 01/14 14:00 Disturbed 3 6 0.0 -2.6 

13-1 .0 .2 

13-2 0 6 0.0 4.1 

13-2 8.0 5.3 

13-3 0 6 0.0 -2.8 

13-3 -.l -1.7 

14-1 Y 01/07 10:30 GWX 10 14 .4 7.4 

14-1 -.6 4.3 

14-2 - 0 14 .2 12.9 

14-2 -.B 12.3 

14-3 0 14 .a 6.6 

14-3 6.3 -5.2 

NUMERICAL DATA - AVERAGES 

Single-channel NLP/ANC algorithm with disk-file input 

Program # FHR-353E 

Noise-f i le parameters: 

A/D gain = 2.0 Sampling frequency = 100.0 

Total seconds in nolse file = 575 

Processing parameters: 

Signal = .I414 +j 0.000 

Nsec = 270 Kout = 2 Nrun =.2 

Nset = 2 F3dB = 1 .a0 Thresh = 1.414 

Set Axis Date Time Comnent NSC NSS Gn Ga Gta Gtf 

1-1 Y 01/14 16:45 GWX, SD 10 17 4.8 4.1 

1-1 .8 1.0 

1-2 0 17 -2.4 3.1 

1-2 1.9 10.2 

1-3 0 17 4.8 4.1 

1-3 .8 1.0 

2-1 Z 01/14 18:00 GWX 18 12 -9.2 3.6 

2- 1 -3.9 1.8 

2-2 0 12 -.9 -.6 

2-2 3.2 3.8 

2-3 0 12 -1.3 -3.8 

2-3 2.6 -3.4 

4-1 X 01/11 05:30 GWX 11 13 -3.6 5.6 

4- 1 .2 3.4 

4-2 0 13 -.2 7.8 

4-2 -.I 4.0 

4-3 0 13 -.3 3.8 

4-3 -.6 5.1 
GWX 16 19 

0 19 

0 19 

12 22 

0 22 

0 22 

GWX 9 12 

0 12 

0 12 

GWX 11 18 

0 18 

0 18 

GWX 528s 40 45 

0 10 

0 10 

G WX 9 11 

0 10 

GWX 3 4 

0 4 

0 4 

Disturbed 3 6 

0 6 

0 6 

GWX 10 14 

0 14 

0 14 
NUMER l CAL DATA - AVERAGES 

Single-channel NLP/ANC algorithm with disk-file 

Program # FHR-353E 

Input 

Noise-file parameters: 

A/D gain = 2.0 Sampling frequency = 100.0 

Total seconds in noise file = 575 

Processing parameters: 

Signal= .I414 +j 0.000 

Nsec = 270 Kout = 2 Nrun = 2 

Nset = 2 F3dB = 1 .a0 Thresh = 0.707 

Set Axis Date Time Comnent NSC NSS Gn Gta Gtf 

16:45 GWX, SD 10 17 -3.1 

-.4 

0 17 -1.9 

5.0 

0 17 11.5 

2.8 

18:00 GWX 10 12 -8.7 

-5.1 

0 12 -1.7 

2.9 

0 12 -1.2 

9.3 

4-1 X 01/11 05:30 GWX 11 13 4.9 

4- 1 3.2 

4-2 0 13 .3 

4-2 1-11 

4-3 0 13 -.0 

4-3 -2.2 

5-1 Y 01/11 06:20 GWX 16 19 7.5 

5- 1 1.7 

5-2 0 19 7.0 

5-2 1.5 

5-3 0 19 3.7 

5-3 .9 

7-1 Z 01/11 07:30 GWX 9 12 -1-8 

7- 1 -2.1 

7-2 0 12 -7.3 

7-2 4.1 

7-3 0 12 -3.3 

7-3 -.4 

8-1 X 01/07 09:40 GWX 11 18 -.3 

8- 1 2.2 

8-2 0 18 -4.6 

8-2 -4 -0 
01/05 11:10 GWX 528s 40 45 

0 10 

0 10 

01/07 11:20 GWX 9 11 

0 10 

0 10 

01/12 19:M GWX 3 4 

0 4 

0 4 

01/14 14:00 Disturbed 3 6 

0 6 

0 6 

01/07 10:M GWX 10. 14 

0 14 

0 14