Some elementary Matlab programs to test your system.

[j_j]

This creates a bit longer than 10 second clip of weird broadband noise. It's value is that it is energetic, and has no zeros or near zeros in its spectrum.
It creates a file called "whooosh.wave".

You may have to swap to wavwrite on octave. The rest should "just work" in octave, although I lack a modern octave any longer, having actual matlab.

---file pgen.m---
clc
clear all
close all

fs=48000;

blen=2^ceil(log2(10*fs)); %next step up from 10 seconds.

x(1:blen,1)=0;
x(round(blen/2),1)=1;

xt=fft(x);

%phi only needs to be 2:blen/2 memory is cheap
phi=rand(blen,1)*2*pi;% many will be unused, so?
phi=phi-sum(phi)/blen; % make DC value zero
phi=phi *5.6 ; %just fiddle for the right spread


[bb,aa]=butter(4, .1); %the fractional lp frequency affects spread too

phi=filter(bb,aa,phi); %limit speed of phase change.

phi(1)=0; % to prevent complex DC part

for ii=2:blen/2
xt(ii)=xt(ii)*( cos(phi(ii)) + i*sin(phi(ii)));
end
xt(blen:-1: (blen/2+2))=conj(xt(2: (blen/2))); % conjugate spectrum for real only output.
%notice 1 and blen/2+1 are unmodified and purely real, not complex

xx=ifft(xt);
t=max(abs(xx));
xx=xx*.9/t;

audiowrite('whooosh.wav',xx,fs,'BitsPerSample', 24);

plot(xx)
figure
plot(abs(fft(xx,2*blen)))
---end pgen.m---

Free for non-commercial use.

 

[sam_adams]

@j_j, I know this will be good. Thanks in advance.

 

[j_j]

What do you do with that strange sound, you ask?

First, if you want to run at 96kHz, change "fs" in the script. You literally need to do nothing more. Yeah, the signal will be half the length. No pain will ensue.

Set up your DAC to loop first channel back to your ADC. Make sure level is never EVER clipping (note, any clipping with this will cause utter havok, so only do that for your own amusement).
Set the same DAC output (tee it or play the mono file into both channels identically) into the device under test. Likewise check for "no clipping".

Now, start recording. Wait a second or 3. Play the file. After the file ends, let it run for at least another 20 seconds or so. You want a length of at least 3 times the entire probe signal length including both 'before' and 'after'

Put the captured files (channel 1 and 2) somewhere. An elementary program to show the results will appear at some point in the not too distant future.

 

[j_j]

file calcimp.m
clc
clear all
close all

[tprob, fs]=audioread('whooosh.wav');

tsig=audioread('filtered.wav');

blen=2*max(length(tprob), length(tsig));

blen=2^ceil(log2(blen)); % if you mess up the file read this will choke most likely.

tprobt=fft(tprob,blen);
tsigt=fft(tsig,blen);

%deconvolve which you can do no zeros in tprobt unless you messed up

rest=tsigt ./ tprobt; % yeah believe it!

fres=abs(rest(1:blen/2+1));

t=max(fres);
fres=fres/t; % set peak gain to zero db
fres=max(fres,10^-6); % log of zero uncool!
fres=20*log10(fres);

f=(0:blen/2)/blen*fs;

plot(f, fres);

impr=ifft(rest);

figure
plot(impr)

audiowrite('impres.wav',impr,fs,'BitsPerSample',32); % Attention, added 3/26/25 to write out impulse response in 32bit float


***end file

Free for non-commercial use.

 

[j_j]

You run that on your original "whoosh" and on both your loop back (measures your system very (cough) accurately as well) and then your device under test.
The result file is named "filtered.wav"


SAME SAMPLING RATE USING THE SAME CLOCK AS THE DAC CLOCK otherwise I have a good idea of "wtf" happens, and it's most assuredly "wtf" on steroids.

For the file I created (dry lab) I ran a preposterous elliptical filter digitally on the "whooosh" and wrote it out as, of course "filtered". Attached frequency response and impuse response. I had to zoo

m in A LOT for the impulse response, of course. You will have to do that too. Fortunately matlab is good at that. If you're using octave you're on your own.

 

[j_j]

Any time you find a showing up in the code it is supposed to be : ( In matlab is fine. In the BB software it turns into a frown.

 

[JP]

Put the code in a code block to prevent that.

 

[j_j]

Put the code in a code block to prevent that.

Now that exists! I didn't know that. I wonder if it's possible to add text file with .m suffix so I could drop the program as a file?

 

[sam_adams]

Relevant Xenforo BBCODEs for posting code:

 

[DonH56]

Now that exists! I didn't know that. I wonder if it's possible to add text file with .m suffix so I could drop the program as a file?

You have to change the extension then upload. You can take e.g. xxx.m, change it to xxx.txt (or one of the other allowed extensions, see below), and upload. Users can download and change the extension back to "m" to run it. That is what I have done in the past.

I have a zillion or so functions (and expect you have many times that) but mostly Mathcad since years ago Matlab more than doubled the cost of annual maintenance, and before they added the Home option, so I switched to Mathcad (which then did the same thing, more than doubling the cost). I follow the IEEE Standards for choosing frequencies, processing the results, and so forth (I had a hand in some of them: 1057, 1241, and 1658 are the main ones for ADCs and DACs).

 

[j_j]

This plots the Hilbert Envelope of an impulse response (or anything else, really):
***envplot.m***

clc
clear all
close all

[x,fs]=audioread('impres.wav'); %ASSuming mono file here

len=length(x); %assuming input from impulse response here. If not, pre and
%post with zeros to the next power of two length, please. At least 10%
%padding, right? set "len" accordingly too. thmink first!

xt=fft(x);
xt(len/2:len)=0; %

xa=ifft(xt); % analytic signal

xenv=abs(xa);% Hilbert envelope

plot(xenv)

***end envplot.m***

You can also log the envelope for dB information, best normalize the peak to 1 first.

You could add this to the calcimp code above at the end, just copy the data before the impulse response is done, adjust 'len' to blen, and it should all work.

If, say, somebody wanted to write code that took the probe file and 2 input files, one loopback and the other the DUT through the same equipment as the loopback AT THE SAME TIME, one could use the peak in the loopback and the peak in the other data to determine latency, for instance. You're welcome to do that. In fact, anyone who wants to dress up this stuff is welcome to do so. You can use it for anything EXCEPT making the code all or part of a commercial product.

 

[j_j]

To show the generality of the envelope plot, I calculated the envelope of "epitaph" by King Crimson. Posted is the first 32k samples.

Quite by accident, it worked for two channels without being asked to.

 

[DonH56]

To show the generality of the envelope plot, I calculated the envelope of "epitaph" by King Crimson. Posted is the first 32k samples.

View attachment 439344Quite by accident, it worked for two channels without being asked to.

A number of times Matlab's innate vector processing has done unexpectedly great things, like handling multiple audio/RF channels or arrays of S-parameters, but it can be frustratingly annoying (maddening more like) when it does not...

 

[j_j]

A number of times Matlab's innate vector processing has done unexpectedly great things, like handling multiple audio/RF channels or arrays of S-parameters, but it can be frustratingly annoying (maddening more like) when it does not...

Which is why I don't trust it a lot.

 

[JP]

This plots the Hilbert Envelope of an impulse response (or anything else, really):
***envplot.m***

clc
clear all
close all

[x,fs]=audioread('impres.wav'); %ASSuming mono file here

len=length(x); %assuming input from impulse response here. If not, pre and
%post with zeros to the next power of two length, please. At least 10%
%padding, right? set "len" accordingly too. thmink first!

xt=fft(x);
xt(len/2:len)=0; %

xa=ifft(xt); % analytic signal

xenv=abs(xa);% Hilbert envelope

plot(xenv)

***end envplot.m***

You can also log the envelope for dB information, best normalize the peak to 1 first.

You could add this to the calcimp code above at the end, just copy the data before the impulse response is done, adjust 'len' to blen, and it should all work.

If, say, somebody wanted to write code that took the probe file and 2 input files, one loopback and the other the DUT through the same equipment as the loopback AT THE SAME TIME, one could use the peak in the loopback and the peak in the other data to determine latency, for instance. You're welcome to do that. In fact, anyone who wants to dress up this stuff is welcome to do so. You can use it for anything EXCEPT making the code all or part of a commercial product.

@j_j you can just wrap the entire block in the code block. The inline code block is per line which is cumbersome.

clc
clear all
close all

[x,fs]=audioread('impres.wav'); %ASSuming mono file here

len=length(x); %assuming input from impulse response here. If not, pre and
%post with zeros to the next power of two length, please. At least 10%
%padding, right?  set "len" accordingly too. thmink first!
xt=fft(x);
xt(len/2:len)=0; %

xa=ifft(xt); % analytic signal

xenv=abs(xa);% Hilbert envelope
plot(xenv)

 

[j_j]

Yeah, I figured that out, after I did the first 3 copy/pastes, so I selected it all, and turned it on. And so it did line by line.
But it was readable. So I let it go.

 

[j_j]

Now, I'm waiting for someone to do a loopback measurement on their stuff. Especially some of the critics of measurements on this site. I want to see what their "significant frequency response errors" are. One trip through the first two programs (and you only have to run the first one once, and if somebody is bored and wants to figure out how to put up the wave file somehow, it's even easier) and you're done.

 

[bmc0]

Now, I'm waiting for someone to do a loopback measurement on their stuff.

Here's a first-generation Focusrite Scarlett 2i2:
Magnitude and phase:

Impulse:


Using REW's standard measurement sweep instead:



Interestingly, the IRs look a little different.

Edit: The difference in the IR shape is due to REW oversampling by default. See post #26.

 

[bmc0]

audiowrite('whooosh.wav',xx,48000,'BitsPerSample', 24);

You may want to change the "48000" to "fs" here.

FWIW, the code appears to work without modification in Octave as long as the "signal" package is loaded (needed for butter()). It plots the impulse response strangely, but changing plot(impr) to plot(real(impr)) fixes that.

 

[j_j]

You may want to change the "48000" to "fs" here.

Done.

FWIW, the code appears to work without modification in Octave as long as the "signal" package is loaded (needed for butter()). It plots the impulse response strangely, but changing plot(impr) to plot(real(impr)) fixes that.

Probably Octave is using single precision, and you're getting residual imaginary parts due to mantissa length. That even happens sometimes in Matlab, but only at crazy lengths.