How to generate uniform distributed samples with given auto-correlation function As I mentioned in the question title, I want to generate specific uniformly distributed samples.
I need them to model a real world scenario. For my real data, I estimated a function, which approximates the auto-correlation function by a e-function. I also modeled a distribution model for my real data.
I can generate samples, which follow my distribution model with the inverse transform sampling.
If i think correct, instead of "standard" uniform distributed samples, which i feed into the inverse transform sampling, i need to adapt those uniform distributed samples, so that their auto-correlation function follows my estimated e-function.
Now the question is, if i am correct, how to adapt my uniform distributed samples, that they correspond to my requirements.
It would be great if, somebody could help.
EDIT 1:
Here is a plot of the auto-correlation function of my real data (different scenarios) and the modeled function (green). The ACF of my uniform distributed samples should follow the green (respectively blue) line.

EDIT 2:
Thank you so much so far for the answers. To make my problem a bit clearer, i added a illustration of the problem. For a simulation, i need to generate samples, which follow my distribution model (with inverse transform sampling a easy task) AND have the same auto-correlation behavior like the real data.
My idea was, that i need to adapt the uniformly distributed samples, which i feed into the inverse transform sampling, so that their ACF follows my model.
The problem is, that i don't know if this idea is correct and that i don't know how to adapt the uniformly distributed samples.


The following spells out the details of the approach proposed in the other answer by @AdamO and in its comments by @LucaCiti.
For $i=1,2,\dots,\infty$, let $|\phi_i|$ denote the probability that $U_t$ takes a value identical to either $U_{t-i}$ or $1 - U_{t-i}$ and let these two possibilities be determined by the sign of $\phi_i$.  Let the remaining fraction
$$
\phi_0=1-\sum_{i=1}^\infty |\phi_i|
$$
denote the probability that $U_t$ takes a uniformly distributed value independent of the history of the process.  Clearly, we must have
$$
0\le \phi_0\le 1. \tag{1}
$$
and
$$
-1\le \phi_i\le 1 \tag{2}
$$
for $i=1,2,\dots,\infty$.
Letting $V_t=U_t-\frac12$ denote the mean-centered process, and using the law of total expectation, we have
\begin{align}
E(V_t|V_{t-1},V_{t-2},\dots)
&=|\phi_1|\operatorname{sgn}(\phi_1) V_{t-1} + |\phi_2|\operatorname{sgn}(\phi_2) V_{t-2} + \dots
\\&=\phi_1 V_{t-1} + \phi_2 V_{t-2} + \dots.
\end{align}
Thus it is immedeately clear that $\phi_1,\phi_2,\dots$ are the coefficients in the  $\operatorname{AR}(\infty)$ representation of the model.
Unlike ordinary ARMA models, constraint (1) and (2) implies that not all positive semi-definite autocovariance functions are possible via this construction, however. For example, if the target autocovariance function is that of a MA(1) model with MA polynomial $1-\theta B$, the infinite AR polynomial would equal
$$
\frac1{1-\theta B}=1+\theta B+\theta^2 B^2+\dots,
$$
and we would have
$$
\sum_{i=1}^\infty|\phi_i|=\sum_{i=1}^\infty |\theta^i|=\sum_{i=1}^\infty |\theta|^i = \frac{|\theta|}{1-|\theta|}.
$$
Combined with (1) this limits possible values of $\theta$ to
$$
|\theta|\le \frac12
$$
and the correlation at lag 1 to
$$
-\frac25\le \rho_1=\frac{\theta}{1+\theta^2}\le\frac25
$$
In contrast, via the copula described in my other answer the correlation at lag 1 is limited to $|\rho_1|<0.4825837$ only .  Semipositive definiteness in itself limits the same correlation to $|\rho_1|\le 1/2$.
