function wavwrite(sig, sf, bits, wavefile, web)

% WAVWRITE Write Microsoft WAVE (".wav") sound file.
%
%  wavwrite(sig, sf, bits, wavfile);
%  wavwrite(sig, sf, wavfile);
%  wavwrite(sig, wavfile);
%   
%        sig : signal
%         sf : sampling freq. (default : 8000Hz)
%       bits : the number of bits per sample (default : 16bits)
%              ( 8, 16, 24 bits)
%    wavfile : file name 
%              (The ".wav" extension is appended if no extension is given.)

if(nargin < 1)
  help wavwrite;
  return;
end

if(nargin == 2)
  wavefile = sf;
  sf = 8000;
  bits = 16;
end
if(nargin == 3);
  wavefile = bits;
  bits = 16;
end

fmt = 1; % fmt=0: digion format, fmt=1: 2^23 format
if fmt
  cnst = 1;
else
  cnst = 2;
end

if bits != 16 && bits != 8 && bits != 24
  error('bits per sample is not WAVE standard!');
end

if ischar(wavefile)
  if isempty(strfind(wavefile,'.wav')) && isempty(strfind(wavefile, '.WAV'))
    wavefile=[wavefile,'.wav'];
  end
  fid=fopen(wavefile, 'wb', 'ieee-le');
else
  fid = 1;
  if nargin == 5
    fprintf(1, 'Content-type: audio/x-wav\n');
  end
end

if fid == -1
  error('Can''t open .WAV file for output!');
end;

format(1) = 1;
[n, m]=size(sig);
invs = 0;
if(n > m)
  invs = 1;
  t = n; n = m; m = t;
end
if(n > 2)
  fprintf(2, 'Caution!! More than 2-channels is not supported\n');
%  n = 2;
end
format(2) = n;
format(3) = sf;
format(4) = bits;

nsamples = m;
if (n > format(2))
  error('Can put a binaural vector into a Wave Mono sound file!'); 
end
if format(4) == 8
  sig=round((sig+1) * 255/2);
  Block = n;
end
if format(4) == 16
%  sig = round((sig+1)*65535/2) - 32768;
  sig = round(sig*32768);  % MATLAB compatible
%  if find(sig > 32767) | find(sig < 32768)
%    fprintf(2, 'Clipping!\n');
%  end
  Block = 2*n;
end
if format(4) == 24
  if fmt
    sig = round(sig*2^23);
  else
    sig = round((sig+1)*(2^24-1)/2) - 2^23;
  end
  Block = 3*n;
end

Riffsize=36+(Block*nsamples);

% write riff chunk
%fwrite(fid,'RIFF','uchar');
fwrite(fid,'RIFF','char');
%fwrite(fid,Riffsize,'ulong');
fwrite(fid,Riffsize,'uint32');
%fwrite(fid,'WAVE','uchar');
fwrite(fid,'WAVE','char');

% write format sub-chunk
%fwrite(fid,'fmt ','uchar');
% fwrite(fid,16,'ulong');
fwrite(fid,'fmt ','char');
fwrite(fid,16,'uint32');

Average = (format(4)/8)*format(3)*format(2);

% fwrite(fid,format(1),'ushort');         % PCM or X-law format 
% fwrite(fid,format(2),'ushort');         % x channel
% fwrite(fid,format(3),'ulong');      	% samples per second
% fwrite(fid,Average,'ulong');      	% average bytes per second
% fwrite(fid,Block,'ushort');        	% block alignment
% fwrite(fid,format(4),'ushort');         % bits per sample

fwrite(fid,format(1),'uint16');         % PCM or X-law format 
fwrite(fid,format(2),'uint16');         % x channel
fwrite(fid,format(3),'uint32');      	% samples per second
fwrite(fid,Average,'uint32');      	% average bytes per second
fwrite(fid,Block,'uint16');        	% block alignment
fwrite(fid,format(4),'uint16');         % bits per sample

% write data sub-chunck
% fwrite(fid,'data','uchar');
% fwrite(fid,(nsamples*Block),'ulong');
fwrite(fid,'data','char');
fwrite(fid,(nsamples*Block),'uint32');

if (format(4)+format(2) == 9)
  fwrite(fid,sig,'uchar');
end

if (format(4)+format(2) == 10)
  if invs == 1
%    fwrite(fid,sig(:,1:2)','uchar');
    fwrite(fid,sig(:,1:2)','char');
  else
%    fwrite(fid,sig(1:2,:),'uchar');
    fwrite(fid,sig(1:2,:),'char');
  end
end

if (format(4)+format(2) == 17)
%  fwrite(fid,sig,'ushort');
%  fwrite(fid,sig,'short');
  fwrite(fid,sig,'int16');
end

if (format(4)+format(2) == 18)
  if invs == 1
%    fwrite(fid,sig(:,1:2)','ushort');
%    fwrite(fid,sig(:,1:2)','short');
    fwrite(fid,sig(:,1:2)','int16');
  else
%    fwrite(fid,sig(1:2,:),'ushort');
%    fwrite(fid,sig(1:2,:),'short');
    fwrite(fid,sig(1:2,:),'int16');
  end
end

if format(4) == 24
  if invs == 1
    sig = sig';
  end
  signm = (sig >= 0);
  sigh = (1-signm).*(255-floor(-(sig+cnst)/2^16)) + signm.*floor(sig/2^16);
  sigm = (1-signm).*(255-floor(mod(-(sig+cnst), 2^16)/2^8)) + signm.*floor(mod(sig, 2^16)/2^8);
  sigl = (1-signm).*(255-mod(-(sig+cnst), 2^8)) + signm.*mod(sig, 2^8);
  fwrite(fid, [reshape(sigl, 1, n*m); reshape(sigm, 1, n*m); reshape(sigh, ...
						  1, n*m)], 'uchar');
end

if ischar(wavefile)
  fclose(fid);
end
