Interpolation to find harmonics from fft

Hey guys, I have graphed the fft of a file and have found the fundamental frequencies. Now I need to find the harmonics of 3 notes whos fundamental frequencies are around 196 Hz, 246.9 Hz, and 329.6 Hz. Below is the graph with the 3 fundamentals and now to find the harmonics I am find the nearest returned peak frequency to each of those locations by using the interp1 function with nearest. However, when I plot the function the values arent at the peaks I want them to be. If someone could help me out that would be greatly appreciated.

This graph shows the 3 fundamental frequencies in the fft.

This graph illustrates my problem where I want the plotted value to be at the peak of X coordinate 489.9 instead of 488.6. The same thing applies to the other points.

%read the audio file x samples, output Fs
[x,Fs]= audioread('Guitar_notes.wav');
%take the fourier transform
xdft = fft(x);
xdft = xdft(1:length(x)/2+1);
% create a frequency vector
freq = 0:Fs/length(x):Fs/2;
%locs is now in frequency 
[pks,locs] = findpeaks(20*log10(abs(xdft)), "MinPeakProminence", 30);
%finds the peaks within the range of the fundamental frequency of the notes
indesiredrange = locs > 150 & locs < 500;
%gets the subsets within range
pks_subset = pks(indesiredrange);
locs_subset = locs(indesiredrange);
figure(1)
semilogx(freq,20*log10(abs(xdft)))
hold on
plot(freq(locs_subset), pks_subset, 'or')
hold off
xlabel('f(Hz)');
title('FFT of signal')

harmonic = [1:5];
%gets the locations in frequency
[peaks,freqs] = findpeaks(20*log10(abs(xdft)),freq);
harmB = 244.3*harmonic;
%finds the frequency peaks nearest to the calculated harmonics of note B(244.3hz)
vqB = interp1(freqs,peaks,harmB, 'nearest');
figure(3)
semilogx(freq,20*log10(abs(xdft)))
hold on
plot(harmB, vqB, 'or')
hold off
xlabel('f(Hz)');
title('FFT of signal');