Sound by Matlab

cosine ramp


% prepare tone
cf = 2000;                  % carrier frequency (Hz)
sf = 22050;                 % sample frequency (Hz)
d = 1.0;                    % duration (s)
n = sf * d;                 % number of samples
s = (1:n) / sf;             % sound data preparation
s = sin(2 * pi * cf * s);   % sinusoidal modulation

% prepare ramp
dr = d / 10;
nr = floor(sf * dr);
r = sin(linspace(0, pi/2, nr));
r = [r, ones(1, n - nr * 2), fliplr(r)];

% make ramped sound
s = s .* r;
sound(s, sf);               % sound presentation
pause(d + 0.5);             % waiting for sound end
plot(s);

amplitude modulation

(1 + m sin 2πgt) sin 2πft
m: modulation index, gt: modulation frequency, ft: carrier frequency


sf = 22050;                        % sample frequency (Hz)
d = 1.0;                           % durati0n (s)
n = sf * d;                        % number of samples

% set carrier
cf = 1000;                         % carrier frequency (Hz)
c = (1:n) / sf;                    % carrier data preparation
c = sin(2 * pi * cf * c);          % sinusoidal modulation

% set modulator
mf = 5;                            % modulator frequency (Hz)
mi = 0.5;                          % modulator index
m = (1:n) / sf;                    % modulator data preparation
m = 1 + mi * sin(2 * pi * mf * m); % sinusoidal modulation

% amplitude modulation
s = m .* c;                        % amplitude modulation

% sound presentation
sound(s, sf);                      % sound presentation
pause(d + 0.5);                    % waiting for sound end
plot(s);

frequency modulation

sin(2πft - β cos 2πgt)
β modulation index, gt: modulation frequency, ft: carrier frequency


sf = 22050;                    % sample frequency (Hz)
d = 1.0;                       % duration (s)
n = sf * d;                    % number of samples

% set carrier
cf = 1000;                     % carrier frequency (Hz)
c = (1:n) / sf;                % carrier data preparation
c = 2 * pi * cf * c;

% set modulator
mf = 5;                        % modulator frequency (Hz)
mi = 0.5;                      % modulator index
m = (1:n) / sf;                % modulator data preparation
m = mi * cos(2 * pi * mf * m); % sinusoidal modulation

% frequency modulation
s = sin(c + m);                % frequency modulation

% sound presentation
sound(s, sf);                  % sound presentation
pause(d + 0.5);                % wating for sound end

filtered noise


% set general variables
sf = 22050;  % sample frequency
nf = sf / 2; % nyquist frequency
d = 1.0;     % duration
n = sf * d;  % number of samples
nh = n / 2;  % half number of samples

% =========================================================================
% select filter
while 1
    clc;
    v = input('lowpass(1), highpass(2), bandpass(3), notch(4) ? ');
    if v >= 1 & v <= 4;
        break
    end
end

% =========================================================================
% set variables for filter
lf = 2000;   % lowest frequency
hf = 4000;   % highest frequency
lp = lf * d; % ls point in frequency domain    
hp = hf * d; % hf point in frequency domain

% design filter
clc;
switch v
    case 1
    a = ['LOWPASS'];
    filter = zeros(1, n);           % initializaiton by 0
    filter(1, 1 : lp) = 1;          % filter design in real number
    filter(1, n - lp : n) = 1;      % filter design in imaginary number
    case 2        
    a = ['HIGHPASS'];
    filter = ones(1, n);            % initializaiton by 1
    filter(1, 1 : hp) = 0;          % filter design in real number
    filter(1, n - hp : n) = 0;      % filter design in imaginary number
    case 3
    a = ['BANDPASS'];
    filter = zeros(1, n);           % initializaiton by 0
    filter(1, lp : hp) = 1;         % filter design in real number
    filter(1, n - hp : n - lp) = 1; % filter design in imaginary number
    case 4
    a = ['NOTCH'];
    filter = ones(1, n);
    filter(1, lp : hp) = 0;
    filter(1, n - hp : n - lp) = 0;
end

% =========================================================================
% make noise
rand('state',sum(100 * clock));  % initialize random seed
noise = randn(1, n);             % Gausian noise
noise = noise / max(abs(noise)); % -1 to 1 normalization

% do filter
s = fft(noise);                  % FFT
s = s .* filter;                 % filtering
s = ifft(s);                     % inverse FFT
s = real(s);
% =========================================================================
% play noise
disp('WHITE noise');
sound(noise, sf);                % playing sound
pause(d + 0.5);                  % waiting for sound end

% play filtered noise
clc;
disp([a, ' noise']);
sound(s, sf);              % playing sound
pause(d + 0.5);                  % waiting for sound end

% =========================================================================
% plot sound
x = linspace(0, d, n);
subplot(2,2,1); plot(x, noise); xlabel('time (s)'); title('sound: noise');               
subplot(2,2,2); plot(x, s); xlabel('time (s)'); title('sound: filtered noise');

% plot Fourier spectrum
x = linspace(0, nf, nh);
t = fft(noise);
t = t .* conj(t);
subplot(2,2,3); semilogy(x, t(1,1:nh) ./ max(t));  xlabel('frequency (Hz)'); title('spectrum: noise'); 
t = fft(s);
t = t .* conj(t);
subplot(2,2,4); semilogy(x, t(1,1:nh) ./ max(t));  xlabel('frequency (Hz)');  title('spectrum: filtered noise'); 

figure(1);

write as wav file

useful when not enough memory for 'wavwrite'


% prepare sample sound
sf = 22050;                   % sample frequency (Hz)
n = floor(sf * 1.0);
t = (1:n) / sf;               % sound data preparation
sL = sin(2 * pi *  500 * t);  % sound left: pure tone
sR = rand(1, n) * 2 - 1;      % sound right: noise
s = [sL; sR];                 % sound stereo

% save as wav file
bit = 16;                                  % bit rate (8 or 16)
byte = bit / 8;
[ch, dl] = size(s);
% header
fid = fopen('sample.wav', 'wb');
fwrite(fid, hex2dec('52'), 'uchar');       % 'R' in ASCII code
fwrite(fid, hex2dec('49'), 'uchar');       % 'I' in ASCII code
fwrite(fid, hex2dec('46'), 'uchar');       % 'F' in ASCII code
fwrite(fid, hex2dec('46'), 'uchar');       % 'F' in ASCII code
fwrite(fid, dl * ch * byte + 36, 'ulong'); % file size from next to EOF
fwrite(fid, hex2dec('57'), 'uchar');       % 'W' in ASCII code
fwrite(fid, hex2dec('41'), 'uchar');       % 'A' in ASCII code
fwrite(fid, hex2dec('56'), 'uchar');       % 'V' in ASCII code
fwrite(fid, hex2dec('45'), 'uchar');       % 'E' in ASCII code
fwrite(fid, hex2dec('66'), 'uchar');       % 'f' in ASCII code
fwrite(fid, hex2dec('6D'), 'uchar');       % 'm' in ASCII code
fwrite(fid, hex2dec('74'), 'uchar');       % 't' in ASCII code
fwrite(fid, hex2dec('20'), 'uchar');       % ' ' in ASCII code
fwrite(fid, bit, 'ulong');                 % quantization bit rate
fwrite(fid, 1, 'ushort');
fwrite(fid, ch, 'ushort');                 % # of channels
fwrite(fid, sf, 'ulong');                  % sample frequency
fwrite(fid, sf * ch * byte, 'ulong');      % data speed
fwrite(fid, ch * byte, 'ushort');          % block size
fwrite(fid, bit, 'ushort');                % bit rate
fwrite(fid, hex2dec('64'), 'uchar');       % 'd' in ASCII code
fwrite(fid, hex2dec('61'), 'uchar');       % 'a' in ASCII code
fwrite(fid, hex2dec('74'), 'uchar');       % 't' in ASCII code
fwrite(fid, hex2dec('61'), 'uchar');       % 'a' in ASCII code
fwrite(fid, dl * ch * byte, 'ulong');      % sound data size (byte)
% data
s=s-min(s(:));                             % minimum value = 0
s=s/max(s(:))*(2^bit)-2^(bit-1);           % -2^(bit-1) - 2^(bit-1)
s=int16(s);
for y = 1:dl
    for x = 1:ch
        fwrite(fid, s(x, y), 'short');     % sound data 
    end
end
fclose(fid);

How to make and play sound with Matlab

pure tone

cosine ramp

amplitude modulation

frequency modulation

filtered noise

write as wav file