English: made with Mathematica using the following commands:

probX[x_, t_, m_, h_, om1_, om2_, x0_] =

Sqrt[m om1 om2^2 /(Pi h (om1^2 Sin[om2 t]^2 + om2^2 Cos[om2 t]^2))]*
 Exp[-(m om1 om2^2 / h)*(x - x0*Cos[om2 t])^2/(om1^2 Sin[om2 t]^2 + 
      om2^2 Cos[om2 t]^2)]

probP[p_, t_, m_, h_, om1_, om2_, x0_] =

Sqrt[om1/(Pi m h (om1^2 Cos[om2 t]^2 + om2^2 Sin[om2 t]^2))]*
 Exp[-(om1/ (m h))*(p + m om1 x0*Sin[om2 t])^2/(om1^2 Cos[om2 t]^2 + 
      om2^2 Sin[om2 t]^2)]

v[x_, m_, om2_] = m om2^2 x^2/2 frame2[t_] :=

Plot[{probX[x, t, 1, 1, 2, 1, -10], probP[x, t, 1, 1, 2, 1, -10], 
  v[x, 1/200, 1]}, {x, -25, 25}, Axes -> {True, False}, 
 Ticks -> None, PlotRange -> {0, 1.2}]

frames2 = ParallelTable[frame2[t], {t, 0, 2*Pi, 0.1}];