How to interpret the results of a 2 sample KS-test

Question

I have some data which I want to analyze by fitting a function to it. To do that, I have two functions, one being a gaussian, and one the sum of two gaussians. To test the goodness of these fits, I test the with scipy's ks-2samp test. The result of both tests are that the KS-statistic is $0.15$, and the P-value is $0.476635$.

There is clearly visible that the fit with two gaussians is better (as it should be), but this doesn't reflect in the KS-test.

I am sure I dont output the same value twice, as the included code outputs the following: (hist_cm is the cumulative list of the histogram points, plotted in the upper frames)

KStest= stats.ks_2samp(hist_cm,hist_fit_cm)
KStest_1gauss = stats.ks_2samp(hist_cm,hist_fit_1gauss_cm)
print KStest
print KStest_1gauss
#output
#(0.15000000000000002, 0.47663525071642981)
#(0.15000000000000002, 0.47663525071642981)

I know the tested list are not the same, as you can clearly see they are not the same in the lower frames.

So i've got two question: Why is the P-value and KS-statistic the same? (this might be a programming question). And how to interpret these values?

edit: I calculate radial velocities from a model of N-bodies, and should be normally distributed. In the figure I showed I've got 1043 entries, roughly between $-300$ and $300$. I then make a (normalized) histogram of these values, with a bin-width of 10. To this histogram I make my two fits (and eventually plot them, but that would be too much code)

bins = np.arange(-300,301,10)
hist, bin_edges = np.histogram(v[i], normed=True, bins=bins)
hist_cm=np.cumsum(hist)
bin_centres = (bin_edges[:-1] + bin_edges[1:])/2
coeff, pcov2 = leastsq(residuals, x0=(0.01,0.,60.,0.01,150.,40.) ,args=(bin_centres, hist))
coeff, pcov2 = leastsq(residuals, x0=(0.01,0.,60.,0.01,150.,40.) ,args=(bin_centres, hist))
hist_fit = gauss2(bin_centres, *coeff)
hist_fit_cm=np.cumsum(hist_fit)
coeff_1gaus, pcov2_1gauss = leastsq(residualsOneGauss, x0=(0.01,0,100), args=(bin_centres, hist))
hist_fit_1gauss = gauss(bin_centres, *coeff_1gaus)
hist_fit_1gauss_cm = np.cumsum(hist_fit_1gauss)

#and the used functions:
def gauss(x, A, mu, sigma):
  return A*np.exp(-(x-mu)**2/(2.*sigma**2))

def gauss2(x,A, mu, sigma, A2, mu2, sigma2):
  if A2<0:
    return 1000
  return A*np.exp(-(x-mu)**2/(2.*sigma**2))+ A2*np.exp(-(x-mu2)**2/(2.*sigma2**2))

def residuals(p, x,y):
  integral = quad( gauss2, -500, 500, args= (p[0],p[1],p[2],p[3],p[4],p[5]))[0]
  penalization = abs(1-integral)*1000
  return y - gauss2(x, p[0],p[1],p[2],p[3],p[4],p[5] ) - penalization

def residualsOneGauss(p,x,y):
  integral = quad( gauss, -500, 500, args= (p[0],p[1],p[2]))[0]
  penalization = abs(1-integral)*1000
  return y - gauss(x, p[0],p[1],p[2]) - penalization