# Returns the value of basis function at depth 'zlvl'

function [ v ] = bfunval( ind, zlvl )

global z;
global nlayers;

  layer = 1;
  z0 = 0.0;

  while ( z(layer) < zlvl )

	z0 = z(layer);

	++layer;

	if ( layer >= nlayers )
		v = 0.0;
		return;
	endif

  endwhile

#zlvl
#layer
#ind
#bfunpoly(layer,ind,0)

  x = (zlvl - z0)/(z(layer) - z0);
  v = polyval(bfunpoly(layer,ind,0),x);

endfunction