// howling simulation program for Scilab
//
//  +----------------------------------------------+
//  |                                              |
//  |        x  +-----+ r  +------+ y  +-------+ m |
//  +--->(+)--->| FIR +--->| GAIN +--->| LIMIT +---+
//        ^     +-----+    +------+    +-------+    
//        |
//        |
//        S

clear;
rand('uniform');
rand('seed', 924857);

SIMLEN=300000;
FIRLEN=1000;                  
DLYLEN=1;           
GAIN=3.15;
LIMIT=1.0;
TRIGLVL=0.00001;

out=zeros(1:SIMLEN);
firbuf=zeros(1:FIRLEN+DLYLEN-1);

h=rand(1:FIRLEN-1);
h=h-sum(h)/(FIRLEN-1);
h=[zeros(1:DLYLEN), h];
hspec=fft([h, zeros(1:2000-FIRLEN)], -1);
hspec=20*log10(abs(hspec)+1e-100);
h=h/(10^(max(hspec)/20));

rand('seed', 71);
m=0.0;
for i=1:SIMLEN
  if i<=FIRLEN
    s=(2.0*rand(1)-1.0)*TRIGLVL;
  else
    s=0.0;
  end
  x=m+s;
  firbuf=[x, firbuf(1:FIRLEN+DLYLEN-2)];
  r=h*firbuf';
  y=r*GAIN;  
  if LIMIT<=abs(y)
    m=LIMIT*sign(y);
  else
    m=y;
  end
  out(i)=y;
end

clf;
plot(out);