Find best fit common probability density function, if already have a probability density function (Python's Fitter does not suffice)

Question

There are questions about finding best-fit p.d.f.s given some data: How to find probability distribution and parameters for real data? (Python 3), https://medium.com/the-researchers-guide/finding-the-best-distribution-that-fits-your-data-using-pythons-fitter-library-319a5a0972e9.

But I already have a p.d.f.

for i in range(N-100,len(normalized_p),10):
       plt.plot(x[i], normalized_p[i], label=f'n={i+1}', marker='o')

and I want to find the best fit model.

To use the Fitter library mentioned above, I ended up just sampling data from my p.d.f.

data = np.random.choice(x[N-1], 100000, p=p[N-1]/sump[N-1])

f = Fitter(data,
           distributions=['gamma',
                          'lognorm',
                          "beta",
                          "burr"])
f.fit()
f.summary()

This barely works; because of the fact that I only sample discrete numbers, things can be slightly "out of alignment" with the bins, leading to weird spikes One solution is to add a little noise to the data so this minute misalignments don't impact the histogram as much.

Is there any better way for me to do this task?