Making quantitative the assumption "the intersection of $L_n$ and $L_m$ cannot be the point $(x_k, f(x_k))$", it is possible to construct such a function $g$ (as pointed out by Jaume, the nonquantitative assumption is not sufficient).
Let us consider the problem in $\mathbb R^n$.

Given a family of indices $I$, let $(x_i)_{i\in I}\subset \mathbb R^n$ be a family of points and let $(L_i)_{i\in I}$ be a family of affine functions. We assume that there is a convex function $f:\mathbb R^n\to\mathbb R$ such that $f(x_i) = L_i(x_i)$ and $f\ge L_i$.

*Quantitative assumption:* There is an $\varepsilon > 0$ such that for any $i\not=j$ we have $L_i(x)\ge L_j(x)$ for any $x\in B(x_i, \varepsilon)$ (when $x = x_i$ this follows from the convexity of $f$) (this assumption is equivalent to the original one if $I$ is finite).

Let $h=\sup_{i\in I} L_i$. The function $h$ is convex and $h\le f < \infty$.
Let $\rho\in C^{\infty}_c(B(0,\varepsilon))$ be a convolution kernel ($\rho\ge 0$ and $\int\rho=1$).
Define $g=h\star\rho$. It is not hard to check that $g$ is smooth, $g(x_i)=L_i(x_i)$ and $g\ge L_i$ for each $i\in I$.

*Relaxing the assumption:* It shall be possible, by appropriately modifying the argument, to prove the result even if the value of $\varepsilon$ is allowed to depend on $i$ (so we would have $\varepsilon_i$), provided that it remains *locally* bounded away from $0$ (i.e., $\inf_{|x_i|<R}\varepsilon_i > 0$).