; 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, data_c a table containing the concentration

; and coordinates of the observation, cell_b the sorted coordinates

; of the sheet piles, test_n a counter to indicate wether it is a

; pre-remediation test or a post remediation test

; the array gives the mass and elevation for the cell

function cont_mass_2, data_c, cell_b, test_n

;info, cell_b

min_north = long(min(cell_b.northing))-1

max_north = long(max(cell_b.northing))+1

min_east = long(min(cell_b.easting))-1

max_east = long(max(cell_b.easting))+1

;need to emphasize data at the same elevation due to

;lenticular nature of soils

zscale_factor=5.0

scaled_z=data_c.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_b.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_b.northing

arrayy=(arrayy-DOUBLE(min_north))

gridy=DOUBLE(gridy)

scaley=gridy/(DOUBLE(max_north-min_north))

;print, scaley

arrayy=arrayy*scaley

;wind=0

;if test_n eq 0 then wind=0;test_n 0 refers to pre remediation

;if test_n eq 1 then wind=1; test_n 1 refers to post remediation

;window,wind,xsize=250,ysize=400

;wshow,wind,1

;plot, arrayx, arrayy, background = 255, color=0

xx=size(data_c)

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,*)=(data_c.boringsarrayx-DOUBLE(min_east))*scalex

borings_d(1,*)=(data_c.boringsarrayy-DOUBLE(min_north))*scaley

borings_d(2,*)=scaled_z

;borings_d(3,*)=alog10(data_c.bchem)

borings_d(3,*)=(data_c.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=2.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

;2.1 then it is possible to calculate the mass of chemical in

;per zgrid thickness of the profile.

mass_d(0,*)=mass_d(0,*)*0.1954

mass_d(1,*)=mass_d(1,*)/(zscale_factor)

;print,mass_d(0,*)

;print,mass_d(1,*)

return, mass_d

info

finish:

end