We explore the interactions between a quark anti-quark pair in a thermal medium based on lattice QCD ensembles with $N_f = 2+1$ dynamical HISQ flavors. Our dataset spans the phenomenologically relevant temperature range between T=140MeV-2GeV based on lattice sizes
$N_\tau=10,12$ and $16$, with an aspect ratio of $N_\sigma/N_\tau=4$.
The peak position $\Omega$ and the width $\Gamma$ of the spectral function of Wilson-line correlators in Coulomb gauge is computed. We assess the information content in the correlation functions and deploy three complementary strategies to reconstruct spectral information: model fits, Pad\'e approximation and the Bayesian BR method. Limitations of each approach are carefully assessed.