function [ v ] = plot_rayleigh_bfun( ind, kreal )

global model nfreq nbfun n2 z lambda bulk myu l2myu

[ jr1 r1j jr2 r2j il1 il2 ] = rayleigh_basis( ind, kreal );

z1 = [ 0:0.1:z(2) ];
f = z1;
g = z1;

jr1
jr2
il1
il2

clf;

f1 = zeros(size(z1));
f2 = zeros(size(z1));
g1 = zeros(size(z1));
g2 = zeros(size(z1));

if ( jr1 > -1 )
	for i=1:length(z1)
		f1(i) += r1j*bfunval( jr1, z1(i) );
		g1(i) += r1j*bfungrad( jr1, z1(i) );
	endfor
endif

if ( jr2 > -1 )
	for i=1:length(z1)
		f2(i) += r2j*bfunval( jr2, z1(i) );
		g2(i) += r2j*bfungrad( jr2, z1(i) );
	endfor

endif

tauz = l2myu(1)*g2 + kreal*lambda(1)*f1;
taux = g1 - kreal*f2;

figure(1);

plot(z1,f1,'r-;r1;',z1,g1,'b-;r1 grad;',z1,f2,'g-;r2;',z1,g2,'k-;r2 grad;');

figure(2);

plot(z1,tauz,'-;tauz;',z1,taux,'-;taux;');

# axis ([0 100 -5 10]);

v = 0;

endfunction