EPA/600/R-10/110
October 2010
AUTOMATED LONG-TERM REMOTE MONITORING OF
SEDIMENT-WATER INTERRACIAL FLUX
-------
AUTOMATED LONG-TERM REMOTE MONITORING
OF
SEDIMENT-WATER INTERRACIAL FLUX
By
Bob K. Lien1
Carl G. Enfield2
1U.S. EPA National Risk Management Research Laboratory
Office of Research and Development
U.S. Environmental Protection Agency
Cincinnati, OH 45268
2Pegasus Technical Services,
Cincinnati, OH 45219
-------
Notice
The U. S. Environmental Protection Agency through its Office of Research and Development
funded the research described here. This project was conducted under an approved Quality
Assurance Project Plan (743-Q2-1). This report has been subjected to the Agency's peer review
and has been approved for publication as an EPA document. Mention of trade names or
commercial products does not constitute endorsement or recommendation for use.
The authors wish to acknowledge George Enfield (Pegasus Technical Service Inc.), Patrick Clark
and Robert Ford (USEPA/ORD) for their assistance in the field activities. Kim McClellan
(USEPA/ORD) is acknowledged for providing editorial review of this report. Steve Acree
(USEPA/ORD-Ada) and Robert Ford (USEPA/ORD) are acknowledged for providing Red Cove
staff gauge data and technical review of this report.
-------
Foreword
The U.S. Environmental Protection Agency (EPA) is charged by Congress with protecting the
Nation's land, air, and water resources. Under a mandate of national environmental laws, the
Agency strives to formulate and implement actions leading to a compatible balance between
human activities and the ability of natural systems to support and nurture life. To meet this
mandate, EPA's research program is providing data and technical support for solving
environmental problems today and building a science knowledge base necessary to manage our
ecological resources wisely, understand how pollutants affect our health, and prevent or
reduce environmental risks in the future.
The National Risk Management Research Laboratory (NRMRL) is the Agency's center for
investigation of technological and management approaches for preventing and reducing risks
from pollution that threaten human health and the environment. The focus of the Laboratory's
research program is on methods and their cost-effectiveness for prevention and control of
pollution to air, land, water, and subsurface resources; protection of water quality in public
water systems; remediation of contaminated sites, sediments and groundwater; prevention
and control of indoor air pollution; and restoration of ecosystems. NRMRL collaborates with
both public and private sector partners to foster technologies that reduce the cost of
compliance and to anticipate emerging problems. NRMRL's research provides solutions to
environmental problems by: developing and promoting technologies that protect and improve
the environment; advancing scientific and engineering information to support regulatory and
policy decisions; and providing the technical support and information transfer to ensure
implementation of environmental regulations and strategies at the national, state, and
community levels.
This publication has been produced as part of the Laboratory's strategic long-term research
plan. It is published and made available by EPA's Office of Research and Development to assist
the user community and to link researchers with their clients.
Sally Gutierrez, Director
National Risk Management Research Laboratory
IV
-------
Abstract
Advective flux across the sediment-water interface is temporally and spatially heterogeneous in
nature. For contaminated sediment sites, monitoring spatial as well as temporal variation of
advective flux is of importance to proper risk management. This project was conducted to
develop a ruggedized advective flux meter capable of unattended long-term remote monitoring
of the sediment-water interfacial flux. The flux meter, based on the heat-pulse technique, is
capable of quantitatively measuring bi-directional seepage flux from 0 to 8 cm/day with the
current configuration. The system has automatic data acquisition, real time data display, and
signal display/analysis capability. A remotely controlled monitoring module was included to
permit system monitoring and modification over the internet. The instrument has undergone
several calibrations to establish relations between the flow rate and the heat-pulse travel time.
For very low flow rates, an alternate approach was adopted based on the ratio of temperature
rises. Using this approach, it is possible to accurately interpret flow rates down to zero flow
conditions. Depending on the magnitude, flow rate can be derived as a function of peak
temperature arrival time or as a ratio of rise in temperature on both sides of the heater. Using
these two methods, reliable measurements can be made over a wide bi-directional range of
flows. In the field operation, the flow across the sediment-water interface is isolated by a
dome chamber and funneled through a flow sensor where the flow rate is measured. The
advective flux through the sediment-water interface, in terms of the vertical Darcy velocity, is
calculated by dividing the flow rate through the flow sensor by the cross-section area of the
dome.
The flux meter underwent a continuous long-term field test from June 2008 to April 2009 in the
Red Cove area of Plow Shop Pond, Fort Devens, Massachusetts. Two flux meters were placed
side by side; the unattended operation ran automatically with a self-sufficient power supply
system. The measurement frequency was changeable through the wireless network. Seepage
flux data was available for download after each measurement cycle. All of the data suggest
very low flow and the majority of the data suggest an upward flow direction from the
groundwaterto the surface water. The average advective flux, in terms of Darcy velocity, for
flux meter unit SN300 and SN303 was 0.44 ± 0.34 mm/day and 0.76 ± 0.33 mm/day,
respectively. Calculated seepage flux based on hydraulic conductivity and gradient data
suggests the average flux was 1.12 ± 0.71 mm/day. During the winter months when the surface
water was frozen, the data superficially suggests there was significant rapid oscillatory pumping
of the ice-covered water during this period of time.
IV
-------
Table of Contents
Notice iii
Foreword iv
Abstract iv
Table of Contents v
List of Figures vii
1 Introduction 1
1.1 Background 2
1.2 Scope and Objectives 3
2 Flux Meter Design 4
2.1 Basic Concept 4
2.2 Design of Measurement System 6
2.3 Flux Meter Mechanicals 10
3 Theory of Operation 13
4 Laboratory Calibrations 15
4.1 Thermistor Calibration 15
4.2 Circuitry Description and Calibration 18
4.3 Flow Sensor Calibration 20
5 Extended Time laboratory Study 24
6 Field Deployment 25
6.1 Solar Power System 29
6.2 Heater Response 31
6.3 Pressure Transducer Response 32
6.4 Flux Meter Response 35
-------
6.4.1 Alternative Approach 38
6.4.2 Flux Meter Results 41
6.5 Flux Estimation using Darcy's Law 46
7 Summary and Discussion 49
8 References 52
Appendix A Visual Basic Code 54
frm Data Acquisition 54
frm Data Analysis 85
frm Differentiates 91
frmMain 98
frmSignal Display 101
Appendix B Circuit Diagram of Measurement System 118
Appendix C Signal Conditioning Circuit Board Layout 120
Appendix D Interface Circuit Board Layout 121
Appendix E Calibration constants of flow sensor SN300, SN303 122
VI
-------
List of Figures
Figure 1 Schematic drawing of the measurement system 5
Figure 2 Field implementation of the advective flux meter 6
Figure 3 Conceptual drawing of the flow sensor. The flow tube is constructed with
Polyetheretherketone (PEEK) material. Note that T4 is on the outside surface of the heater not in
contact with water 7
Figure 4 Block diagram of the measurement system 8
Figure 5 Signal conditioning and data interface boards as used in the field 9
Figure 6 Flow tube and differential pressure transducer with drawings of flow tube case end plates.... 11
Figure 7 Flow meter Dome Chamber. (A) Photo of dome halves separated; (B) Photo of dome as it
would be installed in the field; (C) Assembled unit with the flow tube and the differential pressure
transducer installed inside the chamber (shown from the bottom); (D) Drawing of the dome chamber. A
flange is welded around the bottom cylinder that serves as a gauge for installing the unit to a consistent
depth of 15 cm 12
Figure 8 Thermistors in an aluminum block used for calibration 15
Figure 9 Example results of a typical thermistor calibration. The blue-diamond marks the measured
versus the estimated temperature. The red-cross marks the error (on secondary Y-axis) associated with
each estimated temperature. In this particular example, the root mean squared error of the estimator
was4.5E-3°C 17
Figure 10 Comparative plot of calibration curves for flow sensor unit SN300 and SN303. For the
purpose of distinguishing flow direction, flow toward T5 and T6 is considered positive flow while flow
toward T3 and T2 is considered negative flow 20
Figure 11 Comparison between measured and estimated flows for T2 and T6 using Equation 3. The error
associated with each estimated flow rate is shown on the secondary Y-axis. The RMSE of T2 and T6 flow
estimation was 0.095 and 0.090 ml/min, respectively 21
Figure 12 Comparison between measured and estimated flow for T3 and T5 using Equation 3. The error
associated with each estimated flow rate is shown on the secondary Y-axis. The RMSE of T3 and T5 flow
estimation was 0.099 and 0.112 ml/min, respectively 22
Figure 13 Mean absolute error in defining flow rate as a function of peak temperature arrival time 23
Figure 14 Demonstration of long-term measurement repeatability in the laboratory 24
vii
-------
Figure 15 Aerial photo of the Red Cove area of Plow Shop Pond, Fort Devens, Massachusetts. The
yellow dots mark the location of the instruments; the blue square marks the location of STAFF1 Red
Cove staff gauge 25
Figure 16 Placement of a perforated PVC mat prior to the flux meter installation. The perforated mat
(4' x 4' in size, 1/4" thick with 1/4" hole-diameter) was placed on the sediment with a cutout in the
center for the flux meter 26
Figure 17 Deployed locations [latitude, longitude] of the flux meter and piezometer. Water lily covered
most of the cove at the time of installation 27
Figure 18 Time trend of water temperature and average air temperature at the Red Cove of Plow Shop
Pond, Fort Devens, Massachusetts. Air temperature data were obtained from the Fitchburg Municipal
Airport meteorological station, National Climatic Data Center of the National Oceanic and Atmospheric
Administration 28
Figure 19 Solar power system was installed uphill south of Red Cove 29
Figure 20 Battery Voltage Compared to Sky Cover. Sky cover data was obtained from the Fitchburg
Municipal Airport meteorological station, National Climatic Data Center of the National Oceanic and
Atmospheric Administration 30
Figure 21 Heater Temperature Rise Versus Battery Voltage 31
Figure 22 Comparative plots of water elevations (primary Y-axis) at piezometers PZ5S, PZ5C, PZ5N and
Red Cove STAFF1 gauge station and time-series precipitation (secondary Y-axis). The daily precipitation
data of Red Cove area was obtained from the Fitchburg Municipal Airport meteorological station,
National Climatic Data Center of the National Oceanic and Atmospheric Administration 32
Figure 23 Temporal vertical hydraulic gradients at each piezometer location. Temporal variation of the
three piezometers appeared to follow similar trend. The gradient peaked in August 2008 and then
gradually descended overtime 33
Figure 24 Effect of ambient temperature on the movement of heat-pulse in zero flow condition showing
increasing peak arrival time with decreasing temperature 35
Figure 25 Martin and Lang's Thermal conductivity of water as a function of temperature 36
Figure 26 Ambient temperature effect on estimated flow rate 39
Figure 27 Comparative plot of measured flow rate versus temperature corrected estimated flow rate.
The RMSE of the estimated flow rate was 0.013 ml/min 40
Figure 28 Daily average advective flux data at the Red Cove of Plow Shop Pond, Fort Devens,
Massachusetts. The data were corrected according to the ambient water temperature 41
Figure 29 Histogram of advective flux measured at Red Cove flux meter location SN300 42
viii
-------
Figure 30 Histogram of advective flux measured at Red Cove flux meter location SN303 42
Figure 31 Monthly average advective flux and precipitation 43
Figure 32 Apparent flow reversals during measurement cycle. Winter processed data (A), spring
processed data (B) 45
Figure 33 Falling-head permeameter test to measure sediment hydraulic conductivity 46
Figure 34 Falling head permeameter test at the Red Cove of Plow Shop Pond, Fort Devens,
Massachusetts. Hydraulic conductivity is proportional to the slope of the best fit line 48
Figure 35 Comparison between measured and calculated advective flux. The calculated flux was based
on hydraulic conductivity and vertical hydraulic gradient measurements 48
IX
-------
1 Introduction
Advective flux (also referred to as seepage flux) at the groundwater-surface water interaction
zone is the magnitude and direction of water movement across a sediment-water interface.
The convention is that positive advective flux values indicate flow of groundwater to the
surface water and negative flux values indicate flow from the surface water to the
groundwater. The seepage flux is measured in terms of volume per unit area per unit time
(L3/L2/T), or simply translates to unit length per unit time (L/T).
There are many indirect and direct methods of estimating advective flux at the groundwater-
surface water interaction zone (Kalbus, Reinstorf et al. 2006). The indirect methods include
hydraulic conductivity and hydraulic gradient, heat tracer, hydrographic analysis, artificial
and/or environmental tracers, or hydrogeological mapping. The most common direct method
is the use of seepage flux meters. The basic concept of the flux meter is to cover and isolate an
area of the sediment-water interface with a submerged chamber opened at the bottom to
measure the volume of water exchange between the groundwater and the surface water over
an interval of time. The classic flux meter design of Lee (1977) consists of a cutoff 55 gallon
(200L) drum inserted into the sediment. A stopper with a tube is inserted into a port on the top
of the drum and a plastic bag is attached to the tube with rubber bands. The time interval when
the bag is connected and subsequently disconnected is recorded, as well as the change in the
volume of water in the bag. The basic Lee-type flux meter has its shortcomings (Schincariol and
McNeil 2002), and many factors could affect the accuracy of the flux measurements (Murdoch
and Kelly 2003). One of the major disadvangates of the Lee-type flux meter is that it only
measures the accumulated seepage over the deployed time period, which does not provide any
data on temporal changes of seepage during that time period.
Over the years, various modifications have been made to the basic design of Lee's flux meter.
The same principle of an inverted chamber to isolate and direct seepage flux remains.
However, instead of using a plastic bag or similar collection device, different instruments are
used to continually measure the time-series change in flow automatically. From the literature,
examples of these modified devices include the ultrasonic meter (Paulsen, Smith et al. 2001),
the heat-pulse meter (Taniguchi and Fukuo 1993; Lien 2006; Smith, Herne et al. 2009), the
continuous heat-type meter (Taniguchi and Iwakawa 2001), the dye-dilution meter (Sholkovitz,
Herbold et al. 2003) and the electromagnetic meter (Rosenberry and Morin 2004). The range of
measurable flow rates varies with each device. Most of the devices were deployed for periods
of days or weeks in various subsequent field studies (Taniguchi, Burnett et al. 2003; Taniguchi,
Burnett et al. 2006; Taniguchi, Ishitobi et al. 2007; Swarzenski and Izbicki 2009; Mwashote,
Burnett et al. 2010), with the exception of (Fritz, Mendoza et al. 2009) who deployed an
ultrasonic meter in a river for a period of two months. None of these devices were specifically
-------
developed to include remote data transmittal capability. Although the above mentioned heat-
pulse meters are similar in that they all employ heat-pulse technique, Lien's and Smith's heat-
pulse meters are also capable of measuring bi-directional water flow at the groundwater-
surface water interaction zone.
Temporal hydrological variation exists in nature. For reliable environmental risk management,
it is essential to collect and analyze repeated measurements to evaluate changes in conditions.
Knowledge of the variation in seepage flux is invaluable to proper remedial design at a
contaminated sediment site. The need to develop a flux meter that is sustainable to extensive
environmental exposure is therefore warranted. This document presents the development and
field testing of a ruggedized advective flux meter for measurement of advective water flux
across the sediment-water interface. Evolved from a previously developed prototype by Lien
(2006), the flux meter was enhanced and configured for long-term deployment and remote
data transmittal capabilities.
1.1 Background
Previously, a prototype bi-directional advective flux meter for measuring water flux across the
sediment-water interface was developed by U.S. EPA and field tested at multiple sites (Lien
2006). The prototype consists of a dome chamber opened at the base to isolate and direct flow
through a flow sensor for measurement. The flow sensor employs a principle of heat-pulse
technique (Campbell, Calissendorff et al. 1991; Choi, Lee et al. 1991; Taniguchi and Fukuo 1993;
Paillet, Crowder et al. 1996; Ren, Kluitenberg et al. 2000; Mortensen, Hopmans et al. 2006)
derived from the effect of flowing water velocity on heat flow propagation along a flow tube.
When a temperature gradient exists along the flow path, the heat will be transported in the
flowing water by advection and conduction. The advective component of the heat flow is
related to the flow velocity of the fluid. Empirical relations between flow rate and heat-pulse
peak arrival time were established; flow rate could be derived as a function of peak
temperature arrival time, or as a function of first temporal moment of the heat-pulse. The
complete prototype system was also equipped with automatic data acquisition and signal
display/analysis capability. The prototype flux meter has undergone field tests at three
different settings. The first test site was a shallow turbulent stream at Santo Domingo,
Nicaragua. The second test site was at Grand Calumet River, Hammond, Indiana, which is a slow
moving river with fine texture organic rich sediment. And the third test site was at Lake
Hartwell, South Carolina which is a deeper water reservoir with sandy sediment. The field
deployments were successful in that the flux meter was operational and bi-directional flow
measurement capability was demonstrated. Nevertheless, ideas for further improvement
emerged from lessons learned from earlier work.
-------
The early prototype flux meter had notable limitations:
1. The prototype required constant physical presence of an operator and is not suitable for
unattended long-term monitoring from a remote location;
2. In fast moving water bodies, differential pressures were observed between inside and
outside the dome chamber, but could not be recorded in real-time to correct for the
Bernoulli Effect;
3. It is necessary to separate gas flow from liquid flow due to differences in thermal
conductivity. The prototype had potential problems with gas venting configuration. The
prototype used a vent tube to the surface, which was successful when the surface water
was calm and with little change in the surface water elevation. However, when the
surface water is dynamic unintended flow could be induced;
4. The prototype flux meter was not built to sustain lengthy field deployment.
1.2 Scope and Objectives
The scope of this project is to aim at improving the limitations and shortcomings of the
prototype flux meter and further enhance its capability for long-term operation. The objectives
of this project include:
• Redesign and ruggedize the flux meter mechanicals for long-term deployment;
• Redesign an imbedded low power data collection system for unattended operation, and
capable of being queried and modified through the internet;
• Redesign circuitry to incorporate differential pressure measurements for use in
determining differences in pressure between the inside and the outside of the dome to
correct for the Bernoulli Effect;
• Redesign signal conditioning at the data collection system for optional pressure
transducers (up to four with 0-10 VDC output and 15 VDC excitation) and temperature
transducers (lOkQ thermistors up to two);
• Redesign gas venting mechanism to minimize the interference of gas ebullition on flow
measurements;
• Demonstrate the automated long-term seepage flux measurement at a field location.
-------
2 Flux Meter Design
The design of the meter incorporates the same principle of the previous prototype with
improvements to ensure good operational integrity over an extended period of time.
Significant modifications were made to address power consumption, remote data transmittal,
and long-term environmental exposure issues. The meter was designed to withstand up to 30
feet (9.1 m) of hydraulic pressure. The in situ components of the meter consist of a cylinder
that when pushed into the sediment isolated 2,922.5 cm2 of interfacial sediment area, a dome
attached to the cylinder with a water-tight seal, and a flow tube that allows water to freely flow
between inside the dome and the surface water at the rate equal to the specific discharge rate
across the sediment-water interface.
2.1 Basic Concept
The principles of operation of a groundwater/surface water flow meter were originally
presented by (Zuber 1970). If a known area of the groundwater/surface water interface is
covered by a chamber which is open at the bottom and has a pipe outlet near the top (Figure
1), the following expression can be written to express the flux (q) of water normal to a plane
passing through the open bottom of the container:
where
• v is the velocity of the water in the measuring pipe
• A is the area covered by the chamber
• a is the cross sectional area of the measuring pipe.
Two assumptions must be satisfied for the validity of Equation [1]. The first assumption
requires that area A be tightly covered by the instrument, so the flow through area A is equal to
the flow through the measuring pipe. The second assumption requires a negligible loss of
potential in the pipe compared with the loss along the edge of the cylinder which is tightly
pressed into the bed material. In other words, the presence of the measuring instrument
should not disturb natural flow lines. Taniguchi and Fukuo (1993) concluded that inaccuracy of
the flow measurement caused by pressure drop in the pipe could be minimized by using a large
diameter pipe (i.e. > 0.9 cm ID).
-------
Degassing valve
Measuring Pipe
Cross sectional area (a)
Length (I)
Collection funnel
Cross sectional area (A)
Sediment water interface
Figure 1 Schematic drawing of the measurement system.
A heat-pulse flow sensor is incorporated along the flow tube to measure advective flow rate.
Flow through an isolated section of the sediment-water interface is captured and directed
through the flow sensor. As shown in a conceptual drawing in Figure 2, a typical field
installation would include a data collection/controller system, internet interface, underwater
instrumentation and power source. The device is capable of minimizing interference of gas
ebullition on flow measurements by installing an automatic degassing valve (Plast-0-matic
Automatic Vent Valve DGV050V-PV1) on the dome. The top of the dome is slightly convex
allowing any gas generated from the sediments to escape through the degassing valve rather
than into the flow tube.
A differential pressure transducer (Setra Systems 23010R5PB2FEB 10.5PSID1) is incorporated to
measure the difference in pressure between the inside and the outside of the dome. In high
velocity streams, the orientation of the discharge can create pressure gradients (the Bernoulli
Effect) which would potentially influence flow through the device (Shinn, Reich et al. 2002;
Rosenberry 2008). During installation, when significant stream flows are present, the dome
should be rotated to an optimal orientation where this pressure difference is minimal.
Alternatively, the differential pressure transducer can also be used to determine the difference
in hydraulic head from two piezometers that are installed at different elevations in the
sediment formation.
1 Mention of trade names or commercial products does not constitute endorsement or recommendation for use.
-------
Solar collectors
Data collection system with internet interface
Air water interface
Automatic gas vent valve
Battery backup
Flow sensor
Differential pressure
transducer
Figure 2 Field implementation of the advective flux meter.
2.2 Design of Measurement System
Power for the system operation was provided by a 300 W solar panel with regulator and
battery backup provided by Sun Wise Systems. A pulse of heat is injected into in the center of
the flow sensor using a resistive heater (Figure 3). The heater is a section of 3/8" stainless steel
pipe with a 5 mm stainless steel tubing helix fused to the pipe and filled with the resistive
heating element. The migration of this heat-pulse is monitored using four thermistors (General
Electric BR14KA1003K 10K ohm ± 10% @ 25°C housed in 20 GA 304 Stainless Steel thin wall
hypo tubing, closed end). Two of the thermistors are placed on each side of the heater to
measure the bi-directional movement of the heat-pulse. A fifth thermistor is attached to the
surface of the heater to monitor the temperature rise and fall of a heating cycle. A reference
thermistor measures the temperature of the flow case. The five thermistors shown in Figure 3
are identified as T2, T3, T4, T5 and T6, with T2 being closest to the cable end of the device. This
nomenclature is used throughout this report to identify the individual thermistors. Not shown
in Figure 3 are Tl (a reference thermistor imbedded in the case of the flux meter) and T7 (a
thermistor outside the dome used to measure the ambient water temperature). For the three
remaining channels on the 10-channel counter, Channel 8 is the differential pressure
transducer, Channel 9 and 10 can be connected to either pressure transducers or thermistors
that are connected above the water surface.
-------
1.5 0.5 2.0
Heater
Flow tube
.
O-Ring
Note: Dimensions In inches not to scale
Figure 3 Conceptual drawing of the flow sensor. The flow tube is constructed with
Polyetheretherketone (PEEK) material. Note that T4 is on the outside surface of the heater
not in contact with water.
The system is controlled by an embedded single board PC with Windows XP as the operating
system (WinSystems PPM-GX500-G). Software installed on the single board computer included
Symantec PC Anywhere to provide a user interface between the remote and local computers,
and Tethys Automation Anywhere to automate the remote boot and program initialization
sequence, as well as a Visual Basic code used to operate the device. The application that
operated the device was written in Visual Basic. The code is exhibited in Appendix A. A cellular
modem (Janus Terminus GSM/GPRS Terminal) was used to provide the field modem for the T-
Mobile virtual private network (VPN). This allowed remote communication with the field
instrument from an office environment. The single board computer was interfaced with the
sensors using a Measurement Computing USB 4303 10 channel counter. A counter system was
used to overcome difficulties in transmitting low voltage analog signals over the long cables.
All of the analog signals needed to be converted to variable frequency digital pulse strings for
recording. A simplified block diagram of the temperature measurement and signal conversion
system is shown in Figure 4. There are five different functional blocks in the analog to the
digital circuitry. Regulated constant sources (Burr Brown REF200) feed a fixed 100 u.A current
into each of two thermistors used to measure temperature. The voltage generated is
proportional to the temperature at the thermistors reactive element. The reference thermistor
(Figure 4) is imbedded into the case of the flow sensor. The two voltages are subtracted with
-------
the difference amplifier and the voltage difference is converted to a frequency for transmission
to the frequency counter. A circuit diagram is shown in Appendix B and the printed circuit
board layout is shown in Appendix C. To minimize errors caused by temperature induced
buoyancy gradients, the flow sensor is installed level and the amount of energy injected into
the water is small. The above described circuit board provides eight channels of data and
requires one input power voltage (24 VDC) and two 5 VDC control signals. The USB 4303 used
as the interface between the signal data and the embedded computer has ten channels
available. An additional circuit board was constructed to support the USB 4303 and provide up
to three additional signals. Two of the additional signals are designed for the pressure
transducers, which are excited by a 15 VDC power supply, and one thermistor supplied a 100
u.A current source. The same voltage to frequency circuitry was used as discussed above. The
board layout is shown in Appendix D. Both circuit boards are shown in Figure 5 as they were
configured for the work at Fort Devens. Some of the components on the circuit boards are
labeled. In Figure 5, one of the circuit boards containing the USB 4303 is not fully populated
with components. These components were not installed, because they were not needed at
Fort Devens, and the power consumption needed to be reduced.
Figure 4 Block diagram of the measurement system.
-------
variable Rc&Ktors used
tacalibratoVFCs
Figure 5 Signal conditioning and data interface boards as used in the field.
-------
2.3 Flux Meter Mechanicals
The flow sensor is shown in Figure 6. The flow tube is constructed with Polyetheretherketone
(PEEK) material. The flow sensor is encased in a 304 Stainless Steel cylindrical flow tube case
shown in Figure 6A. In addition to containing the flow tube, the flow case provides a water
tight environment for the thermistor that measures the ambient water temperature and the
connections for the differential pressure transducer. These items are shown in the photo in
Figure 6A. The end plates drawings of the flow tube case are shown in Figure 6B and C. The
flow tube was held in the case using 1" 316 Stainless Steel compression fittings. The end plates
of the flow tube case were sealed with 0-rings. The differential pressure transducer was
encased in PVC. The two halves of the differential pressure transducer case were sealed with a
neoprene gasket. In the photo of Figure 6A, each pressure test port of the flow tube case and
the differential pressure transducer case is pointed out by an arrow. Both cases were pressure-
tested by immersing in water and pressurizing the case to 20 psig. The seals were considered
adequate if no air bubbles were observed.
The dome chamber consists of two separable pieces constructed out of 304 Stainless Steel. The
bottom cylinder is pressed into the sediment formation during installation in the field and an
upper dome houses the flow tube case and other components. The bottom cylinder isolates
2,922.5 cm2 of interfacial sediment area when pushed into the sediments. The dome chamber
also serves as an amplifier to the generally low seepage flow rates. A photo of the two parts
separated is shown in Figure 7A. Between the two halves there is a closed cell PVC gasket
providing a water tight seal. The two halves are held together by lever clamps. A flange is
welded around the lower part that serves as a gauge for installing the unit to a consistent depth
of 15 cm. The radius of the bottom cylinder is 30.5 cm. There are cross supports and a drive
index pin in the center for use in stiff sediments. The upper part of the dome has several ports
to connect the required fittings (Figure 7D). The 1/4" NPT fitting in the top is for the vent
degassing valve. The 2" NPT port in the top serves as an air relief port used during installation
to allow flow through the dome or a location where a standard addition can be performed. On
the left side of the figure there are three ports. The 1" NPT fitting is for the flow tube. The 1/4"
NPSF fitting is for the thermistor that measures the ambient temperature. The 1/4" NPT fitting
is for a pressure port for the differential pressure transducer. On the right side of the figure
there are two more ports. The 2" NPT fitting is for the cable going to the surface and the 1/4"
NPT fitting not labeled is an additional pressure port for the differential pressure transducer.
The unit is shown assembled in Figure 7B and C. Note that the flow tube is positioned inside
the dome for maximum protection from any potential outside destruction.
10
-------
• lurr>
Therm i
Figure 6 Flow tube and differential pressure transducer with drawings of flow tube case end
plates.
11
-------
Figure 7 Flow meter Dome Chamber. (A) Photo of dome halves separated; (B) Photo of dome
as it would be installed in the field; (C) Assembled unit with the flow tube and the differential
pressure transducer installed inside the chamber (shown from the bottom); (D) Drawing of
the dome chamber. A flange is welded around the bottom cylinder that serves as a gauge for
installing the unit to a consistent depth of 15 cm.
12
-------
3 Theory of Operation
The theory of heat flow has been a subject of investigation for centuries. Numerous papers and
books have been written on the subject. Probably the most comprehensive book on the
subject within the last century was Conduction of Heat in Solids written by H.S. Carslaw and J.C.
Jaeger originally published in 1946. When the developed theories are applied to physical
problems it is often necessary to approximate either the initial or boundary conditions. For a
simple solution that can describe the thermograph, one would either want an instantaneous
pulse of heat, a step of heat, or a square wave pulse. None of the normally assumed initial
conditions is a reasonable description of the heat input produced by our heater. As a result,
the shape of the thermograph observed in the physical system is inconsistent with the shape
produced by the common the mathematical developments. To overcome this problem,
Taniguchi and Fukuo, (1993) proposed using the time the peak temperature is detected at the
thermistors to determine the velocity of the heat. They took the derivative of an analytical
solution with ideal boundary conditions. The peak temperature occurs when the derivative
equals zero. Their development yields the equation:
Y2 — 'fY
U2 = ^H^- [2]
max
where U is the velocity, tmax is the peak temperature arrival time measured by the thermistor at
distance x from the heat source, and a is the thermal conductivity of water divided by the
specific heat and density. The water velocity is proportional but not equal to the thermal
velocity. If one can assume there is no thermal gradient in the radial direction (instantaneous
heat transfer in the radial direction) and there is no heat loss from the pipe, then the thermal
velocity will be slower and directionally proportional to the water velocity. The (Taniguchi and
Fukuo 1993) technique is suitable for measuring higher flow rates when the advective process
dominates. It uses the difference in peak arrival times between thermistors to estimate flow
rates. The precise distance between the thermistors must be known.
Lien (2006) proposed using moment analysis to determine the travel time of the heat-pulse. If
there are no other external inputs, the first moment of the thermograph is the heat-pulse peak
arrival time (tp). The time of travel of the heat-pulse between the center of the generated heat-
pulse and a measuring thermistor is equal to the distance between the points divided by the
sum of the diffusive flux and the advective flux, plus the time required for the heat to travel
from the heater to the water. Heat flow is controlled by two primary components. At low flow
rates the predominant component is heat conduction through the media. At high flow rates
13
-------
the primary component is advective flow of the water. The empirical function which describes
the flow of heat though the flow sensor has the form:
Q = ^ - [3]
t
p
where Q is the volumetric flow rate in the pipe and Ci; C2, C3 are calibration constants optimized
by mathematical fitting of data. The calibration constants are presumed to be dependent on
the physical geometry and material properties of each individual device.
14
-------
4 Laboratory Calibrations
There are three aspects of the instrument calibration. In the flux meter, thermistors are used
to measure temperature. Thermistors are non linear temperature dependent resistors and
need to be calibrated when precise measurements are needed. Since the heat-pulse
thermograph is defined by a series of data points with small temperature changes, all of the
thermistors in the flux meter need to be calibrated so small temperature changes can be
detected. Thermistors are analog devices and the signals must be "conditioned" before they
can be converted to digital pulse streams. The signal conditioning circuitry needs to be
calibrated because of the differences in components used in the construction of the circuitry.
Finally, the entire flow sensor needs to be calibrated because of the variability in the machining
and construction of the device.
4.1 Thermistor Calibration
Thermistors are temperature dependent resistors frequently used to measure temperature. All
of the thermistors used are commercially manufactured in stainless steel jackets designed for
continuous emersion in water. The nominal resistance, of all the thermistors, is 10,000 ohms at
25°C. Calibration is performed by measuring the resistance of the thermistor as a function of
temperature over the range of temperatures anticipated under field conditions (0 - 40°C). This
is accomplished by placing the thermistor along with a reference thermometer (Traceable
Digital Thermometer Model 4000CC or equivalent) in an aluminum block (Figure 8), and placing
the entire block in a temperature controlled water bath. The resistance of the thermistor is
Figure 8 Thermistors in an aluminum block used for calibration.
15
-------
measured with a four-digit digital multimeter (Fluke Model 85). A minimum of eight data
points are collected over a temperature range of 0 - 40°C. The least-squares method is then
used to fit the data to the Steinhart-Hart equation (Steinhart and Hart 1968):
= 4+ ?ln(p + ;[ln(p]2 [4]
where
• 7 is the temperature (in degree Kelvin)
• p is the resistance at T (in Ohms)
• A, B, and C are the Steinhart-Hart coefficients which vary depending on the type and
model of thermistor and the temperature range of interest.
The Steinhart-Hart coefficients for each thermistor are stored along with the serial number of
the thermistor. The calibration is considered successful when the mean squared error of the
measurements is less than or equal to l.OE-11, and the difference between measured and
calculated temperature does not exceed 0.5°C. When the mean squared error of the
thermistor is less than l.OE-11, the overall reproducibility is better than ± 0.2°C. Calibration is
performed in batches. The thermistors are then grouped or matched based on the resistance
value at 25°C. All thermistors used in each flux meter are from a same group.
All of the thermistors were calibrated as described above. Most of the thermistors followed the
expected response curve within the accuracy of the reference thermometer. Typical calibration
results are shown in Figure 9. The blue-diamond marks the measured temperature using the
reference thermometer versus estimated temperature using the thermistor. The red-cross
marks the error (secondary Y-axis) associated with each estimated temperature in degrees C.
The error appears to be small but random. In this particular example, the root mean squared
error (RMSE) of the estimator was 4.5E-3°C. There were no problems with the Steinhart-Hart
function. Calibration curves were, however, significantly different. The thermistors were
grouped such that the calculated resistances at 25°C of each calibration curves were as close to
the same value as possible. Thermistors within the same group were used in the same flux
meter.
16
-------
M
20 -
01
to 10
I
UJ 5
•
-
8.M-03
6.0E03
4.0E-03
2.DE-03
- 0-OE+UD of.
.2.0E.Q3 2
UJ
40E 03
-6,06-03
-8.0E-03
-l.OE-02
10
15
25
Measured Temperature (°C)
Figure 9 Example results of a typical thermistor calibration. The blue-diamond marks the
measured versus the estimated temperature. The red-cross marks the error (on secondary Y-
axis) associated with each estimated temperature. In this particular example, the root mean
squared error of the estimator was 4.5E-3°C.
17
-------
4.2 Circuitry Description and Calibration
The circuit diagram for the flux meter sensor, not including the data logger, is shown in
Appendix B. To measure the reference temperature, a non-adjustable precision reference
current source supplies a constant current (100 u.A half of Burr Brown REF 200, U15) to the
thermistor connected to J2. This thermistor is in series with the 1.82 kQ resistor (R32) the
produced voltage is buffered, gain of 1, by an operational amplifier (Burr Brown OPA134, U19)
to reduce the impedance and minimize stray current sinks. The output voltage from the
OPA134, pin 6, is proportional to a reference temperature (the temperature of the flow case)
and used as a reference for all of the temperature measurements in the flow path. The
reference temperature voltage is converted to a frequency by a voltage to frequency converter
U2 (an Analog Devices AD650), the output of U2 is connected to connector Jl pin 5 for
connection to the cable going to the interface board. Temperatures along the flow path are
measurements of difference. A measured temperature is subtracted from the reference
temperature and the difference is amplified by an instrumentation amplifier. For example the
thermistor nearest the cable end of the flow sensor is connected to J3 and terminal 2 is
connected to the negative input of the instrumentation amplifier U3 to make a temperature
measurement in the flow tube. The reference voltage (which is proportional to the case
temperature) is applied to the positive input of one of the instrumentation amplifiers (Burr
Brown INA 103). The output of the instrumentation amplifier is connected to the input of a
voltage to frequency converter (VFC) (Analog Devices AD650, U3). The voltage to frequency
converter circuit is designed to produce 100,000 pulses per second with a 10 VDC input. The
linearity of the circuit with the components used is specified by the manufacturer. The linearity
was verified during the design phase but is not re-checked during calibration. The output of the
VFC is transmitted by cable to the computer interface (Jl).
The VFCs require calibration. The process used for VFC U2 is to substitute a decade resistance
substitution box in place of the Thermistor at Jack J2 (see Appendix B for parts layout). The
resistance, approximately 8180 ohms, is adjusted until 1 VDC is produced at the output (pin 6)
of the OPA134 (U19) then Rl is adjusted until 10,000 Hz is produced at pin 8 of the AD650 (U2).
Once calibrated the adjusting screw on Rl is secured using a small drop of fingernail polish.
A similar approach is used on the circuitry for the temperature outside the dome connected to
J8 (see Appendix B). The decade resistance substitution box is connected to J8. The resistance
should be approximately lOkQ. The decade resistance substitution box is adjusted until 1VDC is
obtained at pin 6 of the OPA 134, U20. Then adjust R22 until 10,000 Hz is obtained at pin 8 of
theAD650, Ull.
The thermistors connected to J3 through J7 are measured on a differential circuit (the voltage
generated by the thermistors connected to J3 through J7 is sequentially subtracted from the
18
-------
reference voltage from OPA134, U19 and then multiplied by ten using an instrumentation
amplifier). Thermistors decrease in resistance with an increase in temperature. The thermistors
connected to J3 through J7 are at or above the temperature of the reference thermistor
connected to J2. Therefore, the output from the instrumentation amplifier will always be
positive. To calibrate the voltage to frequency converters, connect a resistance substitution
box to J2, as before, and adjust the output of the OPA 134, U19to 1 VDC. Short the pins of J3
and adjust R7 for 100,000 Hz at pin 8 of the AD650, U6. Repeat the process for the remaining
four temperature detectors.
The last input on the circuit board is the differential pressure transducer. Supply a known
reference voltage between 1 and 10 VDC to pin 2, and adjust R25 until the output from the
AD650, pin 8, of U12 is equal to 10,000 Hz per volt DC. Calibration is considered successful
when all adjustable frequencies are within 0.5% of the target test voltage or frequency.
On the circuit board there is also a connection for the heater and two digital control lines. One
digital control line operates the heater and the other control line operates the power supplies
used for the remainder of the circuitry. The power supplies are turned off between sampling
events to preserve power.
19
-------
4.3 Flow Sensor Calibration
The flow sensors were calibrated by pumping water through each unit at 7 different flow rates
bi-directionally from 0.27 to 16.7 ml/min (equivalent to Darcy velocities of 0.13 to 8.22 cm/day
in our current flux meter configuration). Calibrations were performed at room temperature.
The flow sensor was placed in a tank of water to minimize short term temperature fluctuation
but no attempt was made to maintain a constant temperature during the calibration process.
The measurements were repeated at least five times at each flow rate to develop the
calibration curves. Data was then fitted to Equation 3 using the least-squares method. For the
purpose of distinguishing flow direction in the flow tube, flow vector is assigned. Flow toward
T5 and T6 direction is considered positive while flow toward T3 and T2 direction is considered
negative flow. In the field the convention is to define positive flow as upward groundwater
discharge to surface water and negative flow as downward flow from surface water to
groundwater.
20
10
*c"
1 H
1
01 0 -
re
-5 -
-10 -
-IB -
-20
-C^T2 SN303
-•-T3SN303
-•-TSSN303
-O-T6 SN303
+ T3SN300
— I— T5 SN300
— X-J6 SN300
10
100 1000
Time of Peak Arrival (sec)
10000
Figure 10 Comparative plot of calibration curves for flow sensor unit SN300 and SN303. For
the purpose of distinguishing flow direction, flow toward T5 and T6 is considered positive
flow while flow toward T3 and T2 is considered negative flow.
20
-------
Figure 10 shows the calibration curves for flow sensor unit SN300 and SN303 (see Appendix E
for the calibration constants). The calibration curves of the two flow sensors are inherently
similar because of the similarity in their geometry and material used. To find out the difference
of the two flow sensors, a series of peak arrival times were applied to each calibration to
calculate flow rates. The results indicate the mean absolute error (MAE) between the two
sensors was 0.29 ml/min. The data points in Figure 10 represent the average of n replicate
measurements at the same flow rate, with n typically equal or greater than five. Taking all
replicate measurements into account the average coefficient of variation was 0.038 for SN300
and 0.037 for SN303. Figure 11 compares the estimated flow rate versus the measured flow
rate for T2 and T6 (the two outer thermistors) using Equation 3. It was perceptible that the
estimated and measured flow rates correlated very well. The error associated with each
estimated flow rate was shown on the secondary Y-axis. The root mean squared error (RMSE)
20 i
- 0.5
-0.5
.E • T2
E A T6
"p + Error T2
"^" + Error T6
2 — Linear (T2)
LU Linear (T6)
Measured Flow Rate {ml/min)
Figure 11 Comparison between measured and estimated flows for T2 and T6 using Equation
3. The error associated with each estimated flow rate is shown on the secondary Y-axis. The
RMSE of T2 and T6 flow estimation was 0.095 and 0.090 ml/min, respectively.
21
-------
20 -
.E • T3
^ A T5
~£ + T5 Error
^ + T3 Error
O
w — Linear (T3)
^ Linear (T5)
-0.5
R =0.9997
-20
- -1
Measured Flow Rate (ml/min)
Figure 12 Comparison between measured and estimated flow for T3 and T5 using Equation 3.
The error associated with each estimated flow rate is shown on the secondary Y-axis. The
RMSE of T3 and T5 flow estimation was 0.099 and 0.112 ml/min, respectively.
of T2 was 0.095 ml/min and the RMSE of T6 was 0.090 ml/min. The magnitude of error
appeared to be greater at lower flow rates. In a similar manner, Figure 12 compares the
measured flow rate and the estimated flow rate forTS and T5, which are the thermistors
closest to the heater. Similarly, the magnitude of error appeared to be greater at lower flow
rates. The RMSE of T3 was 0.099 ml/min and the RMSE forTS was 0.112 ml/min.
Figure 13 shows how to look at the same error data by plotting the MAE as a function of the
absolute value of the flow rate. The magnitude of error emerged to be greater at lower flow
rates. One potential explanation is based on the principals of operation. There are two
components to the heat flow. Advective flow dominates the time of travel of the heat-pulse at
high flow rates and conductive or diffusive flow dominates the time of travel of the heat-pulse
22
-------
at low flows. One might conclude that the majority of the error associated with the calibration
data presented here is in the estimation of the conductive flow. Another potential explanation
of greater errors under low flow situations is the high uncertainty in determining when the
peak temperature arrives. As the flow velocity slows down, the thermograph spreads and
looses its sharp peak. Lack of a sharp peak increases the uncertainty in determining the peak
arrival time. It appears that Figure 13 represents the behavior one normally encounters with
analytical measurements that have a lower limit of detection. Typically for the analytical
measurements, the error increases as the measured parameter approaches the limit of
detection of the analytical device.
I 0-2.2
I 2.2-8.8
8.8-16.7
0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16
Flow Rate (ml/min)
Figure 13 Mean absolute error in defining flow rate as a function of peak temperature arrival
time.
23
-------
5 Extended Time laboratory Study
A laboratory study was performed to evaluate the stability of the measurements over an
extended period of time. Flow sensor SN302, without an external temperature detector or
differential pressure transducer, was used to determine the repeatability of the measurements.
The pump supplied a constant flow through the sensor. Measurements were repeated 296
times at a flow rate of 2.6 ml/min equivalent to Darcy velocities of 1.28 cm/day in our current
configuration. There was no attempt made to control the ambient temperature beyond the
controlled room temperature in the laboratory. Figure 14 demonstrates the long-term
measurement repeatability in the laboratory. The mean travel time to the first thermistor
down gradient from the heater (T5) was 104.8 seconds with a coefficient of variation of 0.013.
The mean travel time to the second thermistor down gradient from the heater (T6) was 188.3
seconds with a coefficient of variation of 0.022. Based on the calibration curve these travel
times would be equivalent to Darcy velocities of 1.36 and 1.24 cm/day with a coefficient of
variation of 0.2. The error is less than 6.5% compared to the known pumping rate. As
mentioned above the error increased with lower flow rates. The errors observed here are not
different than those observed during calibration.
"o"
Si 200
-------
6 Field Deployment
Field testing of the long-term flux meter operation was carried out continuously from June 12,
2008 through April 7, 2009. The site selected was located in the Red Cove area of Plow Shop
Pond adjacent to the Shepley's Hill Landfill Superfund site at Fort Devens, Massachusetts
(Figure 15). Detailed site information was described in the report by Ford, Acree et al. (2009);
as noted, groundwater was discharging to Red Cove under most conditions based on the
potentiometric surface data and seepage flux measurements with the prototype flux meters.
Note that the flux measurements within this report are specific to the area within Red Cove in
which the automated flux meters were installed. The studies presented in this report were not
designed to estimate seepage flux at the scale of Red Cove.
Figure 15 Aerial photo of the Red Cove area of Plow Shop Pond, Fort Devens, Massachusetts.
The yellow dots mark the location of the instruments; the blue square marks the location of
STAFF1 Red Cove staff gauge.
25
-------
Figure 16 Placement of a perforated PVC mat prior to the flux meter installation. The
perforated mat (4' x 4' in size, 1/4" thick with 1/4" hole-diameter) was placed on the
sediment with a cutout in the center for the flux meter.
Sediments in this portion of the pond have a thick organic layer of material developed from
falling leaves and branches. The surface sediment has very little load bearing capacity. To
support the weight of the instrument, a perforated PVC mat (4'x 4', 1/4" thick, 1/4" hole-
diameter, 1/2" centers straight rows) was placed on the sediment with a cutout in the center
for the flux meter (Figure 16). A perforated mat was used rather than a solid mat to minimize
divergence of the water flow lines around the mat that could impact the measurements. While
the mat design was configured to minimize impact to sediment conditions, it should be noted
that the extent of the sediment compaction due to the placement of the mat and how it could
be a factor in the observed discharge rates at this location was not evaluated as part of this
study. Two flux meter units (SN300 & SN303) along with four piezometers (PZ5S, PZ5C, PZ5G &
PZ5N) were installed in Red Cove as shown in Figure 17. The flux meters were deployed side by
side for comparison. The meters were placed on separate pads but the pads were placed
adjacent each other resulting in a center to center spacing of approximately 4 feet (1.2 m). The
flux meter easily penetrated the sediment and was supported by the flange on the cylinder.
The Southern flux meter was SN300 and the Northern meter was SN303. The three
26
-------
piezometers were placed along the west side of the pads and labeled PZ5N, PZ5C and PZ5S for
the North, Center, and South piezometer, respectively. All three piezometers were installed at
a 4 ft (1.2 m) vertical depth below the sediment. Staff gauge well PZ5G (bundled with PZ5C)
was placed above the sediment surface to measure water stage data. Water lilies covered most
of the cove at the time of field installation; the diver had to clear out the underwater
vegetation before installing the flux meter units. The installation took place during the
turnover of the reservoir and the diver reported near zero visibility underwater. There was no
significant horizontal flow in this water body; the Bernoulli Effect was considered negligible.
Thus there was no need to use the differential pressure transducer to assist in the placement of
the dome. The water depth at the time of installation was 3.46 ft (105.5 cm). The instrument
installation was completed on June 12, 2008 and the first data was collected at 9:28 PM
Greenwich Mean Time (GMT). Monitoring continued until April 7, 2009 at 2:13AM GMT.
PZ5N [+42°33'16.Sa",-71°35'41.34"]
SN303 [442°33'16.74", -71°35'41.10"]
PZ5G, PZ5C ]l42°33'16.6S",-71°35'40.86"]
- SN300 [442°33'16.56", -71°35'41.34"]
_•
PZ5S [+42°33'16.56", -71°35'41.10"]
Figure 17 Deployed locations [latitude, longitude] of the flux meter and piezometer. Water
lily covered most of the cove at the time of installation.
Daily meteorological data for the meteorological station at the Fitchburg Municipal Airport
were obtained from the National Climatic Data Center of the National Oceanic and Atmospheric
Administration. It is noted that the station is located approximately twelve miles from the
study site. Although the exact magnitude of meteorological data at the Fitchburg station may
27
-------
be different from that in the Red Cove, the meteorological patterns should be sufficiently
similar to represent the conditions at the site. Figure 18 shows the time trend (June 2008 -
April 2009) of the water temperature measured at the SN300 flux meter location and the
average air temperature at Red Cove. The majority of the time the water and air temperature
at the site was below 10°C during the sampling period, and was way below room temperature
range during the calibration of the flow sensor.
Air Temperature » Water Temperature
01
O.
0)
iO
zs
»_ 15
QJ
3
E 10 -
i)
-5
-10
-15
20
Figure 18 Time trend of water temperature and average air temperature at the Red Cove of
Plow Shop Pond, Fort Devens, Massachusetts. Air temperature data were obtained from the
Fitchburg Municipal Airport meteorological station, National Climatic Data Center of the
National Oceanic and Atmospheric Administration.
28
-------
6.1 Solar Power System
There was no power source available at the site. Therefore, a solar power system was installed
to provide continuous electrical power for the long-term flux meter operation (Figure 19). The
solar panel was installed uphill, south of Red Cove, for a better sunlight exposure. The solar
system was designed for 300W peak power production. The average solar hours in the area is
about 4.2 (Marion and Wilcox 1992), which should yield an annual average production of
approximately 52W. The system would consume an average of 22W, if the flux meter
measurement operated 4 times a day. The battery backup was designed to maintain power for
up to 10 days. Initially there was no monitoring of the power system status after installation.
In early July and late August, the voltage regulator/battery charger abruptly failed. After a
second failure of the regulator, it was decided to begin monitoring the battery voltage. In
Figure 20, a seven point moving average of the battery voltage is plotted along with a seven
point moving average of the sky cover. A moving average was used to smooth out a significant
amount of apparent randomness in the data. On the sky cover scale a zero is full sun and a 10
is 100% cloud cover. Although logic would suggest that when the sky cover went up the battery
voltage might go down, there does not appear to be much of a correlation between the two.
One possible explanation might be related to when the battery voltage was measured and how
Figure 19 Solar power system was installed uphill south of Red Cove.
29
-------
bright the sun was shining at the time of the measurement. The battery voltage was measured
once per sampling cycle. The sampling frequency was not constant throughout the experiment.
The sampling frequency was reduced during the winter to minimize power consumption. It
should be noted that after the battery voltage was monitored there were no more problems
with the regulator and the solar panel performed as expected. From a qualitative point of view
the power produced was somewhat less than anticipated. While no setbacks were
encountered after resolving the problem with the regulator, it was felt that the sampling
frequency needed to be reduced to reduce power demand during the winter months.
30 -
28 -
26
24 -
> 22 -
01
I 20
> 18 -
16 -
14 -
12 -
10
fN
00
rx
CO
fN
CO
CM
r-
-------
6.2 Heater Response
For reliable results the heater should produce a consistent amount of heat at the time of each
measurement cycle. In the laboratory this was not an issue, because a regulated power supply
was used rather than a battery, and very little variability was noticed in the measured peak
temperature of the heater. In the field, the irregular battery voltage seemed to affect the
consistency of the rise in the temperature of the heater. From the data collected between
10/10/2008 and 4/7/2009 (the battery voltage data collection did not begin until October
2008), there is a weak but apparent correlation (coefficient of correlation = 0.34) between the
rise in temperature of the heater and the battery voltage (see Figure 21). In Figure 21, the
correlation was influenced by a small fraction (11%) of data (associated with random dates over
the measurement period), which appear to have a distinctive lower rise in temperature. The
trends in the two groups of data appear to be similar. If we regard the data with a distinctive
lower rise in temperature as outliers, the coefficient of correlation improves significantly to
0.82. This proves that the consistency of the rise in the temperature of the heater was
compromised by the voltage irregularity. Although there was no indication that the correlation
had a bad effect on the Red Cove data quality (based upon evaluating the shape of the
thermograph as well as the peak arrival time), the power supplied to the heater should be
regulated in future designs to prevent any doubts in the field.
u
\r Z-
-------
6.3 Pressure Transducer Response
In the field portion of the research at Red Cove, an attempt was made to monitor the
groundwater-surface water flux. There are several approaches that can be used to estimate
this flux. If one knows the hydraulic conductivity of the sediment one could monitor the
vertical hydraulic gradient and calculate the flux. In a similar manner if you know the flux and
vertical hydraulic gradient you can calculate the hydraulic conductivity. At Red Cove,
piezometers were installed at four feet below the sediment-water interface. Pressure
transducers (SETRA 5241-010wd-G-Wl-2C-300-F) were installed in the piezometers and
gauging station at common depths to record the vertical hydraulic heads electronically. Due to
the faulty voltage regulator of the solar panel, the flux meter system was out of commission
from 7/9/2008 to 7/24/2008 and from 8/26/2008 to 10/9/2008. There was no piezometric
data available during these time periods. In both incidents the system regained its full
functionality after the faulty voltage regulator was replaced.
STAFF! • P2SS PZSC PZSN • Precipitation
LO
66.7
66.5
66.3
> 66.1
12
I
65.7
65.5
«t -
'
•/*
.•
O
SI
CTi
rj
cri"
X
fM
8.0
- 7.0
- 6.0
I
Q.
- 3.0 'D
£
-• 2.D
1.0
0.0
Figure 22 Comparative plots of water elevations (primary Y-axis) at piezometers PZ5S, PZ5C,
PZ5N and Red Cove STAFF1 gauge station and time-series precipitation (secondary Y-axis).
The daily precipitation data of Red Cove area was obtained from the Fitchburg Municipal
Airport meteorological station, National Climatic Data Center of the National Oceanic and
Atmospheric Administration.
32
-------
0.12 -i
£ 0.10 H
•M 0.08 -
_0)
"ro °-06
5
.y 0.04 H
3
ro
•n 0.02 -
03 0.00 -
u
^ -0.02
-0.04
CO
o
fM
rH
"$•
CO
O
O
fM
00
o
o
rM
CO
o
o
I
00
o
o
r^
at
fM
CO
O
O
CM
CT1
O
O
fM
fJl
O
O
^
CO
fM
ft
* PZ5S
A PZ5C
* PZ5N
o
o
fM
Figure 23 Temporal vertical hydraulic gradients at each piezometer location. Temporal
variation of the three piezometers appeared to follow similar trend. The gradient peaked in
August 2008 and then gradually descended over time.
In the field the vertical hydraulic gradient was manually determined three times, once on
6/13/2008, once on 8/20/2008 and once on 4/7/2009. The vertical hydraulic gradient is
determined by the static elevation in the piezometer minus the static elevation of the water
level in the pond, divided by the vertical depth of piezometer screen below sediment. For PZ5S,
the vertical hydraulic gradients were 0.030, 0.033, and 0.030, respectively. For PZ5C, the
vertical hydraulic gradients were 0.044, 0.041, and 0.040, respectively. For PZ5N the the
vertical hydraulic gradients were 0.028, 0.024, and 0.026, respectively. Each time the manual
measurements were made the pressure transducers readings were calibrated accordingly.
Based on the pressure transducer data and field observation, it was determined the PZ5G staff
gauge had physically slipped downward several feet into the sediment during winter time,
which caused the water stage data to drift upward gradually over time. The exact time and rate
33
-------
of the occurrence were unverifiable. As a result, the PZ5G gauge data was proven unreliable
and was excluded for any further use. Instead, the surface water elevation data obtained at the
Red Cove STAFF1 gauge station were used. The water elevation of the piezometer was
determined by the static head in the piezometer plus the surface water elevation at the Red
Cove STAFF1 gauge. Figure 22 shows comparative plots of water elevations at piezometers
PZ5S, PZ5C and PZ5N and STAFF1 gauge and time-series precipitation (secondary Y-axis). The
daily precipitation data of Red Cove area was obtained from the Fitchburg Municipal Airport
meteorological station, National Climatic Data Center of the National Oceanic and Atmospheric
Administration. Temporal and spatial variations were evident as shown in Figure 22. Vertical
hydraulic gradients were calculated for each piezometer location based on the water elevation
data and the vertical depth of the piezometer below the sedimets. All three piezometers were
installed at a vertical depth of 4 feet (1.22 m) below the sediment. Figure 23 shows the
majority of the time the direction of the vertical hydraulic gradients were upward and their
magnitude fluctuated daily. There is no verifiable reason why the vertical hydraulic gradient for
PZ5N during the time span between 7/21/2008 and 9/9/2008 diverges so much from the other
two piezometers. During the other times, all three piezometers seem to show consistent result.
Due to natural heterogeneity, it is not unexpected to see spatial variability in the piezometer
data. Figure 23 illustrates minor spatial heterogeneity does exist in the vertical hydraulic
gradient. Fifty two percent (52%) of the daily piezometric data have a spatial coefficient of
variation smaller than 25%. Overall, 79% of the data have a daily spatial coefficient of variation
smaller than 50%. Temporal variation appeared to follow similar trend for the piezometers.
The vertical hydraulic gradient descended initially, ascended rapidly, peaked in August 2008,
and then gradually descended over time. The gradient climax in August could be a cumulative
hydrological response to the continual rainfall event from July 31 to August 15, 2008 (Figure 22)
with a total precipitation of 10.8 cm. The apparent precipitation event from September 6 to
September 14, 2008 yielded 12.6 cm (Figure 22). Because of the system stoppage, as
mentioned above, there was no piezometric data recorded during that time. The hydrological
response at this time was not verifiable. Nonetheless, Figure 22 clearly shows noticeable
correlation between the rise in the surface water elevation and the precipitation event.
34
-------
6.4 Flux Meter Response
Unexpectedly, the data collected from the field test was not consistent with the calibrations.
The initial crude observation indicated that the advective flux values were extremely low. But
many of the observed peak temperature arrival times in the field were longer than the peak
temperature arrival times of the zero advection (or zero flow) condition observed during the
laboratory test. This suggested that the calibration curves could not adequately describe the
field data when the flow was extremely low. In the literature, the influence of very low
temperature on the seepage meter measurements has been observed (Mwashote, Burnett et
al. 2010). It was suspected that this temperature effect might have also played a role in our
field measurements, given the fact that the average water temperature was 7.8 ± 5.1°C during
the field test, which was much lower than the room temperature condition during calibration.
After completing the field deployment, the flow sensors were retrieved and tested in the
laboratory to investigate the phenomenon. The test was conducted at zero advection while the
temperature of the water bath, where the flow sensor was placed was allowed to drift. Ice was
used to cool the temperature of the water bath initially, and then the water bath was allowed
to gradually return toward room temperature over time. Numerous cycles of measurements at
u
01
D
4-1
TO
OJ
0.
E
.2
18
16 -
14 -
12
10
8
6
4 -
2 -
n
AU
A •
A •
A •
A •
A •
• A
X
T3
AT5
500 1000
Peak Arrival Time (sec)
1500
Figure 24 Effect of ambient temperature on the movement of heat-pulse in zero flow
condition showing increasing peak arrival time with decreasing temperature.
35
-------
zero flow were repeated. The temperature of the water bath and the peak arrival time at each
cycle were recorded. In Figure 24, preliminary results show the effect of ambient temperature
on the travel time of a heat-pulse under a zero flow condition. It is evident that the heat-pulse
movement becomes sluggish as the water temperature decreases. Although it is reasonable to
believe the temperature effect would become much less significant when advection dominates,
this effect could have profound influence on field measurements when advection is low. Figure
24 illustrates the need to understand flow sensor performance based on knowledge of the
temperature-dependent thermal conductivity of water and the PEEK material.
u
Oi
•P
u
O
U
TO
01
1.6 -|
1.58
1.56 -
1.54 -
1.52
1.5 -
1.48 -
1.46
1.44
1.42
1.4
Y = 0.0031 X + 1.40
4- Data
Function
0
10 20 30 40 50 60
70
Temperature ( C)
Figure 25 Martin and Lang's Thermal conductivity of water as a function of temperature.
Martin and Lang (1933) measured the thermal properties of water and demonstrated the
thermal conductivity is temperature dependent (Figure 25). Lower thermal conductivity at
lower temperatures is consistent with longer travel times at lower temperatures. Martin and
Lang's data is almost linear over the temperature range evaluated, which shows decreasing
thermal conductivity with decreasing temperature. This is consistent with our observations of
approximately linear data points in Figure 24 with increasing peak temperature arrival time
when the ambient water temperature decreases. However, the thermal properties of the PEEK
material used to construct the flow tube are not the same as water. Local thermal equilibrium
36
-------
between the liquid and solid phases may not always be maintained. Runyan and Jones (2008)
observed temperature-dependent thermal conductivity of PEEK at very low temperatures
between 0.3 and 4 degree Kelvin. But at the present time, the literature source(s) for the
temperature-dependent thermal conductivity of PEEK in the temperature range of this study
have not been identified. Therefore, we can not use Martin and Lang's data alone to forecast
the flow sensors response. Instead, an alternative approach will be used to establish an
empirical relationship between the variables in the heat-pulse transfer.
37
-------
6.4.1 Alternative Approach
At the location of the automated flux meter installation within Red Cove, the flow was limited
to very low rates. While these low flow rates may not be significant from a total volume of
flow, they may be important from an ecological standpoint for the sediment's flora and fauna,
if the flow contains toxic chemicals. As mentioned above, the errors appear to increase as the
flow rate decreases. At very low flows the temperature peak spreads over time and it becomes
difficult to accurately determine when the peak occurs. As a result, the estimated flow rates
become unreliable at low flows. To improve the resolution, an alternative method similar to
the approach proposed by (Wang, Ochsner et al. 2002; Ochsner, Morton et al. 2005) for soil
water flux determination was considered. At very low flows heat can be detected on both sides
of the heater. At zero flow, if the flow tube is constructed symmetrically, the peak temperature
that is observed should be the same at equal distances from the heater center. Likewise, if
there is a very small flow in one direction the higher peak will occur in the direction of the flow.
For flow rates where the heat can accurately be detected on both sides of the heater, one can
estimate both the magnitude and direction of the flow from the ratio of the temperature of
equally spaced temperature detectors. The form of the equation describing flow (Q) is:
(T "i
i-l [5]
where 3 is a coefficient and TH and TL are the downstream and upstream temperature rises.
Note that (Wang, Ochsner et al. 2002; Ochsner, Morton et al. 2005) define the term 3 as a
function of the thermal conductivity of the porous media, the heat capacity of water and the
distances from the heater to the upstream and the downstream temperature sensors. The
heater and the temperature sensors were three parallel needles aligned in a commom plane,
which is different from the configuration of our heat-pulse device. Because of the different
boundry condition, we can not use the same expression in (Wang, Ochsner et al. 2002;
Ochsner, Morton et al. 2005) to define p. Instead, a mathematical fitting of data was used to
establish the empirical expression of p.
Intuitively, if there is no water flow, the temperature rise in each of the two temperature
detectors spaced at an equal distance from the heat source should be the same. Likewise, if
there is a very small amount of flow, the convection of heat by the flowing water resulted in a
larger temperature increase downstream from the heat source than upstream from the heat
source. With this concept in mind, the ratio of the rise in temperature downstream to the rise
38
-------
in temperature upstream as a function of flow rate was evaluated. The natural log of the ratio
appears to be directly proportional to the flow rate at a constant temperature. In Figure 26,
the response is shown for the data collected at two different temperatures. The temperature
effect is significant. As the flow rate approaches zero, the importance of the thermal
conductivity dominate over the advective processes. Based on the work of Martin and Lang it is
assumed the coefficient 3 in Equation 5 is linearly correlated to the temperature. Using the
data in Figure 26, we determined by mathematical fitting of the data, the linearly temperature
dependent 3 = 0.0301 + 0.0287 T, where T is the temperature in Celsius. Figure 27 shows a
comparative plot of measured flow rate versus temperature corrected estimated flow rate,
when attempting to estimate flow rates up to 0.3 ml/min through the flow tube over a
temperature range of 5 to 15°C. The estimated flow rate follows Q = (0.0301 + 0.0287T) In
(TH/TL), where the flow rate Q is in ml/min and the ratio TH/TL is the temperature rise of the
downstream thermistor divided by the temperature rise of the upstream thermistor. T is the
ambient water temperature in Celsius. The flow direction is determined by the thermistor with
the higher temperature rise. The RMSE of the estimated flow rate was 0.013 ml/min,
equivalent to a Darcy velocity of 0.006 cm/day in our current dome cylinder configuration. This
is adequate for the measurements being made. It is known that over the flow range studied,
temperature affects are significant, and the lower the flow rate the more significant the effect.
Equation 5 should only be applied to low flow situations where advective heat flow no longer
dominates. For the data set collected at Red Cove, Equation 5 modified with a thermally
dependent coefficient has provided a valid way for interpreting the data. The empirical
equation is assumed universally applicable to both flow sensors for the reason of their near
identical geometry in the construction.
0.28
A 6 degree C
• IS degree C
- 0.191 n(TH/TL)
- 0.471 n (TH/TL)
-0.02
U1(T«/TJ
Figure 26 Ambient temperature effect on estimated flow rate.
39
-------
0.35
0.3 -
• 0.25
0.2 H
0.15
I
y = 0.9943x
R2 = 0.9861
0.05 0.1 0.15 0.2 0.25 0.3 0.35
Measured Flow Rate (ml/min)
Figure 27 Comparative plot of measured flow rate versus temperature corrected estimated
flow rate. The RMSE of the estimated flow rate was 0.013 ml/min.
40
-------
6.4.2
Flux Meter Results
Monitoring was continuous at various frequencies over the time period except when there
were problems with the solar panel's voltage regulator. Sampling frequency was
reprogrammable from any remote location with Internet access. Data were downloaded daily
for analysis. During the winter months when the water surface was ice-covered, data could not
be interpreted with certainty.
Temporal advective flux data for SN300 and SN303 excluding the period of January 17 - March
8, 2009 is shown in Figure 28. The daily sampling frequency was not constant; each data point
represents a daily average value. Piezometric data was collected once per measurement cycle.
The field data shows a significant amount of scattering, which was not seen under controlled
laboratory conditions of constant flow. Data in Figure 23 suggests that there were fluctuations
in the vertical hydraulic gradients, which could be the cause of the variability observed in the
measured seepage flux.
ro
I"
X
3
2.0 n
1.5 -
i.o -
0> 0.5
'•8
V
•o o.o -
-0.5
00
§
§
o
00
§
1
00
§
I
O
O
00
o
o
(N
00
o
o
rM
cT
00
o
o
o~
01
o
o
to
rN
OSN303
OSN300
(Tl
O
O
r^
00
IN
o
o
r\i
Figure 28 Daily average advective flux data at the Red Cove of Plow Shop Pond, Fort Devens,
Massachusetts. The data were corrected according to the ambient water temperature.
41
-------
T 120%
-• 100%
-• 80%
-• 60%
40%
•• 20%
-»- Cumulative %
-0.5 0 0.5 1 1.5 2 More
Advective Flux (mm/day)
Figure 29 Histogram of advective flux measured at Red Cove flux meter location SN300.
12Q -
100 -
80 -
60 -
40 -
20 -
h .
-
•
-• m
...m.->
- 120%
- 100%
- 80%
- 60%
- 40%
- 20%
. n*.
• Cumulative %
•0.5 0 0.5 1 1.5 2
Advective Flux (mm/day)
More
Figure 30 Histogram of advective flux measured at Red Cove flux meter location SN303.
42
-------
Due to the natural heterogeneity, it is expected that the vertical hydraulic gradient at the flux
meter locations were, to some extent, different from the piezometer locations. The gradients
at piezometers were low with sporadic negative values at times. All of the data suggest low
flow and the majority of the data suggest an upward flow direction from the groundwaterto
the surface water. Data also superficially suggests a trend that the flow magnitude was
declining over time. The average advective flux, in terms of Darcy velocity, for SN300 was 0.44
± 0.34 mm/day (median = 0.36 mm/day; minimum = -0.34 mm/day; maximum = 1.81 mm/day).
Figure 29 is the histogram of the flux measured at SN300. It is evident that the advective flux
was extremely low with 77% of the measurements less than 0.5 mm/day and 93% less than 1
mm/day. The average advective flux for SN303 was 0.76 ± 0.33 mm/day (median = 0.77
mm/day; minimum = -0.45 mm/day; maximum = 1.89 mm/day). Figure 30 is the histogram of
the advective flux measured at SN303. Flow measured at SN303 was also extremely low but
slightly higher than flow measured at SN300, with 89% less than Imm/day and 97% less than
1.5 mm/day. The data scattered appears to be less in the spring than during the entire time.
The average spring (3/9/2009 - 4/7/2009) advective flux for SN300 and SN303 was 0.42 ± 0.16
mm/day and 0.38 ± 0.16 mm/day, respectively. Spatial heterogeneity in the advective flux does
exist, eventhough the two flux meters were placed only 1.2 m apart center to center.
i.z -
ra" i.o -
•a
"^
t
.§, 0.8 -
X
3
u.
QJ 0-6 -
>
B
<" nyl
> 0.4 -
3
00
> 0.2 -
— .
p
—
—
L
.U "I — I — — I — — I — — I — — I
00 00 00 00 00 00
O O O O O CD
c -i on a. f! >
3 .5 3 a; £! o
-i < u~> O z
r
I
00 O) O} tTi
O O O O
U c J
a t
a; ro Q, >u
LJ ^^ t i ^>
~~
- U.8
- 0.7
,_,
- 0.6 E
u
*JT DSN303
" °'5 .2 DSN300
£ D Precipitation
- 0.4 '£_
01
- 0.3 ^
00
-0.2 I
- 0.1
n r>
~i r u.u
a^
o
\.
CL
Figure 31 Monthly average advective flux and precipitation.
43
-------
Figure 31 illustrates the relationship between monthly average advective flux (primary Y-axis)
and precipitation (secondary Y-axis). From observations, there appears to be a qualitative
association between the flux and the precipitation. Flux was higher during the summer and
early fall when precipitation was also relatively higher. Flux was lower during winter and early
spring when precipitation was lower. The average precipitation in September 2008 was
elevated, but unfortunately no flux data was available during that time due to the faulty solar
power system. The elevated flux values in October could very well be the result of the large
precipitation in September.
There were problems analyzing the winter data, from January 17, 2009 to March 8, 2009, when
uncharacteristic thermographs resembling wave-driven flow were noticed. Both flux meters
experienced the same phenomenon. The wave influenced artifact was similar to the
thermographs observed by Lien (2006) at Lake Hartwell. It was noticed that the surface water
of Red Cove was covered with ice during the winter months. Figure 32A shows a representative
data set for the winter months. Thermograph of the heater and the two thermistors T3 and T5
nearest to the heater are shown. Figure 32A suggests that the flow in the flow tube was
changing direction back and forth; the thermograph of T3 and T5 appears to be conversing with
each other indicating continuous alternation in flow direction. The oscillation frequency was
about 3.5 minutes. It was suspected this oscillating pattern is comparable to the effect of
wave-induced flow observed by (Smith, Herne et al. 2009). Although with waves this is
intuitive, when the surface is frozen there should not be any waves. Nevertheless, if there were
oscillating instability of the overlaying ice mass it could potentially create a pumping effect in
the ice-covered water akin to waves. The data was unexpected and we did not collect pressure
gradient data near the sediment interface to substantiate alternation of the pressure gradients
to be the cause of the oscillating flow. Uncertainty of meaningful data interpretation could be
attributed to the winter data (from January 17, 2009 to March 8, 2009) being omitted from the
data analysis. Figure 32B is the type of smooth thermographs observed after the ice was
thawed and is typical of much of the remainder of the data, which proves the uncharacteristic
winter data was not caused by any instrument malfunction. It is safe to conclude, the flux
meters were functional throughout the field deployment.
44
-------
-1000-O.OS 0 1000
Time (sec)
1.8
1.6
1.4
1.2
1
0.8
0.6
O.4
P
0.2
-
0
2000-0.2
ro
GJ
00
c
TO
(J
OJ
•«-»
2
01
0.
£ Ft Devens
2/9/2009
• T3
• T5
• Heater
(A)
0.18
0.16
0.14
0.12
0.1
0.08
0.06
0.04
0.02
\
1
;
•
•1000 0 1000
Time (sec)
1.8
1.6
1.4
1.2
1
0.8
0,6
0.4
0.2
0
c3
ro
OJ
o
CUD
C
to
£
3
JD
i—
O>
P
Ft
Devens
4/6/2009
o"
P
• T3
. T5
• Heater
2000
(B)
Figure 32 Apparent flow reversals during measurement cycle. Winter processed data (A),
spring processed data (B).
45
-------
6.5 Flux Estimation using Darcy's Law
Another way of estimating seepage flux is to apply Darcy's law (Hubbert 1956) using hydraulic
conductivity and piezometric data. Darcy's law is expressed in the following equation:
9
q= - =
A
L
[6]
where q is flux (L/T), Q is flow rate (L3/~T), A is area (L2), K is hydraulic conductivity (L/T), and
(AH/L) is the vertical hydraulic gradient. To obtain the hydraulic conductivity, a falling-head
permeameter test was conducted on April 8, 2009 adopting the same approach by (Chen 2004;
Genereux, Leahy et al. 2008). A permeameter was installed closely adjacent to the flux meters.
The permeameter was installed on April 7, 2009, which allowed the disturbance of the
formation to settle before taking measurements the next day. Figure 33 is a schematic of the
field permeameter test. The test involves inserting a pipe (permeameter) vertically into the
sediments, filling the pipe with water, and measuring the rate at which the water level inside
the pipe falls; hydraulic conductivity of the section of sediments isolated is then calculated from
this rate.
Schedule 80 PVC Pipe
v
y
u
Surface water
Sediments
Figure 33 Falling-head permeameter test to measure sediment hydraulic conductivity.
46
-------
Chen (2004) proposed a simplified equation in the form of Equation 7 to estimate vertical
hydraulic conductivity Kv:
Kv - --In [7]
where Lv is the length of the sediment column; H0 and Htare the hydraulic heads inside the pipe
corresponding to time t0 and t. Chen (2000) suggested the error in estimate is less than 5%
when Lv/D > 5 with D denoting the inner diameter of the pipe. The error can be reduced by
increasing the ratio of LV/D. Equation 7 could be derived in the form of Equation 8 (Genereux,
Leahy et al. 2008). Hydraulic conductivity Kv is calculated from the slope of the best fit line.
\nHt= - -^t+*H0 [8]
Our permeameter is a schedule 80 PVC pipe 14.4 cm in inner diameter and 305 cm in length.
The bottom edge of the pipe was shaped in an outward bevel to help minimize sediment
compaction during installation. The permeameter was driven 91.4 cm into the sediment with a
hammer. Figure 34 shows the results of a falling-head permeameter test. Hydraulic
conductivity was calculated from the slope of the trend line. The results show at this particular
location that the hydraulic conductivity is 0.031 m/day, which is a typical value for silt.
Temporal advective flux was then calculated based on the hydraulic conductivity and the
average vertical hydraulic gradient of the three piezometers. Similar to the measured values
the calculated flux values were small. The average calculated flux was 1.12 ± 0.71 mm/day
(median = 1.15 mm/day). Comparison between the measured and calculated flux is shown in
Figure 35. Although the permeameter and flux meters were closely adjacent they were not at
the exact same location. Therefore, spatial variability of the flux value should be anticipated.
Despite the inherent spatial heterogeneity the calculated and measured flux were convincingly
comparable, which also confirmed the validity of the two methods. Results from both
approaches suggest a very low flow, and the majority of the data suggests an upward flow
direction from the groundwaterto the surface water. The magnitude of flow compare
favorably in the early summer of 2008 and in the spring of 2009. During the fall and winter
months the calculated flux values were higher than the measured values, but nevertheless
remain in the same order of magnitude as the measured values.
47
-------
(J
4.247
4.246 j
4.245 -
4.244
4.243
4.242 -
4.241 -
4.240
4.239 -
4.238 -
4.237
0
-3.925E-07X +4.245
R2= 0.982
5000
10000
time (sec)
15000
20000
Figure 34 Falling head permeameter test at the Red Cove of Plow Shop Pond, Fort Devens,
Massachusetts. Hydraulic conductivity is proportional to the slope of the best fit line.
3.0 -,
*> 2.5 -
flj
5. 2.0 -
E
E 1.5 -
X
.3 1.0 -
LJ_
> °'5^
QJ 0.0 -
< -0.5
-1.0
00
8
s
00
o
o
00
o
o
00
o
o
r\l
cn
00
o
o
r\l
CTi
O'
00
o
o
00
r-l
rsl
OSN303
OSN300
Calculated
o
o
fN
tNI
CTi
O
o
(N
00
rM
pri
en
O
O
Figure 35 Comparison between measured and calculated advective flux. The calculated flux
was based on hydraulic conductivity and vertical hydraulic gradient measurements.
48
-------
7 Summary and Discussion
In Situ long-term monitoring of seepage flux across the sediment-water interface is necessary,
for two reasons. Seepage flux is temporally heterogeneous in nature, it is central to collect and
analyze repeated measurements to evaluate changes in condition, and the long-term
monitoring provides a key to evaluating in situ performance of the flux meter. For
contaminated sediment sites, monitoring temporal (as well as spatial) variation of advective
flux is essential to proper risk management. This report presents the development and field
testing of a ruggedized advective flux meter capable of unattended remote monitoring of
sediment-water interfacial flux. The flow sensor employs a heat-pulse technique; the bi-
directional feature of the flux meter is realized through the temperature measuring capability
on either side of the flow tube. In the field operation, flow across the sediment-water interface
is funneled through a dome to a flow tube; the rate of water flowing through the flow tube is
measured. The advective flux, in terms of a vertical Darcy velocity, is calculated by dividing the
flow rate through the flow sensor by the cross-sectional area of the dome. The dome serves as
an amplifier to the generally low Darcy flow rates. A remote controlled monitoring system is
included to permit system monitoring and modification over the internet. The sampling
frequency is changeable through the wireless network. Data is downloadable from any remote
location with Internet access after each measurement cycle.
The system has automatic data acquisition, real time data display, and signal display/analysis
capability. The instrument has undergone several calibrations to establish relations between
flow rate and heat-pulse travel time. Under typical condition when advection dominates, time
of peak temperature arrival is inversely propotional to the flow rate. For very low flow rates
when advection no longer dominates, the spread of the heat-pulse made it difficult to
determine precisely when the peak arrived. A sensible approach using temperature rise ratio
was adopted to provide a creditable alternative for determining flow rate. Using this approach
it is possible to accurately interpret flow rates down to zero flow conditions. Depending on the
magnitude, flow rate can be derived as a function of peak temperature arrival time or as a ratio
in the rise of temperature. Using these two methods, reliable measurements can be made over
a wide bi-directional range of flow rates.
The extended laboratory study of flux meter performance demonstrated the ability to repeat
measurements accurately. In the field, unattended remote monitoring of seepage flux
measurement was demonstrated at Fort Devens. The automated operation ran continuously
for ten months with a self-sufficient power supply system. The only problem encountered in
the field was with a faulty solar panel voltage regulator. The system performed well after the
glitches were resolved. Two flux meters were placed side by side in shallow water. They
successfully monitored the water flux from June 2008 to April 2009. The flux meter is capable
49
-------
of quantitatively measuring very low flow rates. All of the data suggest very low flow at the flux
meter locations and the majority of the data suggest an upward flow direction from the
groundwater to the surface water. The average advective flux for flux meter unit SN300 and
SN303 was 0.44 ± 0.34 mm/day and 0.76 ± 0.33 mm/day, respectively. The results show that
temporal flow variation exists over the study period and the flow at the location of the flux
meter was very low. Hydraulic conductivity and piezometric measurements were not part of
the flux meter. They were intended to provide supporting evidence for the anticipated
variability in the flux measurements. Calculated seepage flux based on hydraulic conductivity
and gradient data suggests the average flux was 1.12 ± 0.71 mm/day. Note that these
measurements are not reflective of measured fluxes in other locations of the cove receiving
contaminated groundwater discharge. The studies presented in this report were not designed
to estimate seepage flux at the scale of Red Cove.
The flux meter works well when the water surface is reasonably still. The apparent oscillating
flow during the winter months is perplexing. It was suspected that wave actions due to
oscillating disturbances in ice-covered water were the cause of the flow oscillation during the
winter period. More field studies are needed to determine if the observation made during this
study is a common occurrence, and to determine the cause of the observation. At the present
time there are no reliable data interpretation methods for handling flow oscillation or wave
action. Efforts need to be made in developing techniques that will permit interpreting data that
is generated during apparent flow oscillation. Measuring the differential pressure across the
sediment-water interface may provide data needed to assist in the interpretation of the flow
behavior. The apparent oscillating flow is potentially significant. It would potentially allow a
major cyclic fluid exchange between the sediment and the overlaying water affecting the
biogeochemical processes in permeable beds (Precht and Huettel 2003).
A sharp peak is essential when the analysis is based on the time of the peak arrival. It is difficult
to determine the precise time of the peak arrival when the heat-pulse spreads in a low flow
condition. A few modifications could be made to improve the data interpretation. An
additional temperature sensor could be installed closer to the heat source. This would reduce
the time available for the heat-pulse to spread before it was detected. The length of the heater
could be reduced making the heat source more of a point source. The material of the heater
could be changed to one with a high thermal conductivity and low heat capacity. The duration
of the heat-pulse could be shortened focusing the energy closer to an instantaneous point
source. Sufficient heat is needed to be able to easily detect the change in temperature, but the
amount of heat injected needs to remain small to minimize the buoyancy effects caused by the
heating of the water.
50
-------
Thermistors were used to measure temperature in this device. They provided excellent
resolution and accuracy. However, there was a considerable amount of variability among the
thermistors. The time and labor involved in generating individual calibration curves for each
thermistor was costly. It would be desirable to find an interchangeable temperature detector
that is sufficiently small and have the required sensibility that could be used to replace the
thermistors. Thermocouples are interchangeable, but have low sensitivity and did not offer a
proper signal-to-noise ratio needed for this application. A suitable replacement has not been
identified at this point in time.
51
-------
8 References
Campbell, G. S., C. Calissendorff, et al. (1991). "Probe for measuring soil specific heat using heat-pulse
method." Soil Sci Am J 55: 291-293.
Chen, X. (2004). "Streambed Hydraulic Conductivity for Rivers in South-Central Nebraska." Journal of the
American Water Resources Association (JAWRA) 40(3): 561-573.
Choi, H.-m., K.-B. Lee, et al. (1991). "Measurement of gas flowrate by using a heat pulse." Flow
Measurement and Instrumentation 2(3): 188.
Ford, R. G., S. Acree, et al. (2009). "Devens 2008 Monitoring Update: Arsenic Fate, Transport and
Stability Study, Groundwater, Surface Water, Soil and Sediment Investigation, Fort Devens
Superfund Site, Devens, Massachusetts." EPA/600/R-09/064.
Fritz, B. G., D. P. Mendoza, et al. (2009). "Development of an electronic seepage chamber for extended
use in a river." Ground water 47(1): 136-140.
Genereux, D. P., S. Leahy, et al. (2008). "Spatial and temporal variability of streambed hydraulic
conductivity in West Bear Creek, North Carolina, USA." Journal of Hydrology 358: 332-353.
Hubbert, M. K. (1956). "Darcy's law and the field equations of the flow of underground fluids." Am. Inst.
Min. Met. Petl. Eng. Trans. 207: 222-239.
Kalbus, E., F. Reinstorf, et al. (2006). "Measuring methods for groundwater,surface water and their
interactions: a review" Hydrol. Earth Syst. Sci. Discuss. 3: 1809-1850.
Lee, D. R. (1977). "Device for Measuring Seepage Flux in Lakes and Estuaries." Limnology and
Oceanography 22(1): 140-147.
Lien, B. K. (2006). "Development and Demonstration of a Bi-directional Advective Flux Meter for
Sediment-Water Interface." EPA/600/R-06/122.
Marion, W. and S. Wilcox. (1992). "Solar Radiation Data Manual for Flat-Plate and Concentrating
Collectors." from http://rredc.nrel.gov/solar/pubs/redbook/.
Martin, L. H. and K. C. Lang (1933). "The Thermal Conductivity of Water." Proc. Phys. Soc. 45(4): 523-
529.
Mortensen, A. P., J. W. Hopmans, et al. (2006). "Multi-functional heat pulse probe measurements of
coupled vadose zone flow and transport." Advances in Water Resources 29(2): 250.
Murdoch, L. C. and S. E. Kelly (2003). "Factors affecting the performance of conventional seepage
meters." Water Resources Research 39(6): SWC 2-1 : 2-10.
Mwashote, B. M., W. C. Burnett, et al. (2010). "Calibration and use of continuous heat-type automated
seepage meters for submarine groundwater discharge measurements." Estuarine, Coastal and
Shelf Science 87: 1-10.
Ochsner, T. E., R. Horton, et al. (2005). "Evaluation of the Heat Pulse Ratio Method for Measuring Soil
Water Flux." Soil Sci. Soc. Am. J. 69: 757-765.
Paillet, F., R. Crowder, et al. (1996). "High-Resolution Flowmeter Logging Applications with the Heat-
Pulse Flowmeter." Journal of Environmental and Engineering Geophysics 1(1): 1-11.
Paulsen, R. J., C. F. Smith, et al. (2001). "Development and evaluation of an ultrasonic ground water
seepage meter." Ground Water 39(6): 904-911.
Precht, E. and M. Huettel (2003). "Advective pore-water exchange driven by surface gravity waves and
its ecological implications." Limnol. Oceanogr. 48(4): 1674-1684.
Ren, T., G. J. Kluitenberg, et al. (2000). "Determining soil water flux and pore water velocity by a heat
pulse technique." Soil Sci Am J 64: 552-560.
Rosenberry, D. O. (2008). "A seepage meter designed for use in flowing water." Journal of Hydrology
359(1-2): 118-130.
Rosenberry, D. O. and R. H. Morin (2004). "Use of an electromagnetic seepage meter to investigate
temporal variability in lake seepage." Ground Water 42(1): 68-77.
52
-------
Runyan, M. C. and W. C. Jones (2008). "Thermal conductivity of thermally-isolating polymeric and
composite structural support materials between 0.3 and 4 K." Cryogenics 48(9-10): 448-454.
Schincariol, R. A. and J. D. McNeil (2002). "Errors with small volume elastic seepage meter bags." Ground
Water 40(6): 649-651.
Shinn, E. A., C. D. Reich, et al. (2002). "Seepage Meters and Bernoulli's Revenge." Estuaries and Coasts
25(1): 126-132.
Sholkovitz, E., C. Herbold, et al. (2003). "An automated dye-dilution based seepage meter for the time-
series measurement of submarine groundwater discharge." Limnology and Oceanography-
Methods 1: 16-28.
Smith, A. J., D. E. Herne, et al. (2009). "Wave effects on submarine groundwater seepage measurement."
Advances in Water Resources 32(6): 820-833.
Steinhart, J. S. and S. R. Hart (1968). "Calibration curves for thermistors." Deep Sea Res. 15: 497-503.
Swarzenski, P. W. and J. A. Izbicki (2009). "Coastal groundwater dynamics off Santa Barbara, California:
Combining geochemical tracers, electromagnetic seepmeters, and electrical resistivity."
Estuarine Coastal and Shelf Science 83(1): 77-89.
Taniguchi, M., W. C. Burnett, et al. (2006). "Submarine groundwater discharge measured by seepage
meters in Sicilian coastal waters." Continental Shelf Research 26(7): 835-842.
Taniguchi, M., W. C. Burnett, et al. (2003). "Spatial and temporal distributions of submarine
groundwater discharge rates obtained from various types of seepage meters at a site in the
Northeastern Gulf of Mexico." Biogeochemistry 66(1-2): 35-53.
Taniguchi, M. and Y. Fukuo (1993). "Continuous measurements of ground-water seepage using an
automatic seepage meter." Ground water 31(4): 675-679.
Taniguchi, M., T. Ishitobi, et al. (2007). "Evaluating ground water-sea water interactions via resistivity
and seepage meters." Ground Water 45(6): 729-735.
Taniguchi, M. and H. Iwakawa (2001). "Measurements of submarine groundwater discharge rates by a
continuous heat—type automated seepage meter in Osaka Bay, Japan." Groundwater Hydrology
42(271-277).
Wang, Q., T. E. Ochsner, et al. (2002). "Mathematical analysis of heat pulse signals for soil water flux
determination." Water Resour. Res. 38(6): 1091.
Zuber, A. (1970). "Method for determine leakage velocities through the bottom of reservoirs." Isotope
Hydrology: 761-771.
53
-------
Appendix A Visual Basic Code
frmDataAcquisition
'Option Explicit' last revision 5/19/08 by Carl Enfield
' Purpose To measure temperatures, differential temperatures and differential pressures
' Scans a range of frequencies Input Channels, stores the data in an array
' and output the data to ASCII files.
'send digital output to switch heater ON/OFF
'sends digital signal to control the power supply on the flux meter
'the temperature measurements are made with thermistors. The thermistors
' are supplied a constant current of 100 uA. For channels 1 and 7 the impedance
' is changed by supplying the voltage to the positive input of an operational
' amplifier which serves as a buffer amplifier. The output from this amplifier from channel 1 is
' used as a reference For the other temperature measurements
' the voltage signals are input to the negative input of instrumentation amplifiers
'channel 1: reference junction temperature a 1.82K ohm 1% resistor is in series with this thermistor
'channel 2 - 3: flow temperatures beginning from cable end
'channel 4: heater temperature CA adhesive is initially used to attach the thermistor to the heater.
'Epoxy (JB weld) is then applied for added strength.
'channel 5-6: flow temperature
'channel 7: temperature outside the dome
'channel 8: Differential pressure transducer
'channels 2 through 6 are subtracted from channel 1 and have a gain often all other channels have
'a gain of one.
Const BoardNuml% = 0 ' Board number Site 1
Const BoardNum2% = 1 ' Board number Site 2
'Dim BoardNum2%
Const FirstPoint& = 0 'set first element in buffer to transfer to array
Const Warmup = 12 ' Time in minutes after activating circuit board power before measurements are initiated
Const Gatelnterval& = 100 ' Number of milliseconds used to count pulses and determine frequency
Const Gatelntervall& = 1000 ' Number of milliseconds used to count pulses and determine frequency used on ref
temperature only
'Dim BoardNum%
Dim ActualCounts%(10)
Public baseline 'time used to collect temperature data prior to injecting heat
Dim CBCount%(10)
Dim ChipNum%
Dim counter%, rounds, RateCount
Dim DataValue% ' digital value to write to Auxport
Dim Dirl As String * 32 ' Directory for data storage for Site 1
Dim Dir2 As String * 32 ' Directory for data storage for Site 2
'Dim DirA As String * 32
Dim DomeDia ' diameter of dome in centimeters
Dim fname As String * 50 ' file name for raw data with ext .dat
Dim fnameS As String * 50 ' file name for supplemental raw data with ext .dat
Dim fnameP As String * 50 ' file name for processed output data with extension .csv
Dim fnamel As String * 50 ' file name for processed relative change in temperature
54
-------
Dim fname2 As String * 50 ' file name for data converted to temperature in deg C
Dim FOutDivider&
Dim Freq&(10)
Public heater 'time in seconds heater is turned on
Dim i, k, j, N As Integer ' counter
Dim MemHandle& ' define a variable to contain the handle for
' memory allocated by Windows through cbWinBufAlloc%()
Public pause ' time in seconds before next cycle of data collection allos heat to dissipate
Public CycleTime
Dim SampleTime, startpause
Dim sitel As String * 64 ' a descriptive string describing the site 1
Dim site2 As String * 64 ' a descriptive string describing the site 2
Dim Source&
Public tail ' time in seconds experimental data is collected
Dim CallStop, mytime
Dim finish, timepast, t, timezero
Dim sum(8) As Double
Dim ULStat As Long 'universal library statistic return
Dim SeriaINo
Dim SigSource(lO)
Dim textline
Dim Pass%
Dim PassWarm As Long ' number of passes through pause cycle
Dim Resistor(), clO(), c20(), c30(), cvOQ, cll(), c21(), _
c31(), cvl(), c!3(), c23(), c33(), cv3(), c!4(), c24(), _
c34(), cv4(), SerNoQ, TslopeQ, Tint(), TSnlQ, chla(), _
chlbQ, chlc(), TSn2(), ch2a(), ch2b(), ch2c(), TSn3(), _
ch3a(), ch3b(), ch3c(), TSn4(), ch4a(), Ch4b(), Ch4c(), _
TSn5(), ch5a(), Ch5b(), Ch5c(), TSn6(), Ch6a(), Ch6b(), _
Ch6c(),TSn7(), Ch7a(), Ch7b(), Ch7c()
Private Sub CmdClose_Click()
Call PowerOffHeatOff(BoardNuml%)
lfCheck2Then
'MsgBox ("Box Checked")
CallPowerOffHeatOff(BoardNum2%)
End If
Unload Me
End Sub
Private Sub PowerOffHeatOff(BoardNum%)
' Turn off heater turn off power supply
ULStat = cbDBitOut(Boardl\lum%, AUXPORT, 0, 0) ' turn off power supply flux meter
If ULStat <>0 Then Stop
ULStat = cbDBitOut(BoardNum%, AUXPORT, 1, 0) ' turn off heater
If ULStat <>0 Then Stop
ULStat = cbDBitOut(Boardl\lum%, AUXPORT, 2, 0) ' turn off power auxiliary sensors
If ULStat <>0 Then Stop
' MsgBox (BoardNum%)
55
-------
End Sub
Private Sub PowerOnHeatOff(BoardNum%)
' Turn off heater turn on power supply
ULStat = cbDBitOut(Boardl\lum%, AUXPORT, 0, 1) ' turn on power supply
If ULStatoO Then Stop
ULStat = cbDBitOut(Boardl\lum%, AUXPORT, 1, 0) ' turn off heater
If ULStatoO Then Stop
ULStat = cbDBitOut(Boardl\lum%, AUXPORT, 2, 1) ' turn on power supply auxiliary sensors
If ULStatoO Then Stop
1 MsgBox (BoardNum%)
End Sub
Private Sub cmdStart_Click()
i*********
' Collect Site Specific Screen Data
sitel = CStr(txtSitel.Text)
site2 = CStr(txtSite2.Text)
SeriaINo = CStr(txtSerNo.Text)
SerialNo2 = CStr(txtSerNo2.Text)
baseline = CDbl(txtBaseline.Text) ' time to collect baseline data (min)
heater = CDbl(txtHeater.Text) ' time for heater to be turned on (sec)
tail = CDbl(txtTail.Text) ' time to collect data after heater is turned off (min)
DomeDia = CStr(txtDDia.Text)
'MsgBox ("Read screen data")
Dirl = "C:\Temp\Sitel"
Dir2 = "C:\Temp\Site2"
'fname = Format$(Now, "yyyy") + ".txt"
'myFile= Dir("C:\\Temp\*.*")
'MsgBox ("looked for file in C:\\Temp")
'lf(myFile = "")Then
1 ChDir(Dirl)
' Open fname For Output As #1
1 Write #1, Format("test")
1 Close #1
1 MkDir "C:\\Temp\"
1 ChDir("C:\\Temp\")
' Open fname For Output As #1
1 Write #1, Format("test")
1 Close #1
'End If
'myFile = Dir("C:\\Temp\Sitel\*.*")
'MsgBox ("looked for file in C:\\Temp\Sitel")
'If (myFile = "") Then
' MkDir "C:\\Temp\Sitel"
1 ChDir(Dirl)
' Open fname For Output As #1
56
-------
1 Write #1, Format("test")
1 Close #1
'End If
'myFile = Dir("C:\\Temp\Site2\*.*")
'MsgBox ("looked for file in C:\\Temp\Site2")
'If (myFile = "") Then
' MkDir "C:\\Temp\Site2"
1 ChDir(Dir2)
' Open fname For Output As #1
1 Write #1, Format("test")
1 Close #1
'End If
Pass% = 1
IblTime.Caption = Now 'display time
Labell2.Caption = Pass%
Dirl = "C:\Temp\Sitel"
Dir2 = "C:\Temp\Site2"
ChDir(Dirl)
fnameP = Format$(Now, "mmddhh") + "p.CSV"
Open fnameP For Output As #10
Write #10, Format("DateTime"), Format("baseline"), Format("heater"), Format("Tail"), Format("Serial No."), _
Format("Dome Dia"), Format("Site")
Write #10, Format(Now), Format(baseline, "#0.00"), _
Format(heater, "#0.00"), Format(tail, "#0.00"), Format(SerialNo, "#0"), _
Format(DomeDia, "#0.00"), Format(sitel)
Write #10, Format("DateTime"), Format("PeakT2"), Format("PeakT3"), Format("PeakT5"), Format("PeakT
6"),_
Format("Moment 2"), Format("Moment 3"), Format("Moment 5"), Format("Moment 6"), Format("Ref
Temp"), _
Format("SupTemp"), Format("Diff Pres"), Format("Ch 9"), Format("Ch 10"), _
Format("PeakRefTemp"), Format("PeakTemp2"), Format("PeakTemp3"), Format("PeakTemp4"),
Format("PeakTemp5"), _
Format("PeakTemp6")
Close #10
lfCheck2Then
ChDir(Dir2)
fnameP = Format$(Now, "mmddhh") + "p.CSV"
Open fnameP For Output As #10
Write #10, Format("DateTime"), Format("baseline"), Format("heater"), Format("Tail"), Format("Serial No."),
Format("Dome Dia"), Format("Site")
Write #10, Format(Now), Format(baseline, "#0.00"), _
Formatfheater, "#0.00"), Format(tail, "#0.00"), Format(SerialNo2, "#0"), _
Format(DomeDia, "#0.00"), Format(site2)
Write #10, Format("DateTime"), Format("PeakT2"), Format("PeakT3"), Format("PeakT5"), Format("Peak
T6"),_
Format("Moment 2"), Format("Moment 3"), Format("Moment 5"), Format("Moment 6"), Format("Ref
Temp"),_
57
-------
Format("SupTemp"), Format)" Diff Pres"), Format("Ch 9"), Format("Ch 10"), _
Format("PeakRefTemp"), Format("PeakTemp2"), Format("PeakTemp3"), Format("PeakTemp4")
Format("PeakTemp5"), _
Format("PeakTemp6")
Close #10
End If
Do While CallStop <> -1 'main loop
Labell2. Caption = Pass%
' Re-read timing values each pass so system can be modified during run
baseline = CDbl(txtBaseline.Text) ' time to collect baseline data (min)
heater = CDbl(txtHeater.Text) ' time for heater to be turned on (sec)
CycleTime = CDbl(txtCycleTime.Text) 'total time for a cycle in hours
tail = CDbl(txtTail.Text) ' time to collect data after heater is turned off (min)
SampleTime = baseline * 60 + heater + tail * 60 'total sampling time in seconds
pause = CycleTime * 60 - (SampleTime / 60) - Warmup ' time in minutes
If Check2 Then pause = CycleTime * 60 - (2 * (SampleTime / 60)) - 2 * Warmup
If pause < 0 Then MsgBox ("Insufient Pause time. Pause in minutes = " & pause)
'MsgBox ("pause= " & pause)
' IblSiteNum. Caption = 1
Refresh
Call Warm_Up(BoardNuml%)
Call DataAcquisition(l, Dirl, sitel, BoardNuml%, Pass%, fnameP, SeriaINo, baseline, heater, CycleTime, tail,
SampleTime, DomeDia) ' Collect Data Set for Site 1
' MsgBox (tail)
Call PlotData(fname, fnameP, fnameS, tail, 1, Dirl, Pass%, 1) ' Plot Site 1 Data
lfCheck2Then
' LblSiteNo. Caption = 2
Refresh
Call Warm_Up(BoardNum2%)
Call DataAcquisition(2, Dir2, site2, BoardNum2%, Pass%, fnameP, SerialNo2, baseline, heater, CycleTime, tail,
SampleTime, DomeDia) ' Collect Data Set for Site 1
Call PlotData(fname, fnameP, fnameS, tail, 2, Dir2, Pass%, 2) ' Plot Site 2 Data
End If
Pass% = Pass% + 1
If Pass% > 100 Then
Pass% = 1
ChDir(Dirl)
fnameP = Format$(Now, "mmddhh") + "p. CSV"
Open fnameP For Output As #10
Write #10, Format("DateTime"), Format("baseline"), Format("heater"), Format("Tail"), Format("Serial No."]
Format("Dome Dia"), Format("Site")
Write #10, Format(Now), Format(baseline, "#0.00"), _
Format(heater, "#0.00"), Format(tail, "#0.00"), Format(SerialNo, "#0"), _
Format(DomeDia, "#0.00"), Format(site)
58
-------
Write #10, Format("DateTime"), Format("PeakT2"), Format("PeakT3"), Format("PeakT5"), Format("Peak
T6"),_
Format("Moment 2"), Format("Moment 3"), Format("Moment 5"), Format("Moment 6"), Format("Ref
Temp"), _
Format("SupTemp"), Format("Diff Pres"), Format("Ch 9"), Format("Ch 10"), _
Format("PeakRefTemp"), Format("PeakTemp2"), Format("PeakTemp3"), Format("PeakTemp4"), _
Format("PeakTemp5"), Format("PeakTemp6")
Close #10
lfCheck2Then
ChDir(Dir2)
fnameP = Format$(Now, "mmddhh") + "p.CSV"
Open fnameP For Output As #10
Write #10, Format("DateTime"), Format("baseline"), Format("heater"), Format("Tail"), Format("Serial
No."), Format("Dome Dia"), Format("Site")
Write #10, Format(Now), Format(baseline, "#0.00"), _
Format(heater, "#0.00"), Format(tail, "#0.00"), Format(SerialNo2, "#0"),_
Format(DomeDia, "#0.00"), Format(site2)
Write #10, Format("DateTime"), Format("PeakT2"), Format("PeakT3"), Format("PeakT5"),
Format("PeakT6"),_
Format("Moment 2"), Format("Moment 3"), Format("Moment 5"), Format("Moment 6"), Format("Ref
Temp"), _
Format("SupTemp"), Format("Diff Pres"), Format("Ch 9"), Format("Ch 10"), _
Format("PeakRefTemp"), Format("PeakTemp2"), Format("PeakTemp3"), Format("PeakTemp4"), _
Format("PeakTemp5"), Format("PeakTemp6")
Close #10
End If
End If
txtCycleTime.BackColor = &HFFOO&
Refresh
CallPowerOffHeatOff(BoardNuml%)
lfCheck2Then
'MsgBox ("power off heat off")
CallPowerOffHeatOff(BoardNum2%)
End If
t = 0
timezero = Timer
' MsgBox ("pause =" & pause)
Do While t< pause* 60 'Loop 14
If Timer < timezero Then 'if cross midnight
t = 86400 - timezero + Timer
Else
t = Timer- timezero
End If
IblNxtCycTime.Caption = Format(((pause -1 / 60) + Warmup), "#0.0")
IblTime.Caption = Now 'display time
DoEvents 'Yield to other processes
If CallStop =-1 Then GoTo 44
Loop 'Loop 14
59
-------
1 Pass% = 1
txtCycleTime.BackColor = &HFFFFFF
Loop 'main
44 txtBaseline.BackColor = &HFFFFFF
txtHeater.BackColor = &HFFFFFF
txtTail.BackColor = &HFFFFFF
txtCycleTime.BackColor = &HFFFFFF
cmdStart. Enabled =-1
CallStop = 0
GoTo 45
'Loop
45 End Sub
Private Sub Warm_Up(BoardNum%)
Call PowerOnHeatOff(BoardNum%)
' timezero_main = Timer
' Warm up required before data valid. 20 min determined adequate to obtain reasonably stable values
txtBaseline.BackColor = &HFFFF&
Screen. MousePointer = 11 'hourglass pointer
Refresh
t = 0
timezero = Timer
PassWarm = 0
Do While t < Warmup * 60
PassWarm = PassWarm + 1
'calculate accmulated time
DoEvents
If Timer < timezero Then 'if cross midnight
t = 86400 - timezero + Timer
Else
t = Timer -timezero
End If
If PassWarm = 64000 Then ' reduces the number of times display is refreshed
IblNxtCycTime. Caption = Format((Warmup - (t / 60)), "#0.0") ' warm up time remaining
Refresh
PassWarm = 0
End If
DoEvents
Loop
End Sub
Private Sub DataAcquisition(site#, DirA, siteA, BoardNum%, Pass%, fnameP, SeriaINo, baseline, heater, CycleTime,
tail, SampleTime, DomeDia)
'Dim DirA As String* 32
ChDir(DirA)
60
-------
I*********
' Initialize signal sources
' There are five inputs on each chip so use chip 2 for the second 5 signal sources
SigSource(l) = CTRINPUT1 ' Input one chip 1
SigSource(2) = CTRINPUT2 ' Input two chip 1
SigSource(S) = CTRINPUT3 ' Input three chip 1
SigSource(4) = CTRINPUT4 ' Input four chip 1
SigSource(S) = CTRINPUT5 ' Input Five chip 1
SigSource(6) = CTRINPUT6 ' Input six chip 2
SigSource(7) = CTRINPUT7 ' Input 7 chip 2
SigSource(S) = CTRINPUT8 ' Input 8 chip 2
SigSource(9) = CTRINPUT9 ' Input 9 chip 2
SigSource(lO) = CTRINPUT10 ' Input 10 chip 2
' Collect Screen Input Data
MemHandle& = cbWinBufAlloc(lOO)
If MemHandle& = 0 Then Stop
'write baseline, heat-pulse duration, SeriaINo, and Site description raw data file #1
Date = Now
Time = Now
fnameS = Format$(Date, "mmdd") + Format$(Time, "hh") + "S" + ".dat"
fname = Format$(Date, "mmdd") + Format$(Time, "hh") + ".dat"
'fnamel = Format$(Date, "mmdd") + Format$(Time, "hh") + "l.dat"
'fname2 = Format$(Date, "mmdd") + Format$(Time, "hh") + "2.dat"
Open fname For Output As #1
Write #1, Format(Now), Format(baseline, "#0.00"), _
Format(heater, "#0.00"), Format(tail, "#0.00"), Format(SerialNo, "#0"), _
Format(DomeDia, "#0.00"), Format(siteA)
Close #1
ULStat = cbWinBufFree(MemHandle&) ' Free up memory for use by
If ULStat <> 0 Then Stop ' other programs
' Initialize settings required to scan thermistors
' Chanels 1 through 6 are scanned.
MemHandle& = cbWinBufAlloc(lOO)
If MemHandle& = 0 Then Stop
Screen. MousePointer = 0 'Default Pointer
Refresh
Date = Now
Time = Now
IblSiteNum. Caption = site#
Refresh
61
-------
timezero = Timer
t = 0
k = 0
Do While t (baseline * 60) Then txtBaseline.BackColor = &HFFFFFF
If (t > (baseline * 60) And t <= (baseline * 60 + heater)) Then DataValue% = 1: txtBaseline.BackColor = _
&HFFFFFF: txtHeater.BackColor = &HFFOO&
If t > (baseline * 60 + heater) Then DataValue% = 0: txtHeater.BackColor = &HFFFFFF: txtTail.BackColor =
&HFFOO&
'digital out
ULStat = cbDBitOut(Boardl\lum% AUXPORT, 1, DataValue%)
If ULStatoO Then Stop
DataSet (BoardNum%) ' Collect a set of data
Open fname For Append As #1
'write data in frequency
Write #1, Format(Now), Format(t, "#0.00000"), _
Format(Freq&(l), "#0"), Format(Freq&(2), "#0"),_
Format(Freq&(3), "#0"), Format(Freq&(4), "#0"), _
Format(Freq&(5), "#0"), Format(Freq&(6), "#0")
Close #1
'display data in frequency
For i = !To6
IblADData(i). Caption = Format(Freq&(i), "#0")
Nexti
IblTime.Caption = Now 'display time
IblRemainTime. Caption = Format((SampleTime - 1) / 60, "#0.0")
Refresh
DoEvents
If CallStop = -1 Then GoTo 44
Loop 'Loop 11
'48 Close #1
DataSet2 (BoardNum%) 'collect suplimental data
For i = 6 To 10
IblADData(i). Caption = Format(Freq&(i), "#0")
Next i
Open fnameS For Output As #9
Write #9, Format(Now), _
62
-------
Format(Freq&(7), "#0"), Format(Freq&(8), "#0"), Format(Freq&(9),
Format(Freq&(10), "#0"), Format(SerialNo, "#0"), Format(site)
Close #9
txtTail.BackColor = &HFFFFFF
44 txtBaseline.BackColor = &HFFFFFF
txtHeater.BackColor = &HFFFFFF
txtTail.BackColor = &HFFFFFF
txtCycleTime.BackColor = &HFFFFFF
cmdStart. Enabled =-1
CallStop = 0
GoTo 45
45 End Sub
Private Sub CmdStop_Click()
CallStop = -1
Call PowerOffHeatOff(BoardNuml%)
If Check2Then Call PowerOffHeatOff(BoardNum2%)
Unload Me
End Sub
Private Sub cmdStopConvert_Click()
'MsgBox("cmdStop")
Call PowerOffHeatOff(BoardNuml%)
'MsgBox ("Board loff")
lfCheck2Then
CallPowerOffHeatOff(BoardNum2%)
'MsgBox ("Board 2 off")
End If
ULStat = cbWinBufFree(MemHandle&) ' Free up memory for use by
' other programs
1 If ULStat <>0 Then Stop
End
Unload Me
End Sub
Private Sub Form_Load()
' declare revision level of Universal Library
ULStat = cbDeclareRevision(CURRENTREVNUM)
' Initiate error handling
' activating error handling will trap errors like
' bad channel numbers and non-configured conditions.
' Parameters:
' PRINTALL :all warnings and errors encountered will be printed
' DONTSTOP :if an error is encountered, the program will not stop,
63
-------
' errors must be handled locally
ULStat = cbErrHandling(PRINTALL, DONTSTOP)
If ULStatoO Then Stop
1 If cbErrHandling% is set for STOPALL or STOPFATAL during the program
' design stage, Visual Basic will be unloaded when an error is encountered.
' We suggest trapping errors locally until the program is ready for compiling
' to avoid losing unsaved data during program design. This can be done by
' setting cbErrHandling options as above and checking the value of ULStat%
' after a call to the library. If it is not equal to 0, an error has occurred.
End Sub
Private Sub lnitialize_Click()
Call PowerOffHeatOff(BoardNuml%)
If Check2Then Call PowerOffHeatOff(BoardNum2%)
End Sub
Private Sub lnitlUL(BoardNum%)
' declare revision level of Universal Library
ULStat& = cbDeclareRevision(CURRENTREVNUM)
' Initiate error handling
' activating error handling will trap errors like
' bad channel numbers and non-configured conditions.
' Parameters:
' PRINTALL :all warnings and errors encountered will be printed
' DONTSTOP :if an error is encountered, the program will not stop,
' errors must be handled locally
ULStat& = cbErrHandling(DONTPRINT, DONTSTOP)
If ULStat&<>0 Then Stop
1 If cbErrHandling& is set for STOPALL or STOPFATAL during the program
' design stage, Visual Basic will be unloaded when an error is encountered.
' We suggest trapping errors locally until the program is ready for compiling
' to avoid losing unsaved data during program design. This can be done by
' setting cbErrHandling options as above and checking the value of ULStat%
' after a call to the library. If it is not equal to 0, an error has occurred.
' Initialize the board level features
' Parameters:
' BoardNum& :the number used by CBCONFIG to describe this board
1 ChipNum& :chip to be initialized (1 for CTR5,1 or 2 for CTR10)
1 FOutDivider&:the F-Out divider (0-15)
' Source& :the signal source for F-Out
' Comparel& :status of comparator 1
' Compare2& :status of comparator 2
' TimeOfDay& :time of day mode control
ChipNum% = 1
64
-------
FOutDivider& = 1 ' sets up OSC OUT for 10kHz signal which can
Source& = FREQ1 ' be used as frequency source for this example
Comparel&= Disabled
Compare2& = Disabled
TimeOfDay& = Disabled
ULStat& = cbC9513lnit(Boardl\lum%, ChipNum%, FOutDivider&, Source&, Comparel&, Compare2&,
TimeOfDay&)
If ULStat&<>0 Then Stop
End Sub
Private Sub lnit2UL(BoardNum%)
' Initialize the board level features
' Parameters:
' BoardNum& :the number used by CBCONFIG to describe this board
1 ChipNum& :chip to be initialized (1 for CTR5, lor 2 for CTR10)
1 FOutDivider&:the F-Out divider (0-15)
' Source& :the signal source for F-Out
' Comparel& :status of comparator 1
' Compare2& :status of comparator 2
' TimeOfDay& :time of day mode control
ChipNum% = 2
1 FOutDivider& = 1 ' sets up OSC OUT for 10kHz signal which can
' Source& = FREQ1 ' be used as frequency source for this example
' Comparel& = Disabled
' Compare2& = Disabled
1 TimeOfDay& = Disabled
ULStat& = cbC9513lnit(Boardl\lum%, ChipNum%, FOutDivider&, Source&, Comparel&, Compare2&,
TimeOfDay&)
If ULStat&<>0 Then Stop
End Sub
Private Sub DataSet(BoardNum%)
lnitlUL(BoardNum%)
' First channel, the reference temperature, is counted longer to improve resolution see Gatelntervall& versus
Gatelnterval&
' Resolution is directly dependent on Gatelnterval
ULStat& = cbCFreqln(Boardl\lum%, SigSource(l), Gatelntervall&, CBCount%(l), Freq&(l))
IfULStat&oOThen
ULStat& = cbCFreqln(BoardNum%, SigSource(l), Gatelntervall&, CBCount%(l), Freq&(l))
If ULStat& <> 0 Then MsgBox ("ULStat& = " & ULStat&)
End If
lfCBCount%(l)
-------
Line Input #6, textline
1 = 1 + 1 'I is the number of calibration curves
Loop
Close #6
'MsgBox(l)
'Allocate space to hold arrays place calibration data into arrays
ReDim Resistor(l), clO(l), c20(l), c30(l), cvO(l), cll(l), c21(l), _
c31(l), cvl(l), c!3(l), c23(l), c33(l), cv3(l), c!4(l), c24(l),_
c34(l), cv4(l), SerNo(l), Tslope(l), Tint(l), TSnl(l), chla(l), _
chlb(l), chlc(l), TSn2(l), ch2a(l), ch2b(l), ch2c(l), TSn3(l), _
ch3a(l), ch3b(l), ch3c(l), TSn4(l), ch4a(l), Ch4b(l), Ch4c(l), _
TSn5(l), ch5a(l), Ch5b(l), Ch5c(l), TSn6(l), Ch6a(l), Ch6b(l), _
Ch6c(l), TSn7(l), Ch7a(l), Ch7b(l), Ch7c(l)
' Import calibration data
Open "Calibration_Data(2).csv" For Input As #7
1 = 1
Do While Not EOF(7)
1 MsgBox (I)
Input #7, Resistor(l), clO(l), c20(l), c30(l), cvO(l), cll(l), c21(l),
c31(l), cvl(l), c!3(l), c23(l), c33(l), cv3(l), c!4(l), c24(l), _
c34(l), cv4(l), SerNo(l), Tslope(l), Tint(l), TSnl(l), chla(l), _
chlb(l), chlc(l), TSn2(l), ch2a(l), ch2b(l), ch2c(l), TSn3(l), _
ch3a(l), ch3b(l), ch3c(l), TSn4(l), ch4a(l), Ch4b(l), Ch4c(l), _
TSn5(l), ch5a(l), Ch5b(l), Ch5c(l), TSn6(l), Ch6a(l), Ch6b(l), _
Ch6c(l), TSn7(l), Ch7a(l), Ch7b(l), Ch7c(l)
1 = 1 + 1
Loop
1 MsgBox (I)
Close #7
LblSiteNo. Caption = SiteNo
IblSamRnd. Caption = Pass%
'Determine the number of records in the data file
ChDir(DirA)
'MsgBox (CurDir)
'MsgBox (fname)
Open fname For Input As #1
Nl = 0
N = 0
Do While Not EOF(l)
Line Input #1, textline
N = N + 1 'n is the number of records or values first redord is header not data
Nl = Nl + 1 'Nl is the number of records or values do not reset this counter first recored is header not data
Loop
N = N-1
N1 = N1-1
68
-------
Close #1
lf(Nl-20<0)Then
MsgBox "Not enough data points. Please check Initial Parameters Initial read data."
GoTo 44
End If
'Allocate space to hold arrays
i********
j = 6 'j is the number of channels of data in the raw data file
ReDim timeO(N + 1) As String '+1?
ReDim t(N) As Double
ReDim tb(N)
ReDim ch(j, Nl) As Double
ReDim chb(j, Nl)
Open fname For Input As #5
j=0
k = 0
Do While Not EOF(5)
lfj = OThen
Input #5, timeO(j), bline, heater, xxx, SerNumber, DomeDia, site ' xxx = tail skiped
End If
j=j + l
k = k+l
Input #5, timeO(j), tb(j), chb(l, j), chb(2, j), chb(3, j), chb(4, j), chb(5, j), chb(6, j)
Loop
Close #5
For i = !To6
For ii = ITo k
ch(i, ii) = CDbl(chb(i, ii)) 'convert byte data to double precision
Next ii
Next i
For ii = ITo k
t(ii) = CDbl(tb(ii)) ' convert byte data to double precision
Next ii
' IblSiteD.Caption = site
' MsgBox (" number of data records = " & j)
'13 j=j-l
' Find serial Number
For i = 2 To I - 1
If SerNumber - SerNo(i) = 0 Then GoTo 100
Nexti
If i = I Then SerNo(i - 1) = "default" ' Serial number not detected
i = i-l
'100 IblSerNo. Caption =SerNo(i)
100 ii = i ' ii is used as the record number of the for the calibration data
' MsgBox (ii)
' Convert Data to temperatures
69
-------
Dim temp '
ReDim Res(6, Nl)
For i = 1 To Nl 'Nl is the number of data records
Res(l,i)=ch(l,i)
1 MsgBox (Res(l, i))
Res(2, i) = Res(l, i) - (ch(2, i) / 10)
Res(3, i) = Res(l, i) - (ch(3, i) / 10)
Res(4, i) = Res(l, i) - (ch(4, i) / 10)
Res(5, i) = Res(l, i) - (ch(5, i) / 10)
Res(6, i) = Res(l, i) - (ch(6, i) / 10)
Res(l, i) = ch(l, i) -1820 '1820 is the series resistance on the circuit board added to Channel 1
Next i
' Now convert to temperature and return as the original variable
Fori = lToNl
If Res(l, i)>OThen
temp = chla(ii) + chlb(ii) * Log(Res(l, i) / 10000)
temp = temp + chlc(ii) * ((Log(Res(l, i) / 10000)) A 3)
ch(l, i) = (I/ temp) - 273.15 ' 273.15 convert from Kelvin to Celsus
Else:ch(l, i) = 0
End If
If Res(2, i)>OThen
temp = ch2a(ii) + ch2b(ii) * Log(Res(2, i) / 10000)
temp = temp + ch2c(ii) * ((Log(Res(2, i) / 10000)) A 3)
ch(2, i) = (I/ temp) -273.15
Else:ch(2, i) = 0
End If
If Res(3, i)>OThen
temp = ch3a(ii) + ch3b(ii) * Log(Res(3, i) / 10000)
temp = temp + ch3c(ii) * ((Log(Res(3, i) / 10000)) A 3)
ch(3, i) = (I/ temp) -273.15
Else:ch(3, i) = 0
End If
If Res(4, i)>OThen
temp = ch4a(ii) + Ch4b(ii) * Log(Res(4, i) / 10000)
temp = temp + Ch4c(ii) * ((Log(Res(4, i) / 10000)) A 3)
ch(4, i) = (I/ temp) -273.15
Else:ch(4, i) = 0
End If
If Res(5, i)>OThen
temp = ch5a(ii) + Ch5b(ii) * Log(Res(5, i) / 10000)
temp = temp + Ch5c(ii) * ((Log(Res(5, i) / 10000)) A 3)
ch(5, i) = (I/ temp) -273.15
Else:ch(5, i) = 0
End If
70
-------
If Res(6, i)>OThen
temp = Ch6a(ii) + Ch6b(ii) * Log(Res(6, i) / 10000)
temp = temp + Ch6c(ii) * ((Log(Res(6, i) / 10000)) A 3)
ch(6, i) = (I/ temp) -273.15
Else:ch(6, i) = 0
End If
Nexti
i******
' Convert temperatures from actual temperature in deg C to differential temperature in Deg C
Fori = lToNl
For jj = 2 To 6
ch(jj, i) = ch(jj, i) - ch(l, i)
Next jj
Next i
'MsgBox (ch(l, 1))
'MsgBox (ch(2, 1))
'Open fname2 For Output As #11 'file for temperatures in deg C
1 Fori = lToN
Write #11, Format(t(i), "#0.0000"), Format(ch(l, i), "#0.0000"), _
Format(ch(2, i), "#0.0000"), Format(ch(3, i), "#0.0000"), Format(ch(4, i), "#0.0000"), _
Format(ch(5, i), "#0.0000"), Format(ch(6, i), "#0.0000")
1 Next i
'Close #11
' Determine if sufficient data for possible analysis
lf(t(j-l)-bline*60<0)Then
MsgBox "Not enough data points. Please check Initial Parameters.
GoTo 44
End If
' Analyze base line data and determine slope and intercept of the base line
jj = 0
Dim offset
offset = bline * 60 + ((5 * Gatelnterval& / 1000) + Gatelntervall& / 1000) ' Note:
' Gatelnterval is in milliseconds calculations are in seconds
'MsgBox ("offset = " & offset)
Do
jj = jj + 1 ' JJ is sample number before heater is turned on
Loop Until (t(jj) - offset) >= 0
jj = jj - 2
'MsgBox ("sample number before heater turned on = " & jj)
'calculate sums used to calculate baseline averages, slopes, and intercepts
For j = 1 To 6
sumx(j) = 0
71
-------
sumy(j) = 0
sumxy(j) = 0
sumxx(j) = 0
Next]
' Determind average values
For ct = 1 To jj
For j = 1 To 6
sumx(j) = sumx(j) + t(ct)
sumy(j) = sumy(j) + ch(j, ct)
Next]
Next ct
For j = 1 To 6
avgx(j) = sumx(j)/jj
avgy(j) = sumy(j)/jj
Next j
'MsgBox ("avgx(4) =" & avgx(4))
'MsgBox(avgx(2))
'MsgBox ("avgy(4)= " & avgy(4))
Dim vail As Double
Dim va!2 As Double
For ct = 1 To jj
For j = 1 To 6
vail = (t(ct) - avgx(j)) * (ch(j, ct) - avgy(j))
va!2 = (t(ct) - avgx(j)) A 2
sumxy(j) = sumxy(j) + vail
sumxx(j) = sumxx(j) + va!2
Nextj
Next ct
' Calculate Slopes and Intercepts
For j = 1 To 6
slope(j) = sumxy(j) / sumxx(j)
slnt(j) = avgy(j) - (slope(j) * avgx(j))
Nextj
Correct for baseline drift
Fori = lToN
For j = 1 To 6
ch(j, i) = ch(j, i) - (slope(j) * t(i) + slnt(j))
Nextj
Next i
Open fnamel For Output As #11 'file with change in temperature adjusted for base line slope
For i = 1 To N
Write #11, Format(t(i), "#0.00"), Format(ch(l, i), "#0.00"), _
Format(ch(2, i), "#0.00"), Format(ch(3, i), "#0.00"), Format(ch(4, i), "#0.00"), _
Format(ch(5, i), "#0.00"), Format(ch(6, i), "#0.00")
Next i
Close #11
72
-------
I********
' find time at maximum signal
' the time is calculated as max time - base line time - (pulse time/2)
' This will be one of the metrics evaluated to determine travel time
Dim jmax(6) As Integer 'sample number at max temperature
For i = !To6
tempmax(i) = ch(i, jj + 1) ' only consider samples collected after baseline data
Nexti
For i = !To6
For j = jj + 2 To N
If (ch(i, j) - tempmax(i) > 0) Then
jmax(i)=j
tempmax(i) = ch(i, j)
End If
Nextj
Nexti
tmax(l) = t(jmax(l)) + Gatelntervall& / 2000
For i = 2 To 6
If tempmax(i) < 0.075 Then ' 0.075 considered minimum peak to be used eliminates some of the noise
tmax(i) = bline * 60 + (heater / 2)
Else: tmax(i) = t(jmax(i)) + (Gatelnterval& / 1000 * (i - 1)) + Gatelnterval& / 2000 + Gatelntervall& / 1000
End If
Next i
'MsgBox ("tempmax(4) = " & tempmax(4))
'MsgBox ("tmax(4) = " & tmax(4))
'Begin Moment analysis - Data will be extrapolated based on a lograthmetic
' function based on a sum of squares fit to the data between 80% of tempmax
' and next occurance of <=0.04deg C or the end of file which ever comes first
' first determine sample number correlated to 80%of tempmax to perform moment analysis
' there should be only one peak and it should occure in the first 40% of the data set
' if peak temperature occures after the 40% of the data points do not make moment analysis
Dim jbegin(6) As Integer 'jbegin is the sample number for the first sample to use in generating fit
Dim jend(6) As Integer 'jend is the last sample to use
For i = 2 To 6
If (jmax(i) + 20 - Nl) > 0 Then 'See if there is sufficient data to calculate moments
GoTo 200
End If
j = jmax(i)
Do Until ((tempmax(i) * 0.8) - ch(i, j) > 0) 'determine sample number for first sample
' to use in interpolating and extrapolating
' moments 80% of peak was qualitatively
' estimated as being sufficient to be away
' from the curvature at the peek
j = j + l
If (j - N = 0) Then
73
-------
jbegin(i) = j
jend(i) = j
GoTo 90
End If
' jend(i) = jend(i) + 1
Loop
jbegin(i)=j
lf(jbegin(i) + 20-Nl)>OThen
GoTo 90
End If
If j>= Ml Then GoTo 89
Do Until (0.05 - ch(i, j) > 0) 'determine last sample number considered valid for
'interpolating or extrapolating moments noise was measured
'at +- 0.005 deg C. data has already been adjusted for baseline
'offset. Resolution is 0.005 deg C.
If j>= Ml Then GoTo 89
Loop
89 jend(i) = j-l
90 Next i
' Fit log transformed data to a streight line first determine if there is enough data
' there should be a minimum of minobs observations to calculate the extrapolation function
Dim sumyy(6)
Dim avgxx(6)
Dim avgyy(6)
Dim Mslope(6) As Double
Dim Mint(6) As Double
Dim minobs As Integer 'minimum number of observations between jbegin(i) and jend(i)
'to permit moment analysis 20 is considered minimum
Dim minsignal As Single 'minimum peak signal in volts for analysis
minobs = 20
minsignal = 0.05
Fori = 2To6 'initialize values
sumxx(i) = 0
sumyy(i) = 0
sumxy(i) = 0
avgxx(i) = 0
avgyy(i) = 0
Mslope(i) = 0
Mint(i)=0
Nexti
For i = 2 To 6
If ((jend(i) - jbegin(i) - minobs) > 0) And (tempmax(i) - minsignal > 0) Then
For j = jbegin(i) To (jend(i) - 1)
kk = j+l
avgxx(i) = avgxx(i) + t(kk)
If ch(i, kk) < 0 Then ch(i, j) = 1E-50
74
-------
avgyy(i) = avgyy(i) + Log(ch(i, kk))
Next]
avgxx(i) = avgxx(i) / (jend(i) - jbegin(i))
avgyy(i) = avgyy(i) / (jend(i) - jbegin(i))
End If
Nexti
For i = 2 To 6
If ((jend(i) - jbegin(i) - minobs) > 0) And (tempmax(i) - minsignal > 0) Then
For j = jbegin(i) To (jend(i) - 1)
kk = j+l
If ch(i, kk) < 0 Then ch(i, kk) = 1E-50
sumxx(i) = sumxx(i) + (t(kk) - avgxx(i)) A 2
sumyy(i) = sumyy(i) + (Log(ch(i, kk)) - avgyy(i)) A 2
sumxy(i) = sumxy(i) + (t(kk) - avgxx(i)) * (Log(ch(i, kk)) - avgyy(i))
Next]
Mslope(i) = sumxy(i) / sumxx(i)
Mint(i) = avgyy(i) - (Mslope(i) * avgxx(i))
End If
Next i
'Calculate zero and first moments
Dim Momzero(6)
Dim Moml(6)
Dim cl, c2
DimTl,T2
Dim timecor 'time correction such that zero time is at the end of base line measurements
Dim jstart 'sample number to begin moment calculations this is to remove switching noise
timecor = bline * 60
For i = 2 To 6
If ((jend(i) - jbegin(i) - minobs) < 0) Or (tempmax(i) - minsignal < 0 Or Mslope(i) > 0) Then
Moml(i) = "N.A."
GoTo 300
Else
Momzero(i) = 0
Moml(i) = 0
ch(i,0) = 0
t(0) = 0
Start moment analysis beginning with sample after heat-pulse
Correct time scale to begin half way through the heat-pulse in moment calculations
First step is to determine sample number after end of heat-pulse
Do Until (t(j) - bline * 60 - heater / 2) > 0
j=j + l
Loop
jstart = j
For j = jstart To jbegin(i) 'first part of curve
Momzero(i) = Momzero(i) + ((t(j) - t(j - 1)) * (ch(i, j) + ch(i, j - 1))) / 2
Moml(i) = Moml(i) + (1 / 6) * ((ch(i, j - 1) * (t(j - 1) - timecor) + ch(i, j) *
(t(j)- timecor)) *(t(j)-t(j-l)) + (ch(i,j-l)_
+ ch(i, j)) * ((t(j) - timecor) * (t(j) - timecor) - (t(j - 1) - timecor) * _
75
-------
(t(j -1) - timecor)))
Next]
'end of first part of curve begin analysis for the curve fit portion of the curve
cl = ch(i,j)
Tl = t(j)
j = l
Do Until (cl - 0.001 * tempmax(i) < 0)
T2 = Tl + 1
c2 = Exp(Mslope(i) * (T2) + Mint(i))
Momzero(i) = Momzero(i) + (cl + c2) / 2
Moml(i) = Moml(i) + (1 / 6) * ((cl * (Tl - timecor) + c2 * (T2 - timecor)) * _
(T2 - Tl) + _
(cl + c2) * ((T2 - timecor) * (T2 - timecor) - (Tl - timecor) * (Tl - timecor)))
T1 = T2
cl = c2
Loop
End If
300 Next i
Dim Mon2, Mom3, Mom4, Mom5, Mom6
'MsgBox(Moml(l))
If ((jend(2) - jbegin(2)) > minobs) And (tempmax(2) - minsignal > 0 And lsNumeric(Moml(2)) _
And Momzero(2) > 0 And lsNumeric(Momzero(2))) Then
1 lblMom2.Caption = Format(Moml(2) / Momzero(2), "#0.0")
Mom2 = Moml(2) / Momzero(2)
Else
1 lblMom2.Caption = Format("N.A.")
Mom2=0
End If
If ((jend(3) - jbegin(3)) > minobs) And (Not Not tempmax(3) - minsignal > 0 And lsNumeric(Moml(3))
And Momzero(3) > 0 And lsNumeric(Momzero(3))) Then
1 IblMomS.Caption = Format(Moml(3) / Momzero(3), "#0.0")
Mom3 = Moml(3) / Momzero(3)
Else
1 IblMomS.Caption = Format("N.A.")
Mom3=0
End If
If ((jend(5) - jbegin(5)) > minobs) And (tempmax(S) - minsignal > 0 And lsNumeric(Moml(4)) _
And Momzero(5) > 0 And lsNumeric(Momzero(5))) Then
1 IblMomS.Caption = Format(Moml(5) / Momzero(5), "#0.0")
Mom5 = Moml(5) / Momzero(5)
Else
1 IblMomS.Caption = Format("N.A.")
Mom5=0
End If
If ((jend(6) - jbegin(6)) > minobs) And (tempmax(6) - minsignal > 0 And lsNumeric(Moml(6)) _
And Momzero(6) > 0 And lsNumeric(Momzero(6))) Then
1 IblMomS.Caption = Format(Moml(6) / Momzero(6), "#0.0")
Mom6 = Moml(6) / Momzero(6)
Else
1 IblMomS.Caption = Format("N.A.")
Mom6 = 0
76
-------
End If
' calculate flows
Dim CL 'Confidence limit
CL = 2.326
' Channel 2
If ((jend(2) - jbegin(2) > minobs) And (tempmax(2) - minsignal > 0) And _
(tmax(2) - bline * 60 - 1.5 * heater > 0) And lsNumeric(Moml(2))) Then
Moml(2) = (Moml(2) / Momzero(2))
Q2 = (clO(ii) - (c30(ii) * Moml(2))) / (Moml(2) + c20(ii))
Q2 = (clO(ii) - (c30(ii) * tmax(2) - (bline * 60 + heater / 2))) / _
(tmax(2) - (bline * 60 + heater / 2) + c20(ii))
'MsgBox(Moml(l))
'MsgBox(clO(ii))
Else
Q2 = 0
End If
' Channel 3
If (jend(3) - jbegin(3) > minobs) And (tempmax(S) - minsignal > 0) And (tmax(3) - bline * 60 - 1.5 * heater > 0 .
And lsNumeric(Moml(3)))Then
Moml(3) = (Moml(3) / Momzero(3))
Q3 = (cll(ii) - c31(ii) * Moml(3)) / (Moml(3) + c21(ii))
Q3 = (cll(ii) - c31(ii) * tmax(3) - (bline * 60 + heater / 2)) / (tmax(3) - (bline * 60 + heater / 2) + c21(ii))
Else
Q3 = 0
End If
' Channel 5
If (jend(5) - jbegin(5) > minobs) And (tempmax(S) - minsignal > 0) And _
(tmax(5) - bline * 60 - 2 * heater > 0 And lsNumeric(Moml(5))) Then
Moml(5) = (Moml(5) / Momzero(5))
Q5 = (c!3(ii) - c33(ii) * Moml(5)) / (Moml(5) + c23(ii))
Q5 = (c!3(ii) - c33(ii) * tmax(5) - (bline * 60 + heater / 2)) / (tmax(5) - (bline * 60 + heater / 2) + c23(ii))
'MsgBox(cl3(ii))
'MsgBox(c33(ii))
'MsgBox(c23(ii))
IblQS.Caption = Format(Q3, "#0.000")
Else
Q5 = 0
End If
' Channel 6 Captions
If (jend(6) - jbegin(6) > minobs) And (tempmax(6) - minsignal > 0) And _
(tmax(6) - bline * 60 - 1.5 * heater > 0 And lsNumeric(Moml(6))) Then
Moml(6) = (Moml(6) / Momzero(6))
Q6 = (c!4(ii) - c34(ii) * Moml(6)) / (Moml(6) + c24(ii))
Q6 = (c!4(ii) - c34(ii) * tmax(6) - (bline * 60 + heater / 2)) / (tmax(5) - (bline * 60 + heater / 2) + c24(ii))
'MsgBox(Q4)
Else
Q6 = 0
End If
' Calculate Darcy velocity
77
-------
I*********
If DomeDia > 0.1 Then 'Data exists to calculate flux
If Q2 <> 0 Then
V2 = Q2 * 60 * 24 / (3.1415 * DomeDia * DomeDia / 4)
Else
V2 = 0
End If
If Q3 <> 0 Then
V3 = Q3* 60*24/(3.1415* DomeDia* DomeDia/ 4)
Else
VI = 0
End If
If Q5 <> 0 Then
V5 = Q5* 60*24/(3.1415* DomeDia* DomeDia/ 4)
Else
V5 = 0
End If
If Q6 <> 0 Then
V6 = Q6* 60*24/(3.1415* DomeDia* DomeDia/ 4)
Else
V6 = 0
End If
End If
'MsgBox ("checks")
' set up arrays for ploting
200 Dim Graph() As Single
Dim x As Integer
Dim myXarrayO As Double Time
Dim myYArrayOQ As Double
Dim myYArraylQ As Double
Dim myYArray2() As Double
Dim myYArraySQ As Double
Dim myYArray4() As Double
Dim myYArraySQ As Double
ReDim myXarray(N)
ReDim myYArrayO(N)
ReDim myYArrayl(N)
ReDim myYArray2(N)
ReDim myYArrayS(N)
ReDim myYArray4(N)
ReDim myYArrayS(N)
'Generate some x y data.
myXarray(O) = t(l)
myYArrayO(O) = ch(l, 1)
78
-------
'MsgBox (ch(2,1))
myYArrayl(O) = ch(2,1)
myYArray2(0) = ch(3,1)
myYArrayS(O) = ch(4,1)
myYArray4(0) = ch(5,1)
myYArrayS(O) = ch(6,1)
'MsgBox N
For x = 1 To N -1
myXarray(x) = t(x) - bline * 60 'value for X-axis
myYArrayO(x) = ch(l, x) 'value for Y-axis Ref Temp
myYArrayl(x) = ch(2, x) ' Ch 2 cable end temp
myYArray2(x) = ch(3, x) 'Ch 3 temp
myYArrayS(x) = ch(4, x) ' Heater Temp
myYArray4(x) = ch(5, x)
myYArrayS(x) = ch(6, x) 'Discharge end temp
Nextx
With TChartl
.AddSeriesscPoint
.Series(0).AddArray UBound(myYArrayO), myYArrayO(), myXarray() ' ref temp
.Series(0).XValues.TempValue = True
.Series(l).AddArray UBound(myYArrayl), myYArrayl(), myXarray() ' Cable end Temp
.Series(2).AddArray UBound(myYArray2), myYArray2(), myXarrayQ
.Series(3).AddArray UBound(myYArray3), myYArray3(), myXarray() ' Heater
.Series(4).AddArray UBound(myYArray4), myYArray4(), myXarrayQ
.Series(5).AddArray UBound(myYArray5), myYArray5(), myXarray()
End With
Refresh
DoEvents
'MsgBox ("check 4")
' Screen. MousePointer = 1
'Convert Suplimential Data
Dim timeOl
Dim ch7, ch8, ch9, ch!9
Open fnameS For Input As #11
Input #11, timeOl, ch7, ch8, ch9, chlO, SeriaINo, site
Close #11
' MsgBox (ch7)
' Convert Channel 7 to temperature
lfch7>OThen
temp = Ch7a(ii) + Ch7b(ii) * Log(ch7) / 10000
temp = temp + Ch7c(ii) * ((Log(ch7)) / 10000 A 3)
ch7 = (I/ temp) -273.15 ' 273.15 convert from Kelvin to Celsus
Else: ch7 = 0
79
-------
End If
' write output file with processed data
'ChDirDirl
Open fnameP For Append As #8
Write #8, Format$(timeO(l)), Format(tmax(2) - (bline * 60+ heater/ 2), "#0.00"), Format(tmax(3) - (bline *
60 + heater/ 2), "#0.00"), _
Format(tmax(5) - (bline * 60 + heater / 2), "#0.00"), Format(tmax(6) - (bline * 60 + heater / 2), "#0.00"), _
Format(Moml(2), "#0.00"), Format(Moml(3), "#0.00"), _
Format(Moml(5), "#0.00"), Format(Moml(6), "#0.00"), _
Format(avgy(l), "#0.00"), Format(ch7, "#0.00"), Format(ch8/ 10000, "#0.00"), _
Format(ch9 / 10000, "#0.00"), FormatfchlO / 10000, "#0.00"), _
Format(tempmax(l), "#0.00"), Format(tempmax(2), "#0.00"), Format(tempmax(3), "#0.00"), _
Format(tempmax(4), "#0.00"), Format(tempmax(5), "#0.00"), Format(tempmax(6), "#0.00")
Close #8
Refresh
DoEvents
TChartS Peak Arrival Times
i********
'Determine number of records No header in file
i******
' determine the number of lines in the data table
Open fnameP For Input As #6
1 = 0
Do While Not EOF(6)
Line Input #6, textline
1 = 1 + 1 'I is the number of lines of data
Loop
Close #6
Dim tmax2(), tmax3(), tmax5(), tmax6(), Mom2d(), Mom3d(), MomSdQ, Mom6d(), Reft(), SupTQ, DifP(),
chan9(),_
chanlO(), txt(34)
Dim DaTime() As Date
ReDim tmax2(l), tmax3(l), tmax5(l), tmax6(l), Mom2d(l), Mom3d(l), Mom5d(l), Mom6d(l), Reft(l), SupT(l),
DifP(l), chan9(l), chanlO(l) ' read in data
ReDim DaTime(l)
i = 4
Open fnameP For Input As #6
Do Until i = 1 + 1
lfi = 4Then
Input #6, txt(l), txt(2), txt(3), txt(4), txt(5), txt(6), txt(7), txt(8), txt(9), txt(10), txt(ll), txt(12), _
txt(13), txt(14), txt(15), txt(16), txt(17), txt(18), txt(19), txt(20), txt(21), _
txt(22), txt(23), txt(24), txt(25), txt(26), txt(27), txt(28), txt(29), _
txt(30), txt(31), txt(32), txt(33), txt(34), _
timeO(i - 4), tmax2(i - 4), tmax3(i - 4), tmax5(i - 4), tmax6(i - 4), Mom2d(i - 4), Mom3d(i - 4), _
Mom5d(i - 4), Mom6d(i - 4), Reft(i - 4), SupT(i - 4), DifP(i - 4), chan9(i - 4), chanlO(i - 4), _
80
-------
tempmax(l), tempmax(2), tempmax(3), tempmax(4), tempmax(5), tempmax(6)
Elself i > 4 Then
Input #6, timeO(i - 4), tmax2(i - 4), tmax3(i - 4), tmax5(i - 4), tmax6(i - 4), Mom2d(i - 4), Mom3d(i - 4),
Mom5d(i - 4), Mom6d(i - 4), Reft(i - 4), SupT(i - 4), DifP(i - 4), chan9(i - 4), chanlO(i - 4), _
tempmax(l), tempmax(2), tempmax(3), tempmax(4), tempmax(5), tempmax(6)
End If
DaTime(i - 4) = timeO(i-4)
Loop
Close #6
MsgBox (I)
MsgBox (DaTime(l))
MsgBox (tmax6(l))
' set up arrays for ploting
Dim Graph3() As Single
Dim myXarray2() As Double
Dim myYArray20() As Double
Dim myYArray21() As Double
Dim myYArray22() As Double
Dim myYArray23() As Double
Dim myYArray24() As Double
Dim myYArray25() As Double
Dim myYArray26() As Double
Dim myYArray27() As Double
Dim myYArray28() As Double
Dim myYArray29() As Double
Dim myYArray210() As Double
Dim myYArray211() As Double
Dim myYArray212() As Double
ReDim
ReDim
ReDim
ReDim
ReDim
ReDim
ReDim
ReDim
ReDim
ReDim
ReDim
ReDim
ReDim
ReDim
myXarray2(l)
myYArray20(l)
myYArray21(l)
myYArray22(l)
myYArray23(l)
myYArray24(l)
myYArray25(l)
myYArray26(l)
myYArray27(l)
myYArray28(l)
myYArray29(l)
myYArray210(l)
myYArray211(l)
myYArray212(l)
1 tmax2 time to peak temp
1 tmaxS
1 tmaxS
1 tmaxS
1 Mom2d first temporal moment
1 MomSd
1 MomSd
1 MomSd
1 Reft Reference temperature
1 SupT Extra temperature measurement Outside/Inside Dome
' DifP Differential pressure
1chan9
1chan10
'Generate some x y data.
81
-------
myXarray2(0) = DaTime(O)
myYArray20(0) = tmax2(0)
'MsgBox (ch(2,1))
myYArray21(0) = tmax3(0)
myYArray22(0) = tmax5(0)
myYArray23(0) = tmax6(0)
If lsNumeric(Mom2d(0)) Then
myYArray24(0) = Mom2d(0)
Else: myYArray24(0) = Empty
End If
If lsNumeric(Mom3d(0)) Then
myYArray25(0) = Mom3d(0)
Else: myYArray25(0) = Empty
End If
If lsNumeric(Mom5d(0)) Then
myYArray26(0) = Mom5d(0)
Else: myYArray26(0) = Empty
End If
If lsNumeric(Mom6d(0)) Then
myYArray27(0) = Mom6d(0)
Else: myYArray27(0) = Empty
End If
myYArray28(0) = Reft(O)
myYArray29(0) = SupT(O)
myYArray210(0) = DifP(O)
myYArray211(0) = chan9(0)
myYArray212(0) = chanlO(O)
'MsgBox N
For i = OTo I
myXarray2(i) = DaTime(i) 'value for X-axis
myYArray20(i) = tmax2(i) 'time to peek channel 2
myYArray21(i) = tmax3(i) '
myYArray22(i) = tmax5(i) '
myYArray23(i) = tmax6(i) '
If lsNumeric(Mom2d(i)) Then
myYArray24(i) = Mom2d(i)
Else: myYArray24(i) = Empty
End If
If lsNumeric(Mom3d(i)) Then
myYArray25(i) = Mom3d(i)
Else: myYArray25(i) = Empty
End If
If lsNumeric(Mom5d(i)) Then
myYArray26(i) = Mom5d(i)
Else: myYArray26(i) = Empty
End If
If lsNumeric(Mom6d(i)) Then
myYArray27(i) = Mom6d(i)
Else: myYArray27(i) = Empty
End If
82
-------
myYArray28(i) = Reft(i)
myYArray29(i) = SupT(i)
myYArray210(i) = DifP(i)
myYArray211(i) = chan9(i)
myYArray212(i) = chanlO(i)
Next i
With TChartS
.AddSeries scPoint
.Series(0).AddArray UBound(myYArray20), myYArray20(), myXarray2()'
.Series(0).XValues.DateTime = True
.Series(l).AddArray UBound(myYArray21), myYArray21(), myXarray2() '
.Series(l).XValues.DateTime = True
.Series(2).AddArray UBound(myYArray22), myYArray22(), myXarray2()
.Series(2).XValues.DateTime = True
.Series(3).AddArray UBound(myYArray23), myYArray23(), myXarray2()
.Series(3).XValues.DateTime = True
End With
Refresh
DoEvents
With TChart4
.AddSeries scPoint
.Series(0).AddArray UBound(myYArray24), myYArray24(), myXarray2()'
.Series(0).XValues.DateTime = True
.Series(l).AddArray UBound(myYArray25), myYArray25(), myXarray2() '
.Series(l).XValues.DateTime = True
.Series(2).AddArray UBound(myYArray26), myYArray26(), myXarray2()
.Series(2).XValues.DateTime = True
.Series(3).AddArray UBound(myYArray27), myYArray27(), myXarray2()
.Series(3).XValues.DateTime = True
End With
Refresh
DoEvents
With TChartS
.AddSeries scPoint
.Series(0).AddArray UBound(myYArray28), myYArray28(), myXarray2()'
.Series(0).XValues.DateTime = True
.Series(l).AddArray UBound(myYArray29), myYArray29(), myXarray2() '
.Series(l).XValues.DateTime = True
.Series(2).AddArray UBound(myYArray210), myYArray210(), myXarray2()
.Series(2).XValues.DateTime = True
End With
Refresh
DoEvents
With TChartS
.AddSeries scPoint
.Series(0).AddArray UBound(myYArray211), myYArray211(), myXarray2()'
.Series(0).XValues.DateTime = True
.Series(l).AddArray UBound(myYArray212), myYArray212(), myXarray2()
.Series(l).XValues.DateTime = True
End With
83
-------
Refresh
DoEvents
GoTo 45
44 Unload Me
' Screen.MousePointer = 1
45 End Sub
84
-------
frmDataAnalysis
' last revised 04/20/05 by Carl Enfield
Private Sub Dir2_Change()
File2.Path = Dir2
End Sub
Private Sub File2_Click()
Dim textline
Dim dt, arm, yl, y2 As Double
Dim timeOQ As String
Dim timel As String
Dim chl(), ch2(), ch3(), ch4(), ch5(), ch6(), ch7(), ch8()
Dimt()
Dim al, a2, a3, a4, a5, a6 As Double
Dim ml, m2, m3, m4, m5, m6 As Double
Dim cl() As Double
Dim c2() As Double
Dim c3() As Double
Dim c4() As Double
Dim c5() As Double
Dim c6() As Double
Dimxl() As Double
Dim x2() As Double
Dimx3() As Double
Dim x4() As Double
Dimx5() As Double
Dim x6() As Double
Dim ct, ipt, pt_odd, pt_half As Integer
Dim bl, si As Single
Dim b2, s2 As Single
Dim b3, s3 As Single
Dim b4, s4 As Single
Dim b5, s5 As Single
Dim b6, s6 As Single
Dim dl, d2, d3, d4, d5, d6 As Double
Dim tmaxl, tmax2, tmaxS, tmax4, tmaxS, tmaxS As Double
Dim xlmax, x2max, xSmax, x4max, xSmax, xSmax
Dim pt
Dim baseline$, heater
pt = CDbl(txtPtavg.Text) 'points average
If pt = 0 Then MsgBox "Invalid Points Average Number.": GoTo 44
pt_odd = pt Mod 2
If pt_odd <> 1 Then MsgBox "Please enter an odd number for points average.": GoTo 44
pt_half=pt/2-0.5
baseline$ = 10#
heater = 5#
85
-------
heater = heater + baseline$
'ChDir "c:\temp\"
ChDirDir2
Open File2 For Input As #1
N = 0
Do While Not EOF(l)
Line Input #1, textline
N = N + 1
Loop
Close #1
k = N - pt + 1
g = k+ pt_half
ReDim xl(g), x2(g), x3(g), x4(g), x5(g), x6(g)
ReDim timeO(N + 1) As String
ReDim chl(N), ch2(N), ch3(N), ch4(N), ch5(N), ch6(N), t(N), ch7(N), ch8(N)
ReDim cl(N), c2(N), c3(N), c4(N), c5(N), c6(N)
Open File2 For Input As #3
Input #3, timel
Close #3
If (Year(timel)) > 2004 Then GoTo 3
If (Year(timel) >= 2004 And Month(timel) > 4) Then GoTo 3
If Month(timel) = 4 And Day(timel) >= 26 Then GoTo 3 Else GoTo 4
3 Open File2 For Input As #2
j = 0
N = 0
Do While Not EOF(2)
lfj = OThen
Input #2, timeO(j), baseline, heater
End If
j = j+l
N = N + 1
Input #2, timeO(j), t(j), chl(j), ch2(j), ch3(j), ch4(j), ch5(j), ch6(j), ch7(N), ch8(N)
'baseline = 10 sec.
'heater = 3 sec.
'If (t(j) > 10 And t(j) <= 13) Then GoTo 13 'skip data collected while heater was on
Loop
GoTo 5
4 Open File2 For Input As #2
j = 0
Do While Not EOF(2)
j = j+l
14 Input #2, timeO(j), t(j), chl(j), ch2(j), ch3(j), ch4(j), ch5(j), ch6(j), ch7(j), ch8(j)
'If t(j) > 10 And t(j) < 14 Then GoTo 14 'skip data collected while heater was on
86
-------
Loop
5 Close #2
If t(j) <= baseline$Then MsgBox "Not enough data points. Please check file.": GoTo 44
ct = l
Do While t(ct) <= baseline$
sl = sl + chl(ct)
s2 = s2 + ch2(ct)
s3 = s3 + ch3(ct)
s4 = s4 + ch4(ct)
s5 = s5 + ch5(ct)
's6 = s6 + ch6(ct)
ct = ct + 1
Loop
bl = sl/(ct-l)
b2 = s2 / (ct -1)
b3 = s3 / (ct -1)
b4 = s4 / (ct -1)
b5 = s5 / (ct -1)
'b6 = s6/(ct-l)
For i = 1 To j
dl = chl(i)-bl
cl(i) = dl
d2 = ch2(i)-b2
c2(i) = d2
d3 = ch3(i)-b3
c3(i) = d3
d4 = ch4(i) - b4
c4(i) = d4
d5 = ch5(i)-b5
c5(i) = d5
1 d6 = ch6(i)-b6
1 c6(i) = d6
Next i
'take pt points average
For m = pt_half + 1 To k + ptjialf
xl(m) = 0
x2(m) = 0
x3(m) = 0
x4(m) = 0
x5(m) = 0
1 x6(m)=0
For z = m - pt_half To m + pt_half
xl(m) = xl(m) + cl(z)
x2(m) = x2(m) + c2(z)
x3(m) = x3(m) + c3(z)
x4(m) = x4(m) + c4(z)
x5(m) = x5(m) + c5(z)
x6(m) = x6(m) + c6(Z)
87
-------
Nextz
xl(m) = xl(m)/pt
x2(m) = x2(m)/pt
x3(m) = x3(m)/pt
x4(m) = x4(m)/ pt
x5(m) = x5(m)/pt
x6(m) = x6(m)/pt
Next m
'find t at maximum signal
xlmax = xl(l)
x2max = x2(l)
xSmax = x3(l)
x4max = x4(l)
xSmax = x5(l)
x6max = x6(l)
tmaxl = t(l)
tmax2 = t(l)
tmax3 = t(l)
tmax4 = t(l)
tmax5 = t(l)
tmax6 = t(l)
For i = 1 + pt_half To k + pt_half
If (xl(i) > xlmax) Then
xlmax = xl(i)
tmaxl = t(i)
End If
If (x2(i) > x2max) Then
x2max = x2(i)
tmax2 = t(i)
End If
If (x3(i) >x3max)Then
xSmax = x3(i)
tmaxS = t(i)
End If
If (x4(i) >x4max)Then
x4max = x4(i)
tmax4 = t(i)
End If
If (x5(i) >x5max)Then
xSmax = x5(i)
tmax5 = t(i)
End If
If (x6(i) >x6max)Then
x6max = x6(i)
tmaxS = t(i)
End If
Nexti
al = 0
ml = 0
88
-------
For i = 1 To k -1
yl = xl(i)
y2 = xl(i + 1)
dt = t(i + l)-t(i)
al = al+((yl + y2)/2)*dt
arm = (t(i) * 1 + t(i + 1) * 1) / 2 'moment arm
ml = ml + (((yl + y2) / 2) * dt) * arm
Nexti
a2=0
m2 = 0
For i = 1 To k -1
yl = x2(i)
y2 = x2(i + 1)
dt = t(i + l)-t(i)
a2 = a2 + ((yl + y2)/2)*dt
arm = (t(i)*l + t(i + l)* I)/2
m2 = m2 + (((yl + y2) / 2) * dt) * arm
Next i
a3=0
m3 = 0
For i = 1 To k -1
yl = x3(i)
y2 = x3(i + 1)
dt = t(i + l)-t(i)
a3 = a3 + ((yl + y2)/2)*dt
arm = (t(i)*l + t(i + l)* I)/2
m3 = m3 + (((yl + y2) / 2) * dt) * arm
Next i
a4 = 0
m4 = 0
For i = 1 To k -1
yl = x4(i)
y2 = x4(i + 1)
dt = t(i + l)-t(i)
a4 = a4+((yl + y2)/2)*dt
arm = (t(i)*l + t(i + l)* I)/2
m4 = m4 + (((yl + y2) / 2) * dt) * arm
Next i
a5=0
m5 = 0
For i = 1 To k -1
yl = x5(i)
y2 = x5(i + 1)
dt = t(i + l)-t(i)
a5 = a5 + ((yl + y2)/2)*dt
arm = (t(i)*l + t(i + l)* I)/2
m5 = m5 + (((yl + y2) / 2) * dt) * arm
Next i
89
-------
1 a6=0
1 m6 = 0
1 Fori = lTok-l
1 yl = x6(i)
1 y2 = x6(i + 1)
1 dt = t(i + 1) - t(i)
1 a6 = a6 + ((yl + y2)/2)*dt
1 arm = (t(i)*l + t(i+ !)*!)/2
1 m6 = m6 + (((yl + y2)/2)*dt)*arm
1 Nexti
Textl.Text = Format(al, "#.00") 'zero moment
Text2.Text = Format(a2, "#.00")
TextS.Text = Format(a3, "#.00")
Text4.Text = Format(a4, "#.00")
TextS.Text = Format(a5, "#.00")
1 TextS.Text = Format(a6, "#.00")
Text7.Text = Format(ml / al, "#.00") 'normalized first moment
TextS.Text = Format(m2 / a2, "#.00")
Text9.Text = Format(m3 / a3, "#.00")
TextlO.Text = Format(m4/ a4, "#.00")
Textll.Text = Format(m5 / a5, "#.00")
1 Textl2.Text = Format(m6 / a6, "#.00")
TextlS.Text = Format(tmaxl, "#.00") 'tmax
Textl4.Text = Format(tmax2, "#.00")
TextlS.Text = Format(tmax3, "#.00")
TextlS.Text = Format(tmax4, "#.00")
Textl7.Text = Format(tmax5, "#.00")
1 TextlS.Text = Format(tmax6, "#.00")
44 End Sub
Private Sub CmdClose_Click()
Unload Me
End Sub
Private Sub Form_Load()
Dir2.Path = "c:\temp"
File2.Path = Dir2
End Sub
90
-------
frmDifferentialPres
'Last Revision 11/26/07 by Carl Enfield
' The purpose is to display differential pressure between the surrounding water and
' inside the dome.
'Const BoardNum% = 0 ' Board number
Const FirstPoint& = 0 'set first element in buffer to transfer to array
Const Gatelntervall& = 500 ' Number of milliseconds used to count pulses and determine frequency used on ref
temperature only
Dim CBCount%(10)
Dim ChipNum%
Dim counter%, rounds, RateCount
Dim DataValue% ' digital value to write to Auxport
Dim DomeDia ' diameter of dome in centimeters
Dim fnamePres As String * 50 ' file neame for raw data with ext .dat
Dim FOutDivider&
Dim Freq&(10)
Dim i, k, j, N As Integer ' counter
Dim MemHandle& ' define a variable to contain the handle for
' memory allocated by Windows through cbWinBufAlloc%()
Public CycleTime
Dim SampleTime, startpause
Dim SigSource(S)
Dim ULStat As Long 'universal library styatistic return
Dim site As String * 64 ' a descriptive string describing the site
Dim CallStop, mytime
Dim ActualCounts%(8)
Dim Logil, Logi2 As Boolean
Private Sub cmdStart_Click()
i********
' Determine which device Site 1 or Site 2
lfChecklAndCheck2Then
MsgBox ("Check only one box")
Call CmdStop_Click
End If
Logil = Checkl
Logi2 = Check2
'MsgBox (Logil)
'MsgBox (Logi2)
'MsgBox (Logil + Logi2)
If Logil + Logi2 = 0 Then
MsgBox ("Check one site")
Call CmdStop_Click
End If
If Checkl Then BoardNum% = 0
If Check2 Then BoardNum% = 1
91
-------
I*******
' Turn off heater turn on power supply
ULStat = cbDBitOut(Boardl\lum%, AUXPORT, 0, 1) ' turn on power supply
If ULStatoO Then Stop
ULStat = cbDBitOut(Boardl\lum%, AUXPORT, 1, 0) ' turn off heater
If ULStatoO Then Stop
1 MsgBox (BoardNum%)
' Initialize signal sources
' There are five inputs on each chip so use chip 2 for the second 5 signal sources
SigSource(S) = CTRINPUT8 ' Input 8 chip 2 Pressure transducer
' Collect Screen Input Data
MemHandle& = cbWinBufAlloc(lOO)
If MemHandle& = 0 Then Stop
site = CStr(txtSite.Text) ' site description
'write baseline, heat-pulse duration, SeriaINo, and Site description raw data file #1
i******
Date = Now
Time = Now
fnamePres = "c:\temp\" + Format$(Date, "mmm dd ") + Format$(Time, "hh mm ss") + "Pres" + ".dat"
Open fnamePres For Output As #1
Write #1, Format(Now), Format(site)
Close #1
MemHandle& = cbWinBufAlloc(lOO)
If MemHandle& = 0 Then Stop
ULStat = cbDBitOut(BoardNum%, AUXPORT, 0, 1) ' turn on power supply
If ULStatoO Then Stop
Screen. MousePointer = 0 'Default Pointer
Refresh
Date = Now
Time = Now
timezero = Timer
t = 0
k = 0
Do While CallStop o-l 'main loop
'calculate accmulated time
If Timer < timezero Then 'if cross midnight
t = 86400 - timezero + Timer
Else
t = Timer -timezero
End If
i*******
92
-------
'control when to activate digital output for rempte power supply and heater
DataSet ' Collect a data point
Open fnamePres For Append As #1
'write data in frequency
Write #1, Format(Now), Format(t, "#0.00000"), _
Format(Freq&(8), "#0")
Close #1
'display data in frequency
IblADData. Caption = Format(Freq&(8), "#0")
PlotData (fnamePres)
Refresh
DoEvents
If CallStop = -1 Then GoTo 44
Loop 'Loop 11
44 End Sub
Private Sub CmdStop_Click()
CallStop = -1
ULStat = cbDBitOut(0, AUXPORT, 0, 0) ' turn off power supply
If ULS tat <>0 Then Stop
ULStat = cbDBitOut(0, AUXPORT, 1, 0) ' turn off heater
If ULS tat <>0 Then Stop
ULStat = cbDBitOut(l, AUXPORT, 0, 0) ' turn off power supply
If ULS tat <>0 Then Stop
ULStat = cbDBitOut(l, AUXPORT, 1, 0) ' turn off heater
If ULS tat <>0 Then Stop
Unload Me
End Sub
Private Sub DataSetQ
InitlUL
lnit2UL
i = 8
ULStat& = cbCFreqln(Boardl\lum%, SigSource(i), Gatelntervall&, CBCount%(i), Freq&(i))
If ULStat&<>0 Then Stop
lfCBCount%(i)
-------
ULStat& = cbDeclareRevision(CURRENTREVNUM)
' Initiate error handling
' activating error handling will trap errors like
' bad channel numbers and non-configured conditions.
' Parameters:
' PRINTALL :all warnings and errors encountered will be printed
' DONTSTOP :if an error is encountered, the program will not stop,
' errors must be handled locally
ULStat& = cbErrHandling(PRINTALL, DONTSTOP)
If ULStat&<>0 Then Stop
1 If cbErrHandling& is set for STOPALL or STOPFATAL during the program
' design stage, Visual Basic will be unloaded when an error is encountered.
' We suggest trapping errors locally until the program is ready for compiling
' to avoid losing unsaved data during program design. This can be done by
' setting cbErrHandling options as above and checking the value of ULStat%
' after a call to the library. If it is not equal to 0, an error has occurred.
' Initialize the board level features
' Parameters:
' BoardNum& :the number used by CBCONFIG to describe this board
1 ChipNum& :chip to be initialized (1 for CTR5, lor 2 for CTR10)
1 FOutDivider&:the F-Out divider (0-15)
' Source& :the signal source for F-Out
' Comparel& :status of comparator 1
' Compare2& :status of comparator 2
' TimeOfDay& :time of day mode control
ChipNum% = l
FOutDivider& = 1 ' sets up OSC OUT for 10kHz signal which can
Source& = FREQ1 ' be used as frequency source for this example
Comparel& = Disabled
Compare2& = Disabled
TimeOfDay& = Disabled
ULStat& = cbC9513lnit(Boardl\lum%, ChipNum%, FOutDivider&, Source&, Comparel&, Compare2&,
TimeOfDay&)
If ULStat&<>0 Then Stop
End Sub
Private Sub lnit2UL()
' Initialize the board level features
' Parameters:
' BoardNum& :the number used by CBCONFIG to describe this board
1 ChipNum& :chip to be initialized (1 for CTR5, lor 2 for CTR10)
1 FOutDivider&:the F-Out divider (0-15)
' Source& :the signal source for F-Out
' Comparel& :status of comparator 1
' Compare2& :status of comparator 2
' TimeOfDay& :time of day mode control
ChipNum% = 2
94
-------
1 FOutDivider&= 1 ' sets up OSC OUT for 10kHz signal which can
' Source& = FREQ1 ' be used as frequency source for this example
' Comparel& = Disabled
' Compare2& = Disabled
1 TimeOfDay& = Disabled
Source& = FREQ1 ' be used as frequency source for this example
ULStat& = cbC9513lnit(Boardl\lum%, ChipNum%, FOutDivider&, Source&, Comparel&, Compare2&,
TimeOfDay&)
If ULStat&<>0 Then Stop
End Sub
Private Sub PlotData(fnamePres)
' Screen. MousePointer = 11
Dim textline
Dim timeO() As String
Dim timel As String 'used to store time file created needed to determine age of file and be compatable with
older data
Dim chb() 'dynamic array to hold analog byte data stored on .dat file
Dim ch() As Double 'dynamic array to hold double precision data from .dat file
Dim tb() ' dynamic array to hold elapsed time data in byte format
Dim t() As Double 'dynamic array to hold elapsed time data
Dim ct As Integer ' , pt_odd, ptjialf As Integer
Dim tmax(6) As Double
Dim tempmax(6) As Double
Dim SerNumber
' Dim tail
Dim N, j, i, I, m, g, k, z, ii, jj, kk As Integer
Dim site As String * 64
Dim dat_type As Boolean
'Determine the number of records in the data file
i*********
ChDir ("C:\Temp")
Open fnamePres For Input As #1
Nl = 0
N = 0
Do While Not EOF(l)
Line Input #1, textline
N = N + 1 'n is the number of records or values first redord is header not data
Nl = Nl + 1 'Ml is the number of records or values do not reset this counter first recored is header not data
Loop
N = N-1
N1 = N1-1
Close #1
i********
'Allocate space to hold arrays
ReDim timeO(N + 1) As String '+1?
ReDim t(N) As Double
95
-------
ReDimtb(N)
ReDim ch(Nl) As Double
ReDim chb(Nl)
Open fnamePres For Input As #5
j=0
k = 0
Do While Not EOF(5)
lfj = OThen
Input #5, timeO(j), site ' xxx = tail skiped
End If
j=j+l
k = k + l
Input #5, timeO(j), tb(j), chb(j)
Loop
Close #5
For ii = ITo k
ch(ii) = CDbl(chb(ii)) 'convert byte data to double precision
Next ii
For ii = ITo k
t(ii) = CDbl(tb(ii)) ' convert byte data to double precision
Next ii
' set up arrays for ploting
200 Dim Graph() As Single
Dim x As Integer
Dim myXarrayO As Double Time
Dim myYArrayOQ As Double
ReDim myXarray(N)
ReDim myYArrayO(N)
'Generate some x y data.
myXarray(O) = t(l)
myYArrayO(O) = ch(l)
'MsgBox N
For x = 1 To N -1
myXarray(x) =t(x)
myYArrayO(x) = ch(x)
Next x
'value for X-axis
'value for Y-axis Ref Temp
With TChartl
.AddSeriesscPoint
96
-------
.Series(0).AddArray UBound(myYArrayO), myYArrayO(), myXarray()' ref temp
.Series(0).XValues.TempValue = True
End With
Refresh
DoEvents
'MsgBox ("check 4")
' Screen.MousePointer = 1
45 End Sub
97
-------
frmMain
Option Explicit
Private Sub MDIForm_Load()
frmMain.Show
frmMain.AutoShowChildren =True
frmDataAcquisition.Show
'frmDataDisplay.Show
'frmDifferentialPres.Show
frmDataAcquisition.Initialize = True
End Sub
Private Sub mnuFOpen_Click()
dlgCommon.FILTER = "All Files (*.*) | *.* | Data Files(*.dat) | *.dat"
dlgCommon.ShowOpen
On Error GoTo errorhandler
If FileLen(dlgCommon.Filename) > 65000 Then
MsgBox "The file is too large to open."
Exit Sub
End If
'frmDataDisplay.Show
frmDifferentialPres.Show
OpenFile (dlgCommon.Filename)
errorhandler: Exit Sub
End Sub
Private Sub mnuHAbout_Click()
MsgBox "Flux Meter Data Acquisition Program" _
&vbCrLf_
& vbCrLf & "Beta Version 5.01" _
&vbCrLf_
&vbCrLf&"CarlGEnfield"_
& vbCrLf & "Enfield Environmental Research" _
& vbCrLf & "5017 Park Ridge Ct" _
& vbCrLf & "West Chester, OH 45069, USA" _
& vbCrLf _
& vbCrLf & "Bob K Lien" _
& vbCrLf & "U.S. EPA"_
& vbCrLf &"NRMRL"_
& vbCrLf & "Cincinnati, Ohio" _
, vblnformation
End Sub
Private Sub mnuMAInScan_Click()
frmDataAcquisition.Show
frmDataAcquisition.SetFocus
End Sub
98
-------
Private Sub mnuMReference_Click()
MsgBox "Measurement Computing. Universal Driver Library: User's Guide"
& vbCrLf
End Sub
Private Sub mnuTAnalysis_Click()
frmDataAnalysis.Show
End Sub
Private Sub mnuTDisplay_Click()
frmSignalDisplay.Show
End Sub
Private Sub mnuTPres_Click()
frmDifferentialPres.Show
End Sub
Private Sub mnuWArrangelcons_Click()
Me.Arrange vbArrangelcons
End Sub
Private Sub mnuWCascade_Click()
Me.Arrange vbCascade
End Sub
Private Sub mnuWTileHorizontal_Click()
Me.Arrange vbTileHorizontal
End Sub
Private Sub mnuWTileVertical_Click()
Me.Arrange vbTileVertical
End Sub
Private Sub mnuFCIose_Click()
On Error GoTo errorhandler
ActiveForm.CmdClose.Value = True
errorhandler: Exit Sub
End Sub
Private Sub mnuFExit_Click()
Unload Me
End Sub
Private Sub mnuFPrint_Click()
Dim NumCopies, i
dlgCommon.CancelError = True
On Error GoTo errorhandler
If ActiveForm Is Nothing Then Exit Sub
99
-------
dlgCommon.Flags = 0
dlgCommon.ShowPrinter
NumCopies = dlgCommon.Copies
For i = ITo NumCopies
Screen. ActiveForm.PrintForm
Next i
errorhandler: Exit Sub
End Sub
Private Sub mnuFPrintSetup_Click()
dlgCommon.CancelError = False
dlgCommon.Flags = cdlPDPrintSetup
dlgCommon.ShowPrinter
End Sub
Sub OpenFile(Filename)
Dim flndex As Integer
On Error Resume Next
' Open the selected file.
Open Filename For Input As #1
If Err Then
MsgBox "Can't open file:" + Filename
Exit Sub
End If
' Change the mouse pointer to an hourglass.
Screen. MousePointer = 11
' Change the form's caption and display the new text.
'frmDataDisplay.Caption = UCase(FileName)
'frmDataDisplay.txtNote.Text = StrConv(lnputB(LOF(l), 1), vbUnicode)
Close #1
' Reset the mouse pointer.
Screen.MousePointer = 0
End Sub
100
-------
frmSignal Display
'last revised 05/19/08 by Carl Enfield
Option Explicit
Private Sub CmdClose_Click()
Unload Me
End Sub
Private Sub Dirl_Change()
Filel.Path = Dirl
End Sub
Private Sub Form_Load()
Dirl. Path = "c:\temp"
Filel. Path = Dirl
End Sub
Private Sub Commandl_Click()' Display Information
Screen.MousePointer = 11
Dim textline
Dim timeOQ As String
Dim timea As String
Dim ch7, ch8, ch9, chlO As Double
Dim timel As String 'used to store time file created needed to determine age of file and be compatable with
older data
Dim chb() 'dynamic array to hold analog byte data stored on .dat file
Dim ch() As Double 'dynamic array to hold analog data stored on .dat file
Dim t() As Double 'dynamic array to hold elapsed time data
Dim tb() 'dynamic array to hold byte elapsed time data
Dim ct As Integer ' , pt_odd, pt_half As Integer
Dim tmax(6) As Double
Dim tempmax(6)
Dim bline
Dim heater
Dim SerNumber
Dim tail
Dim N, Nl, j, i, I, m, g, k, z, ii, jj, kk As Integer
Const Gatelnterval& = 200
Const Gatelntervall& = 1000
Dim CalResistor
Dim Resistor))
Dim DomeDia
Dim site As String * 64
Dim QO, Ql, Q3, Q4, QOmin, Qlmin, QSmin, Q4min, QOmax, Qlmax, QSmax, Q4max
Dim VO, VI, V3, V4
Dim dat_type As Boolean
Dim sumx(6) As Double
Dim sumy(6) As Double
Dim sumxx(6) As Double
Dim sumxy(6) As Double
Dim avgx(6) As Double
101
-------
Dim avgy(6) As Double
Dim slope(6) As Double
Dim slnt(6) As Double
Dim File2 As String
ChDirDirl
' Initialize all of the calculated outputs other than peak arrival time to "N.A.'
' to make sure only recomputed values are displayed if this isn't done
' it is possibe to display a value from a previous data set
lblMom2.Caption = Format("N.A.")
IblMomS.Caption = Format("N.A.")
IblMomS.Caption = Format("N.A.")
IblMomS. Caption = Format("N.A.")
QO = "N.A."
QOmin = "N.A."
QOmax = "N.A."
IblQO.Caption = Format(QO)
IblQOmin. Caption = Format(QOmin)
IblQOmax.Caption = Format(QOmax)
Q1 = "N.A."
Qlmin = "N.A."
Qlmax = "N.A."
IblQl.Caption = Format(Ql)
IblQlmin.Caption = Format(Qlmin)
IblQlmax. Caption = Format(Qlmax)
Q3 = "N.A."
Q3min = "N.A."
Q3max = "N.A."
IblQS.Caption = Format(QS)
IblQSmin.Caption = Format(Q3min)
IblQSmax.Caption = Format(Q3max)
Q4 = "N.A."
Q4min = "N.A."
Q4max = "N.A."
lblQ4.Caption = Format(Q4)
lblQ4min.Caption = Format(Q4min)
lblQ4max.Caption = Format(Q4max)
VO = "N.A."
VI = "N.A."
V3 = "N.A."
V4 = "N.A."
IblVO.Caption = Format(VO)
IblVl.Caption = Format(Vl)
IblVS.Caption = Format(VS)
lblV4.Caption = Format(V4)
' Import calibration data
i******
ChDir ("C:\Program Files\FluxMeter\")
102
-------
' determine the number of calibration curves in the data table
Open "Calibration_Data(2).csv" For Input As #6
1 = 0
Do While Not EOF(6)
Line Input #6, textline
1 = 1 + 1 'I is the number of calibration curves
Loop
Close #6
'MsgBox (I)
'Allocate space to hold arrays place calibration data into arrays
ReDim Resistor(l), clO(l), c20(l), c30(l), cvO(l), cll(l), c21(l), _
c31(l), cvl(l), c!3(l), c23(l), c33(l), cv3(l), c!4(l), c24(l),_
c34(l), cv4(l), SerNo(l), Tslope(l), Tint(l), TSnl(l), chla(l), _
chlb(l), chlc(l), TSn2(l), ch2a(l), ch2b(l), ch2c(l), TSn3(l), _
ch3a(l), ch3b(l), ch3c(l), TSn4(l), ch4a(l), Ch4b(l), Ch4c(l), _
TSn5(l), ch5a(l), Ch5b(l), Ch5c(l), TSn6(l), Ch6a(l), Ch6b(l), _
Ch6c(l), TSn7(l), Ch7a(l), Ch7b(l), Ch7c(l)
' Import calibration data
Open "Calibration_Data(2).csv" For Input As #7
1 = 1
Do While Not EOF(7)
1 MsgBox (I)
Input #7, Resistor(l), clO(l), c20(l), c30(l), cvO(l), cll(l), c21(l),
c31(l), cvl(l), c!3(l), c23(l), c33(l), cv3(l), c!4(l), c24(l), _
c34(l), cv4(l), SerNo(l), Tslope(l), Tint(l), TSnl(l), chla(l), _
chlb(l), chlc(l), TSn2(l), ch2a(l), ch2b(l), ch2c(l), TSn3(l), _
ch3a(l), ch3b(l), ch3c(l), TSn4(l), ch4a(l), Ch4b(l), Ch4c(l), _
TSn5(l), ch5a(l), Ch5b(l), Ch5c(l), TSn6(l), Ch6a(l), Ch6b(l), _
Ch6c(l), TSn7(l), Ch7a(l), Ch7b(l), Ch7c(l)
1 = 1 + 1
Loop
1 MsgBox (I)
Close #7
ChDirDirl
'Determine the number of records in the data file
Open Filel For Input As #1
Nl = 0
N = 0
Do While Not EOF(l)
Line Input #1, textline
N = N + 1 'n is the number of records or values
Nl = Nl + 1 'Nl is the number of records or values do not reset this counter
Loop
N = N-1
103
-------
N1=N1-1
Close #1
1 MsgBox ("Number of points = "& Nl)
lf(Nl-20<0)Then
MsgBox "Not enough data points. Please check file."
GoTo 44
End If
'Allocate space to hold arrays
i********
j = 6 'j is the number of channels of data in the raw data file
ReDim timeO(N + 1) As String '+1?
ReDimt(Nl)
ReDim tb(Nl)
ReDim ch(j, Nl) As Double
ReDim chb(j, Nl)
Open Filel For Input As #5
j=0
N = 0
Do While Not EOF(5)
lfj = OThen
Input #5, timeO(j), bline, heater, tail, SerNumber, DomeDia, site
End If
j=j + l
N = N + 1
Input #5, timeO(j), tb(j), chb(l, j), chb(2, j), chb(3, j), chb(4, j), chb(5, j), chb(6, j)
t(j) = CDbl(tb(j))
ch(l,j) = CDbl(chb(l,j))
ch(2,j) = CDbl(chb(2,j))
ch(3,j) = CDbl(chb(3,j))
ch(4,j) = CDbl(chb(4,j))
ch(5,j) = CDbl(chb(5,j))
ch(6,j) = CDbl(chb(6,j))
'MsgBox ("t(j) = "&t(j))
Loop
Close #5
Tori = lTo6
' Forii = lTok
' ch(i, ii) = CDbl(chb(i, ii)) 'convert byte data to double precision
MsgBox("ch(i,ii) = "&ch(i, ii))
' Next ii
'Nexti
'For ii = ITo k
' t(ii) = CDbl(tb(ii)) ' convert byte data to double precision
1 MsgBox("t(ii) = "
'Next ii
IblSiteD. Caption = site
'MsgBox (ch(l, 1))
'13 j=j-l
104
-------
' Find serial Number
For i = 2 To I - 1
If SerNumber - SerNo(i) = 0 Then GoTo 100
Next i
If i = I Then SerNo(i - 1) = "default" ' Serial number not detected
i = i-l
100 IblSerNo. Caption = SerNo(i)
ii = i ' ii is used as the record number of the for the calibration data
' MsgBox (ii)
' Convert Data to temperatures
Dim temp '
ReDim Res(6, Nl)
For i = 1 To Nl 'Nl is the number of data records
Res(l, i) = ch(l, i)
Res(2, i) = Res(l, i) - (ch(2, i) / 10)
Res(3, i) = Res(l, i) - (ch(3, i) / 10)
Res(4, i) = Res(l, i) - (ch(4, i) / 10)
Res(5, i) = Res(l, i) - (ch(5, i) / 10)
Res(6, i) = Res(l, i) - (ch(6, i) / 10)
Res(l, i) = ch(l, i) - 1820
Nexti
' Now convert to temperature and return as the original variable
Fori = lToNl
If Res(l, i)>OThen
temp = chla(ii) + chlb(ii) * Log(Res(l, i) / 10000)
temp = temp + chlc(ii) * ((Log(Res(l, i) / 10000)) A 3)
ch(l, i) = (1 / temp) - 273.15 ' 273.15 convert from Kelvin to Celsus
Else:ch(l, i) = 0
End If
If Res(2, i)>OThen
temp = ch2a(ii) + ch2b(ii) * Log(Res(2, i) / 10000)
temp = temp + ch2c(ii) * ((Log(Res(2, i) / 10000)) A 3)
ch(2, i) = (I/ temp) -273.15
Else:ch(2, i) = 0
End If
If Res(3, i)>OThen
temp = ch3a(ii) + ch3b(ii) * Log(Res(3, i) / 10000)
temp = temp + ch3c(ii) * ((Log(Res(3, i) / 10000)) A 3)
ch(3, i) = (I/ temp) -273.15
Else:ch(3, i) = 0
End If
If Res(4, i)>OThen
temp = ch4a(ii) + Ch4b(ii) * Log(Res(4, i) / 10000)
temp = temp + Ch4c(ii) * ((Log(Res(4, i) / 10000)) A 3)
ch(4, i) = (I/ temp) -273.15
Else:ch(4, i) = 0
105
-------
End If
If Res(5, i)>OThen
temp = ch5a(ii) + Ch5b(ii) * Log(Res(5, i) / 10000)
temp = temp + Ch5c(ii) * ((Log(Res(5, i) / 10000)) A 3)
ch(5, i) = (I/ temp) -273.15
Else:ch(5, i) = 0
End If
If Res(6, i)>OThen
temp = Ch6a(ii) + Ch6b(ii) * Log(Res(6, i) / 10000)
temp = temp + Ch6c(ii) * ((Log(Res(6, i) / 10000)) A 3)
ch(6, i) = (I/ temp) -273.15
Else:ch(6, i) = 0
End If
Nexti
i*****
' get supplemental data
File2 = Replace(Filel, ".dat", "S.dat")
1 MsgBox (File2)
Open File2 For Input As #1
Input #1, timea, ch7, ch8, ch9, chlO
Close #1
lfch7>OThen
temp = Ch7a(ii) + Ch7b(ii) * Log(ch7 / 10000)
temp = temp + Ch7c(ii) * ((Log(ch7 / 10000)) A 3)
ch7 = (1 / temp) - 273.15 ' 273.15 convert from Kelvin to Celsus
Else:ch7 = 0
End If
' Convert suplimental frequencies from ch8 through ch 10 to voltages
ch8 = ch8 / 10000
ch9 = ch9 / 10000
ch 10 = ch 10/10000
' Determine if sufficient data for possible analysis
14 If (t(Nl - 1) - bline * 60 < 0) Then
MsgBox ("t(Nl-l) = " & t(Nl - 1))
MsgBox "Not enough data points. Check Initial Parameters."
GoTo 44
End If
' Convert temperatures from actual temperature in deg C to differential temperature in Deg C
Fori = lToNl
For jj = 2 To 6
ch(jj, i) = ch(jj, i) - ch(l, i)
106
-------
Next jj
Next i
' Analyze base line data and determine slope and intercept of the base line
Dim offset
offset = bline * 60 + ((5 * Gatelnterval& / 1000) + Gatelntervall& / 1000) ' Note:
' Gatelnterval is in milliseconds calculations are in seconds
Do
jj = jj + 1 ' JJ is sample number before heater is turned on
Loop Until (t(jj) - offset) >= 0
jj = jj - 2
'calculate sums used to calculate baseline averages, slopes, and intercepts
For j = 1 To 6
sumx(j) = 0
sumy(j) = 0
sumxy(j) = 0
sumxx(j) = 0
Nextj
' Determind average values
For ct = 1 To jj
For j = 1 To 6
sumx(j) = sumx(j) + t(ct)
sumy(j) = sumy(j) + ch(j, ct)
Nextj
Next ct
1 ct = ct + 1
1 Loop Until (t(ct) - bline * 60) > 0
For j = 1 To 6
avgx(j) = sumx(j)/jj
avgy(j) = sumy(j)/jj
Nextj
Dim vail As Double
Dim va!2 As Double
For ct = 1 To jj
For j = 1 To 6
vail = (t(ct) - avgx(j)) * (ch(j, ct) - avgy(j))
va!2 = (t(ct) - avgx(j)) A 2
sumxy(j) = sumxy(j) + vail
sumxx(j) = sumxx(j) + va!2
Nextj
Next ct
Calculate Slopes and Intercepts
For j = 1 To 6
slope(j) = sumxy(j) / sumxx(j)
slnt(j) = avgy(j) - (slope(j) * avgx(j))
107
-------
Next]
IblTemp. Caption = Format(avgy(l), "#0.0")
LblBLSIope.Caption = Format(slope(l) * 3600, "O.OE-00"
LblBLInt.Caption = Format(slnt(l), "#0.0")
'Correct for baseline drift
For i = 1 To N
For j = 1 To 6
ch(j, i) = ch(j, i) - (slope(j) * t(i) + slnt(j))
Next]
Next i
' find time at maximum signal
' the time is calculated as max time - base line time - (pulse time/2)
' This will be one of the metrics evaluated to determine travel time
Dim jmax(6) As Integer 'sample number at max temperature
For i = !To6
tempmax(i) = ch(i, jj + 1)
Next i
For i = !To6
For j = jj + 2 To N
If (ch(i, j) - tempmax(i) > 0) Then
tempmax(i) = ch(i, j)
jmax(i)=j
End If
Next j
Nexti
tmax(l) = t(jmax(l)) + Gatelntervall& / 2000
For i = 2 To 6
If tempmax(i) < 0.075 Then ' 0.075 considered minimum peak to be used eliminates some of the noise
tmax(i) = bline * 60 + (heater / 2)
Else: tmax(i) = t(jmax(i)) + (Gatelnterval& / 1000 * (i - 1)) + Gatelnterval& / 2000 + Gatelntervall& / 1000
End If
Next i
' display the peak arrival time
Ibltmaxl.Caption = Format(tmax(2) - (bline * 60 + heater / 2), "#0.00")
Ibltmax2. Caption = Format(tmax(3) - (bline * 60 + heater / 2), "#0.00")
IbltmaxS. Caption = Format(tmax(4) - (bline * 60 + heater / 2), "#0.00")
lbltmax4.Caption = Format(tmax(5) - (bline * 60 + heater / 2), "#0.00")
IbltmaxS. Caption = Format(tmax(6) - (bline * 60 + heater / 2), "#0.00")
'Begin Moment analysis - Data will be extrapolated based on a lograthmetic
' function based on a sum of squares fit to the data between 80% of tempmax
' and next occurance of <= 0.04 deg C or the end of file which ever comes first
' first determine sample number correlated to 80%of tempmax to perform moment analysis
' there should be only one peak and it should occure in the first 40% of the data set
' if peak temperature occures after the 40% of the data points do not make moment analysis
108
-------
I********
Dim jbegin(6) As Integer 'jbegin is the sample number for the first sample to use in generating fit
Dim jend(6) As Integer 'jend is the last sample to use
If ((0.4 * Nl) - jmax(l) < 0) And ((0.4 * Nl) - jmax(2) < 0) And ((0.4 * Nl) - jmax(3) < 0) _
And ((0.4 * Nl) - jmax(4) < 0) And ((0.4 * Nl) - jmax(5) < 0) Then
lblMom2. Caption = Format("N.A.")
IblMomS. Caption = Format("N.A.")
IblMomS. Caption = Format("N.A.")
IblMomS. Caption = Format("N.A.")
GoTo 200
End If
For i = 2 To 6
If (jmax(i) + 20 - Nl) >OThen 'See if there is sufficient data to calculate moments
GoTo 90 ' changed from 200 2/24/08
End If
j=jmax(i)
Do Until ((tempmax(i) * 0.8) - ch(i, j) > 0) 'determine sample number for first sample
' to use in enterpolating and extrapolating
' moments 80% of peak was qualitatively
' estimated as being sufficient to be away
' from the curvature at the peek
j=j + l
If (j - N = 0) Then
jbegin(i)=j
jend(i)=j
GoTo 90
End If
' jend(i) = jend(i) + 1
Loop
jbegin(i)=j
lf(jbegin(i) + 20-Nl)>OThen
GoTo 90
End If
j=j + l
lfj>=N!ThenGoTo89
Do Until (0.05 - ch(i, j) > 0) 'determine last sample number considered valid for
'interpolating or extrapolating moments noise was measured
'at +- 2 mv. data has already been adjusted for baseline
'offset
lfj>=N!ThenGoTo89
Loop
89 jend(i)=j-l
90 Next i
' Fit log transformed data to a streight line first determine if there is enough data
' there should be a minimum of minobs observations to calculate the extrapolation function
Dim sumyy(6)
Dim avgxx(6)
Dim avgyy(6)
109
-------
Dim Mslope(6) As Double
Dim Mint(6) As Double
Dim minobs As Integer 'minimum number of observations between jbegin(i) and jend(i)
'to permit moment analysis 20 is considered minimm
Dim minsignal As Single 'minimum peak signal in volts for analysis
minobs = 20
minsignal = 0.05
Fori = 2To6 'initialize values
sumxx(i) =0
sumyy(i) = 0
sumxy(i) = 0
avgxx(i) = 0
avgyy(i) = 0
Mslope(i) = 0
Mint(i) = 0
Nexti
For i = 2 To 6
If ((jend(i) - jbegin(i) - minobs) > 0) And (tempmax(i) - minsignal > 0) Then
For j = jbegin(i) To (jend(i) - 1)
kk = j + 1
avgxx(i) = avgxx(i) + t(kk)
If ch(i, kk) < 0 Then ch(i, kk) = 1E-50
avgyy(i) = avgyy(i) + Log(ch(i, kk))
Next]
avgxx(i) = avgxx(i) / (jend(i) - jbegin(i))
avgyy(i) = avgyy(i) / (jend(i) - jbegin(i))
End If
Nexti
For i = 2 To 6
If ((jend(i) - jbegin(i) - minobs) > 0) And (tempmax(i) - minsignal > 0) Then
For j = jbegin(i) To (jend(i) - 1)
kk = j+l
If ch(i, kk) < 0 Then ch(i, kk) = 1E-50
sumxx(i) = sumxx(i) + (t(j) - avgxx(i)) A 2
sumyy(i) = sumyy(i) + (Log(ch(i, kk)) - avgyy(i)) A 2
sumxy(i) = sumxy(i) + (t(kk) - avgxx(i)) * (Log(ch(i, kk)) - avgyy(i))
Next]
Mslope(i) = sumxy(i) / sumxx(i)
Mint(i) = avgyy(i) - (Mslope(i) * avgxx(i))
End If
Nexti
'Calculate zero and first moments
Dim Momzero(6)
Dim Moml(6)
Dim cl, c2
Dim Tl, T2
Dim timecor 'time correction such that zero time is at the end of base line measurements
Dim jstart 'sample number to begin moment calculations this is to remove switching noise
timecor = bline * 60
For i = 2 To 6
110
-------
If ((jend(i) - jbegin(i) - minobs) < 0) Or (tempmax(i) - minsignal < 0 Or Mslope(i) > 0) Then
Moml(i) = "N.A."
GoTo 300
Else
Momzero(i) =0
Moml(i) = 0
ch(i, 0)=0
t(0)=0
Start moment analysis beginning with sample after heat-pulse
Correct time scale to begin half way through the heat-pulse in moment calculations
First step is to determine sample number after end of heat-pulse
Do Until (t(j) - bline * 60 - heater / 2) > 0
j=j + l
Loop
jstart = j
For j = jstart To jbegin(i) 'first part of curve
Momzero(i) = Momzero(i) + ((t(j) - t(j -1)) * (ch(i, j) + ch(i, j -1))) / 2
Moml(i) = Moml(i) + (1 / 6) * ((ch(i, j -1) * (t(j -1) - timecor) + ch(i, j) * _
(t(j)-timecor)) *(t(j)-t(j-l)) + (ch(i,j-l)_
+ ch(i, j)) * ((t(j) - timecor) * (t(j) - timecor) - (t(j -1) - timecor) * _
(t(j-l)-timecor)))
Next]
'end of first part of curve begin analysis for the curve fit portion of the curve
d = ch(i,j)
Tl = t(j)
j = l
Do Until (cl - 0.001 * tempmax(i) < 0)
T2 = Tl + 1
c2 = Exp(Mslope(i) * (T2) + Mint(i))
Momzero(i) = Momzero(i) + (cl + c2) / 2
Moml(i) = Moml(i) + (1 / 6) * ((cl * (Tl - timecor) + c2 * (T2 - timecor)) * _
(T2-Tl) + _
(cl + c2) * ((T2 - timecor) * (T2 - timecor) - (Tl - timecor) * (Tl - timecor)))
T1 = T2
cl = c2
Loop
End If
300 Next i
'MsgBox(Moml(l))
If ((jend(2) - jbegin(2)) > minobs) And (tempmax(2) - minsignal > 0 And lsNumeric(Moml(2))) Then
lblMom2.Caption = Format(Moml(2) / Momzero(2), "#0.0")
Else
lblMom2. Caption = Format("N.A.")
End If
If ((jend(3) - jbegin(3)) > minobs) And (Not Not tempmax(S) - minsignal > 0 And lsNumeric(Moml(3))) Then
IblMomS.Caption = Format(Moml(3) / Momzero(3), "#0.0")
Else
111
-------
IblMomS. Caption = Format("N.A.")
End If
If ((jend(5) - jbegin(5)) > minobs) And (tempmax(S) - minsignal > 0 And lsNumeric(Moml(5))) Then
IblMomS. Caption = Format(Moml(5)/ Momzero(5), "#0.0")
Else
IblMomS. Caption = Format("N.A.")
End If
If ((jend(6) - jbegin(6)) > minobs) And (tempmax(6) - minsignal > 0 And lsNumeric(Moml(6))) Then
IblMomS. Caption = Format(Moml(6) / Momzero(6), "#0.0")
Else
IblMomS. Caption = Format("N.A.")
End If
' generate labels for flow
DimCL
Dim maxtime
Dim CLV
CL = 2.326
' Channel 2 Captions
If ((jend(2) - jbegin(2) > minobs) And (tempmax(2) - minsignal > 0) And _
(tmax(2) - bline * 60 - 1.5 * heater > 0) And lsl\lumeric(Moml(2))) Then
maxtime = tmax(2) - (bline * 60 + heater / 2)
CLV = CL * cvO(ii)
QO = (clO(ii) - (c30(ii) * maxtime)) / (maxtime + c20(ii))
QOmin = (clO(ii) - (c30(ii) * maxtime * (1 + CLV))) /_
(maxtime *(! + CLV) + c20(ii))
QOmax = (clO(ii) - (c30(ii) * maxtime * (1 - CLV))) /_
(maxtime *(1- CLV) + c20(ii))
IblQO.Caption = Format(QO, "#0.00")
IblQOmax.Caption = Format(QOmax, "#0.00")
IblQOmin. Caption = Format(QOmin, "#0.00")
Else
QO = "N.A."
QOmin = "N.A."
QOmax = "N.A."
IblQO.Caption = Format(QO)
IblQOmin. Caption = Format(QOmin)
IblQOmax.Caption = Format(QOmax)
End If
' Channel 3 Captions
If (jend(3) - jbegin(3) > minobs) And (tempmax(3) - minsignal > 0) And (tmax(3) - bline * 60 - 1.5 * heater > 0
And lsl\lumeric(Moml(3)))Then
maxtime = tmax(3) - (bline * 60 + heater / 2)
CLV = CL * cvl(ii)
Ql = (cll(ii) - (c31(ii) * (maxtime))) / (maxtime + c21(ii))
Qlmin = (cll(ii) - (c31(ii) * maxtime * (1 + (CLV)))) / _
(maxtime *(! + CLV) + c21(ii))
Qlmax = (cll(ii) - c31(ii) * (maxtime) * (1 - (CLV))) / _
(maxtime *(1- CLV) + c21(ii))
IblQl.Caption = FormatfQl, "#0.00")
IblQlmax.Caption = Format(Qlmax, "#0.00")
112
-------
IblQlmin. Caption = Format(Qlmin, "#0.00")
Else
Q1 = "N.A."
Qlmin = "N.A."
Qlmax= "N.A."
IblQl.Caption = Format(Ql)
IblQlmin. Caption = Format(Qlmin)
IblQlmax. Caption = Format(Qlmax)
End If
Channel 5 Captions
If (jend(5) - jbegin(5) > minobs) And (tempmax(5) - minsignal > 0) And _
(tmax(5) - bline * 60 - 2 * heater > 0 And lsNumeric(Moml(5))) Then
maxtime = tmax(5) - (bline * 60 + heater / 2)
CLV = CL * cv3(ii)
Q3 = (c!3(ii) - (c33(ii) * maxtime)) / (maxtime + c23(ii))
QSmin = (c!3(ii) - c33(ii) * maxtime * (1 + CLV)) / _
(maxtime *(! + CLV) + c23(ii))
QSmax = (c!3(ii) - c33(ii) * maxtime * (1 - CLV)) / _
(maxtime *(1- CLV) + c23(ii))
IblQS.Caption = Format(Q3, "#0.00")
IblQSmax.Caption = Format(Q3max, "#0.00")
IblQSmin. Caption = Format(Q3min, "#0.00")
Else
Q3 = "N.A."
QSmin = "N.A."
QSmax = "N.A."
IblQS.Caption = Format(QS)
IblQSmin. Caption = Format(QSmin)
IblQSmax.Caption = Format(QSmax)
End If
Channel 6 Captions
If (jend(6) - jbegin(6) > minobs) And (tempmax(6) - minsignal > 0) And _
(tmax(6) - bline * 60 - 1.5 * heater > 0 And lsNumeric(Moml(6))) Then
maxtime = tmax(6) - (bline * 60 + heater / 2)
CLV = CL * cv4(ii)
Q4 = (c!4(ii) - (c34(ii) * maxtime)) / (maxtime + c24(ii))
Q4min = (c!4(ii) - c34(ii) * (maxtime) * (1 + CLV)) /_
(maxtime *(! + CLV) + c24(ii))
Q4max = (c!4(ii) - c34(ii) * (maxtime) * (1 - CLV)) /_
(maxtime *(1- CLV) + c24(ii))
lblQ4.Caption = Format(Q4, "#0.00")
lblQ4max.Caption = Format(Q4max, "#0.00")
lblQ4min. Caption = Format(Q4min, "#0.00")
Else
Q4 = "N.A."
Q4min = "N.A."
Q4max = "N.A."
lblQ4.Caption = Format(Q4)
lblQ4min. Caption = Format(Q4min)
lblQ4max.Caption = Format(Q4max)
End If
' Generate Darcy velocity lables
113
-------
I*********
If DomeDia > 0.1 Then 'Data exists to calculate flux
lfQOo"N.A."Then
VO = QO * 60 * 24 / (3.1415 * DomeDia * DomeDia / 4)
IblVO.Caption = Format(VO, "#0.0")
Else
VO = "N.A."
IblVO.Caption = Format(VO)
End If
lfQlo"N.A."Then
VI = Ql * 60 * 24 / (3.1415 * DomeDia * DomeDia / 4)
IblVl.Caption = FormatfVl, "#0.0")
Else
V1 = "N.A."
IblVl.Caption = Format(Vl)
End If
lfQ3o"N.A."Then
V3 = Q3 * 60 * 24 / (3.1415 * DomeDia * DomeDia / 4)
IblVS.Caption = FormatfVS, "#0.0")
Else
V3 = "N.A."
IblVS.Caption = Format(VS)
End If
lfQ4o"N.A."Then
V4 = Q4 * 60 * 24 / (3.1415 * DomeDia * DomeDia / 4)
lblV4.Caption = Format(V4, "#0.0")
Else
V4 = "N.A."
lblV4.Caption = Format(V4)
End If
Else
VO = "N.A."
V1 = "N.A."
V3 = "N.A."
V4 = "N.A."
IblVO.Caption = Format(VO)
IblVl.Caption = Format(Vl)
IblVS.Caption = Format(VS)
lblV4.Caption = Format(V4)
End If
'MsgBox ("checks")
' set up arrays for ploting
200 Dim Graph() As Single
Dim x As Integer
Dim myXarrayO As Double Time
Dim myYArrayO() As Double
Dim myYArraylO As Double
Dim myYArray2() As Double
Dim myYArrayS() As Double
114
-------
Dim myYArray4() As Double
Dim myYArraySQ As Double
Re Dim
Re Dim
Re Dim
Re Dim
Re Dim
Re Dim
Re Dim
myXarray(N)
myYArrayO(N)
myYArrayl(N)
myYArray2(N)
myYArrayS(N)
myYArray4(N)
myYArrayS(N)
'Generate some x y data.
myXarray(O) = t(l)
myYArrayO(O) = ch(l, 1)
'MsgBox (ch(2,1))
myYArrayl(O) = ch(2,1)
myYArray2(0) = ch(3,1)
myYArrayS(O) = ch(4,1)
myYArray4(0) = ch(5,1)
myYArrayS(O) = ch(6,1)
'MsgBox N
For x = 1 To N -1
myXarray(x) = t(x) - bline * 60 'value for X-axis
myYArrayO(x) = ch(l, x)
myYArrayl(x) = ch(2, x)
myYArray2(x) = ch(3, x)
myYArrayS(x) = ch(4, x)
myYArray4(x) = ch(5, x)
myYArrayS(x) = ch(6, x)
Next x
'value for Y-axis RefTemp
' Ch 2 cable end temp
'Ch 3 temp
' Heater Temp
'Discharge end temp
With TChartl
.AddSeries scPoint
.Series(0).AddArray UBound(myYArrayO),
.Series(0).XValues.TempValue = True
.Series(l).AddArray UBound(myYArrayl),
.Series(2).AddArray UBound(myYArray2),
.Series(3).AddArray UBound(myYArrayS),
.Series(4).AddArray UBound(myYArray4),
.Series(5).AddArray UBound(myYArray5),
End With
'MsgBox ("check 4")
Screen. MousePointer = 1
myYArrayO(), myXarray()' ref temp
Cable end TEp
myYArrayl(), myXarray()
myYArray2(), myXarray()
myYArray3(), myXarray()
myYArray4(), myXarray()
myYArray5(), myXarray()
Heater
' write output file with processed data
115
-------
I***********
ChDirDirl
Dim fname2 ' file name for output data
'MsgBox ("timeO(O) = " & timeO(O))
fname2 = Format$(timeO(0), "mm dd hh mm") + "a.CSV"
Open fname2 For Output As #8
Write #8, Format("baseline"), Format("heater"), Format("Serial No."), Format("Dome D"), Format("site")
Write #8, Formatfbline, "#0.0"), _
Format(heater, "#0.0"), Format(SerNumber, "#0"), _
Format(DomeDia, "#0.0000"), Format(site)
Write #8, Format("Channel No"), Format("Flow Rate ml/min"), Format("Lower 98% Limit"), _
Format("Upper 98%Limit"), Format("Darcy Velocity cm/day")
Write #8, Format("2"), Format(QO, "#0.00"), Format(QOmin, "#0.00"), Format(QOmax, "#0.00"), Format(VO,
"#0.00")
Write #8, Format("3"), Format(Ql, "#0.00"), Format(Qlmin, "#0.00"), Format(Qlmax, "#0.00"), Format(Vl,
"#0.00")
Write #8, Format("5"), Format(Q3, "#0.00"), Format(Q3min, "#0.00"), Format(Q3max, "#0.00"), Format(V3,
"#0.00")
Write #8, Format("6"), Format(Q4, "#0.00"), Format(Q4min, "#0.00"), Format(Q4max, "#0.00"), Format(V4,
"#0.00")
Write #8, Format("Peak arrival time in seconds from cable end")
Write #8, Format(tmax(2) - (bline * 60 + heater / 2), "#0.0"), Format(tmax(3) - (bline * 60 + heater / 2), "#0.0")
Format(tmax(4) - (bline * 60 + heater / 2), "#0.0"), Format(tmax(5) - (bline * 60 + heater / 2), "#0.0"), _
Format(tmax(6) - (bline * 60 + heater / 2), "#0.0")
Write #8, Format("Supplemental Data Temperature Outside Dome Ch7, Ch8, Ch9, ChlO")
Write #8, Format("Ch7 looks for calibration coefficients while Ch8 - chlO are equivalent voltages")
Write #8, Format("Ch7"), Format("Ch8"), Format("Ch9"), Format("ChlO")
Write #8, Format(ch7, "#0.0"), Format(ch8, "#0.0"), Format(ch9, "#0.0"), FormatfchlO, "#0.0")
Write #8, Format("Maximum differential temperature or reference temperature c")
Write #8, Format("Ch 1 Ref"), Format("Ch 2"), Format("Ch 3"), Format("Ch 4"), Format("Ch 5"), _
Format("Ch 6"), Format("Ref Temp C")
Write #8, Format(tempmax(l), "#0.00"), Format(tempmax(2), "#0.00"), Format(tempmax(3), _
"#0.00"), Format(tempmax(4), "#0.00"), Format(tempmax(5), "#0.00"), Format(tempmax(6), "#0.00"), _
Format(avgy(l), "#0.")
Write #8, Format("Normalized First Moments Calculated beginning at the end of the base line samples")
If ((jend(2) - jbegin(2)) > minobs) And (tempmax(2) - minsignal > 0 And _
Moml(2) <> 0 And lsNumeric(Moml(2))) Then
Write #8, Format("Ch 2"), Format(Moml(2) / Momzero(2), "#0.0")
Else
Write #8, Format("Ch 2"), Format("N.A.")
End If
If ((jend(3) - jbegin(3)) > minobs) And (tempmax(S) - minsignal > 0 And _
Moml(3) <> 0 And lsNumeric(Moml(3))) Then
Write #8, Format("Ch 3"), Format(Moml(3) / Momzero(3), "#0.0")
Else
Write #8, Format("Ch 3"), Format("N.A.")
End If
If ((jend(4) - jbegin(4)) > minobs) And (tempmax(4) - minsignal > 0 And _
Moml(4) <> 0 And lsNumeric(Moml(4))) Then
Write #8, Format("Ch 4"), Format(Moml(4) / Momzero(4), "#0.0")
Else
116
-------
Write #8, Format("Ch 4"), Format("NA")
End If
If ((jend(5) - jbegin(5)) > minobs) And (tempmax(S) - minsignal > 0 And _
Moml(5) <> 0 And lsNumeric(Moml(5))) Then
Write #8, Format("Ch 5"), Format(Moml(5) / Momzero(5), "#0.0")
Else
Write #8, Format("Ch 5"), Format("N.A.")
End If
If ((jend(6) - jbegin(6)) > minobs) And (tempmax(6) - minsignal > 0 And _
Moml(6) <> 0 And lsNumeric(Moml(6))) Then
Write #8, Format("Ch 6"), Format(Moml(6) / Momzero(6), "#0.0")
Else
Write #8, Format("Ch 6"), Format("N.A.")
End If
Write #8, Format("Times in seconds from beginning of heat-pulse and Rel Temperatures in Deg C")
Write #8, Format("E.T."), Format("T-Ref-l"), Format("T-2"), Format("T-3"), Format("T-4"), Format("T-5"),
Format("T-6")
For x = 1 To N -1
Write #8, Format(myXarray(x), "#0.000"), Format(myYArrayO(x), "#0.000"), Format(myYArrayl(x), "#0.000"),
Format(myYArray2(x), "#0.000"), Format(myYArray3(x), "#0.000"), Format(myYArray4(x), "#0.000"), _
Format(myYArray5(x), "#0.000")
Nextx
Close
Refresh
DoEvents
GoTo 45
44 Unload Me
Screen.MousePointer = 1
45 End Sub
Private Sub TeeCommander2_OnEditedChart()
End Sub
Private Sub Mom4_Click()
End Sub
117
-------
Appendix B Circuit Diagram of Measurement System
••.•::
r
• •
. .
. -, -•
... ., . *.. »
w.v*»- -
:-• ' "
• 01 Ht
n, -
. . .
- i
•
118
-------
119
-------
Appendix C Signal Conditioning Circuit Board Layout
120
-------
Appendix D Interface Circuit Board Layout
• jn*' vf^; _y\ .? *—?»
dlrSaR
121
-------
Appendix E Calibration constants of flow sensor SN300, SN303
SN300
SN303
T2
T3
T5
T6
T2
T3
T5
T6
Ci
-674.98
-380.23
370.28
562.70
Ci
-603.84
-420.85
379.77
695.78
C2
-14.27
-14.19
-17.20
-25.58
C2
-24.61
-17.45
-17.03
-14.35
C3
0.081
-0.15
0.46
-0.12
C3
0.11
-0.96
0.69
-0.10
122
------- |