I'm trying to perform a relatively simple task - I have a digital filter implemented on an embedded system, and I'd like to generate the data to describe its magnitude and phase response vs. frequency.
I have been able to successfully generate the magnitude data using a number of methods, but am having trouble with the phase.
In octave, I've generated the magnitude and frequency response plots using the function freqz as follows:
freqz_plot(W_hz, H, 0)
In this case, W_hz is a vector of frequency data in Hz, and H is the computed response of the filter (i.e., the complex value evaluated at each frequency).
I get the following plots, which look like what I expect for the chosen filter design:

However, if I manually compute the phase data using the following formula:
phase_deg = rad2deg(angle(H));
I get this plot instead, full of large discontinuities (forgive the lack of axis labels, I was rushing, its phase in deg vs. frequency in hz):

I can't figure out the math regarding how to correctly eliminate the discontinuities and reproduce the plots that freqz generates. I can see that some of the discontinuities are 180deg and some are 90deg, but I can't figure how to correctly determine the size/direction of the discontinuities and account for it properly in the calculation so that I can generate the data myself.
Thanks!
Managed to figure it out - the answer lies in the explanation of the unwrap function, which successfully generates the expected continuous phase data:
Q = unwrap(P) unwraps the radian phase angles in a vector P. Whenever the jump between consecutive angles is greater than or equal to π radians, unwrap shifts the angles by adding multiples of ±2π until the jump is less than π
In usage, for example:
phase_deg_2 = rad2deg(unwrap(angle(H)));
where H is the computed frequency response of the filter as described above