Using C within Python code for faster fourier transform?

I have a function that computes the fourier transform of an incoming audio stream, and then identifies what musical notes are most dominant within the audio. Everything is working (though I'm sure my code may be janky in many ways). The main thing I'd like to fix is the speed of the function. I'm using a chunk size of 2^14, which makes the function quite accurate, but also way slower than I'd like. I was wondering if there is any way to implement a C fft function that could run faster, but still have it run within my Python code. I'm also open to any other suggestions on how to get this function to run faster.

Here is the current function:

def pitch_calculations(stream, CHUNK, RATE):
    # Read mic stream and then call struct.unpack to convert from binary data back to floats
    data = stream.read(CHUNK, exception_on_overflow=False)
    dataInt = np.array(struct.unpack(str(CHUNK) + 'h', data))
    
    # Apply a window function (Hamming) to the input data
    windowed_data = np.hamming(CHUNK) * dataInt

    # Using numpy fast Fourier transform to convert mic data into frequencies
    fft_result = np.abs(fft.fft(windowed_data, threads=6)) * 2 / (11000 * CHUNK)
    fft_result_copy = fft_result
    freqs = np.fft.fftfreq(len(windowed_data), d=1.0 / RATE)
    
    fft_result = fft_result[:int(len(fft_result)/2)]
    freqs = freqs[:int(len(freqs)/2)]
    
    # Find the indices of local maxima in the frequency spectrum using scipy.signal find_peaks
    localmax_indices = find_peaks(fft_result, height=0.04)[0]
    
    # Get the magnitudes of the local maxima
    strong_freqs = fft_result[localmax_indices]
    
    # Sort the magnitudes in descending order
    sorted_indices = np.argsort(strong_freqs)[::-1]
    
    # Calculate the pitches from the peaks, calculate how many cents sharp/flat
    pitches = []
    
    for index in sorted_indices:
        pitches.append(abs(freqs[localmax_indices[index]]))
        
    note_list = [calculate_tuning(pitch)[0] for pitch in pitches]
    cent_list = [calculate_tuning(pitch)[1] for pitch in pitches]
    note_list = order_pitches(note_list)
    note_list = remove_duplicates(note_list)
    
    
    return note_list, cent_list

As others have said, numpy already performs the underlying FFT computation with optimized native code.

The real problem is that a transform size of 2^14 is a nontrivial cost by any implementation. Suggestions to speed it up:

• Reduce the audio sample rate and/or the transform size, if possible. 2^14 seems very high for determining musical note frequencies. It might be enough to downsample the audio to, say, 24 kHz (e.g. with librosa.resample) and use a transform size of 4096, for a spectrum with frequency resolution of 5.9 Hz per bin. Then use interpolation to estimate peaks with sub-bin percision.
• Use float32 instead of float64: Before the transform, use windowed_data = windowed_data.astype(np.float32) to cast to float32.
• Use real-to-complex FFTs with np.fft.rfft instead of the standard complex-to-complex transform. This might save near a factor 2 in computation time.
• Process a batch of windows at a time, instead of one by one. Particularly, a batch of FFTs is in some cases more efficiently computed than computing FFTs individually. Use for instance np.fft.rfft(x, axis=0) with x being a 2D array to compute the real-to-complex FFT of each column of x.