Relimit pyplot.psd y-axis plot based on frequency range

So I'm trying to find the psd of various signals sampled at different frequencies and compare them in the 0-50Hz range

I have this code:

for i, th in enumerate(ths): 
    NFFT = 2**(int(np.log2(th.size)) - 8)
    ax.psd(th, Fs = freq[i], NFFT = NFFT, scale_by_freq = False)
ax.set_xlabel("freq (Hz)")
ax.set_xlim(0,50)
#ax.relim()              #These
#ax.autoscale()          #Didn't work

What I am trying to do:

Since the signals are low-pass filtered below 50Hz, any signal with a sampling frequency higher than 100Hz has a part that drops to very low psd. So my ylim are very large, and the part of the psd I want to see is not easily visible. I can manually set_ylim every time I call the function, but I'd like some way to automatically set ylim to fit the psd, however every manual method I've found on SO assumes I know the y-values already. But psd does a lot of calculations on the backend that I don't know how to access.

Is there any way to either:

read the y-values out of a psd plot to an array, or
somehow refit the y-axis based on the current x-axis in pyplot?

Solution

You can use the return values of matplotlib.pyplot.psd. It returns

Pxx (1D array) - The values for the power spectrum P_{xx} before scaling (real valued)
freq (1D array) - The frequencies corresponding to the elements in Pxx.
line (Line2D) - The line created by this function. Only returned if return_line is True.

In this case I believe the frequency values could be used to programmatically set the y-limits appropriately, like this:

mx = []; mn = []
for i, th in enumerate(ths): 
    NFFT = 2**(int(np.log2(th.size)) - 8)
    val, freqs = ax.psd(th, Fs = freq[i], NFFT = NFFT, scale_by_freq = False)
    mx.append(val[freqs <= 50].max())
    mn.append(val[freqs <= 50].min())
ax.set_xlabel("freq (Hz)")
ax.set_xlim(0,50)
ax.set_ylim(np.rint(10 * np.log10(np.array(mn).min())) - 2, 
            np.rint(10 * np.log10(np.array(mx).max())) + 2 )