; this function is designed to return an array
; of the mass concentrations/ft depth of a
; given chemical for a given cell
; three variables are passed
; to the procedure, c_id, s_time, and chemical
; the array gives the mass and elevation for the cell
function cont_mass, c_id, s_time, chemical
;remember to put each variable in quotes
; c_id = a one character designation for the cell
; s_time = sampling time pre or post
; chemical = chemical abreviation which you want to
; determine the contaminant mass
env_handle=ODBC_INIT()
;connect to database
access_id=ODBC_CONNECT(env_handle,'MS Access 97 Database','Admin/')
cell = string(c_id, Format='(A1)')
;print, cell
;print, s_time
;print, chemical
;
;data required to make the sample have any value
;one needs to know at a minimum;
;field_vial_water_mecl_wt, and field_vial_water_mecl_solid_wt
;we will assume that concentrations that are not reported
;are nulls can be considered to be 0.0001 a number above
;zero was selected to permit doing a log transformation
;
data_sel = "SELECT northing,easting,elevation,sample_top,sample_bottom"+$
" ,lab_vial_water_wt, lab_vial_water_mecl_wt, "+$
" field_vial_water_mecl_wt,field_vial_water_mecl_solid_wt,"+$
" wet_wt_solid,dry_wt_solid, sample_category"
data_sel = data_sel + ", " + chemical
data_sel = data_sel + " FROM tblObjects_locations "+$
" INNER JOIN (tblCore_samples INNER JOIN tblCore_sample_chemistry "+$
" ON tblCore_samples.sample = tblCore_sample_chemistry.sample) ON (tblObjects_locations.key = "+$
" tblCore_samples.tblkey) AND (tblObjects_locations.key = tblCore_sample_chemistry.tblkey) "+$
" where "+$
" elevation is not null and field_vial_water_mecl_wt is not null"+$
" and field_vial_water_mecl_solid_wt is not null"
data_sel = data_sel + ' and tblCore_sample_chemistry.cell_id = ' + "'" + cell +"'"
data_sel = data_sel + ' and sample_category = '
data_sel = data_sel + "'"+s_time+"'"
;print, data_sel
borings=odbc_sql(access_id, data_sel)
on_error_GOTO, finish
if N_Elements(borings.northing) lt 2 then goto, finish
;put borings data into arrays
boringsarrayx=borings.easting
boringsarrayy=borings.northing
boringsarrayz=borings.elevation-((borings.sample_top + borings.sample_bottom)/2)
boringsarraylww= borings.lab_vial_water_wt
boringsarraylwmw = borings.lab_vial_water_mecl_wt
boringsarrayfwmw = borings.field_vial_water_mecl_wt
boringsarrayfwmsw = borings.field_vial_water_mecl_solid_wt
boringsarraywws = borings.wet_wt_solid
boringsarraydws = borings.dry_wt_solid
;
; assume laboratory weights determine amount of mecl present
; and field weights determine the amount of soil and moisture in soil
;
;
;The chemical concentration is in the 13 field of data so
;use the tag index to eatract the correct field for chemical
;concentration (field-1 or 12)
;
bchem=borings.(12)
;print,bchem(1),boringsarraylww(1),boringsarraylwmw(1),$
;boringsarrayfwmw(1),boringsarrayfwmsw(1),$
;boringsarraywws(1),boringsarraydws(1)
;determine number of observations
xx=size(borings)
lastrow=xx(1)
for i = 0,lastrow - 1 do begin
;print, i
if bchem(i) eq 0.0 then bchem(i)= 0.000100
;print,bchem(i),boringsarraylww(i),boringsarraylwmw(i),boringsarrayfwmsw(i)
if boringsarraylwmw(i) eq 0.0 $
and boringsarraydws(i) ne 0.0 then begin
bchem(i)=$
bchem(i)*(boringsarrayfwmw(i) - $
boringsarraylww(i)) / ((boringsarrayfwmsw(i) - $
boringsarrayfwmw(i))*boringsarraydws(i) $
/boringsarraywws(i))
goto, jump1
endif $
else if boringsarraylwmw(i) eq 0.0 then begin
bchem(i)=$
bchem(i)*(boringsarrayfwmw(i) - $
boringsarraylww(i)) / ((boringsarrayfwmsw(i) - $
boringsarrayfwmw(i))*1.1)
goto, jump1
endif $
else if boringsarraydws(i) ne 0.0 then begin
bchem(i)=$
bchem(i)*(boringsarraylwmw(i) - $
boringsarraylww(i)) / ((boringsarrayfwmsw(i) - $
boringsarrayfwmw(i))*boringsarraydws(i) $
/boringsarraywws(i))
goto, jump1
endif else begin
bchem(i)=$
bchem(i)*(boringsarraylwmw(i) - $
boringsarraylww(i)) / ((boringsarrayfwmsw(i) - $
boringsarrayfwmw(i))*1.1)
endelse ;end of else clause
jump1: ;print,bchem(i)
endfor
;
; this should have created an array of values in bchem which
; is all of the masses on a soil weight basis for the given
; sampling period and given cell for the contaminant of choice
;now we need to
;determine extents of the cell by getting the coordinates of
;the cell we will determine the max and min northing and
;extend to the next whole integer so that we can grid the
;data
;
;
; get the extents of the cell such that an enclosing
; polygon can be drawn
;
data_sel = 'SELECT northing,easting,elevation,object_id,object_type'
data_sel = data_sel + ' FROM tblObjects_locations WHERE '
data_sel = data_sel + ' cell_id = ' + "'" + cell +"' and"
data_sel = data_sel + ' northing is not null and'
data_sel = data_sel + ' easting is not null and'
data_sel = data_sel + " object_type = 'c'"
data_sel = data_sel + ' ORDER BY object_id'
;print, data_sel
cell = odbc_sql(access_id,data_sel)
min_north = long(min(cell.northing))-1
max_north = long(max(cell.northing))+1
min_east = long(min(cell.easting))-1
max_east = long(max(cell.easting))+1
;need to emphasize data at the same elevation due to
;lenticular nature of soils
zscale_factor=5.0
scaled_z=boringsarrayz*zscale_factor
min_elev = long(min(scaled_z))-zscale_factor
max_elev = long(max(scaled_z))+zscale_factor
gridx=50
gridy=50
gridz=nint(max_elev-min_elev)
;info, ival
;info, cell,gridx,gridy,gridz
arrayx=cell.easting
;print, arrayx
arrayx=(arrayx-DOUBLE(min_east))
gridx=DOUBLE(gridx)
scalex=gridx/(DOUBLE(max_east-min_east))
;print, scalex
arrayx=arrayx*scalex
;print, arrayx
arrayy=cell.northing
arrayy=(arrayy-DOUBLE(min_north))
gridy=DOUBLE(gridy)
scaley=gridy/(DOUBLE(max_north-min_north))
;print, scaley
arrayy=arrayy*scaley
wind=0
if s_time eq 'pre' then wind=0
if s_time eq 'post' then wind=1
window,wind,xsize=250,ysize=400
wshow,wind,1
plot, arrayx, arrayy, background = 255, color=0,$
Title=s_time
xx=size(borings)
lastrow=xx(1)
;create an empty array that is (4,n) in size
borings_d=dblarr(4,lastrow)
;place the data into the array remember to start with zero
borings_d(0,*)=(boringsarrayx-DOUBLE(min_east))*scalex
borings_d(1,*)=(boringsarrayy-DOUBLE(min_north))*scaley
borings_d(2,*)=scaled_z
;borings_d(3,*)=alog10(bchem)
borings_d(3,*)=(bchem)
;print,borings_d(0,*)
x=borings_d(0,*)
y=borings_d(1,*)
oplot, x,y, Psym=5, color=0
;print, data
;Scale the data since gridding is just locations
ival=(grid_4D(borings_d,gridx,gridy,gridz,ORDER=4.0, $
ZMIN=min_elev, ZMAX=max_elev))
;print, ival(*,*,8)
;
;calculate the mass per sq ft in each layer
;
;create an empty array that is (2,gridz) in size
mass_d=dblarr(2,gridz)
for i = 0,gridz - 1 do begin
data= ival(*,*,i)
a=data(POLYFILLV(arrayx,arrayy,gridx,gridy))
;print, a
b=data(a)
;b is a scaled value need to put it back in real spatial
;units
b=b/(scalex*scaley)
mass_d(0,i)=total(b)
mass_d(1,i)=FLOAT(min_elev+i)
;print, mass_d(0,i), mass_d(1,i)
endfor
;The mass at this point is actually the average
;concentration calculated per unit thickness. Each layer is one
;gridz thick if we assume that the bulk density of the soil is
;1.5 then it is possible to calculate the mass of chemical in
;per zgrid thickness of the profile.
mass_d(0,*)=mass_d(0,*)*0.1395
mass_d(1,*)=mass_d(1,*)/(zscale_factor)
ODBC_DISCONNECT, access_id
;print,mass_d(0,*)
;print,mass_d(1,*)
return, mass_d
info
finish:
end