I am using sf.plot_quad for the same circuit. However, I get different y-axis values depending on the range of x-axis values over which I produce the plots. Below is a short code snippet that should hopefully allow you to easily reproduce this on your end. The circuit is simple enough (3 Gaussian gates on a single qumode) that the correct output state can be computed by hand. In case this is useful, I have done this calculation and I can confirm that the correct output state is being computed. In other words, results.state.means() is returning the correct values.

```
#
# code snippet begins
#
import strawberryfields as sf
from strawberryfields import ops
import numpy as np
# set hbar
sf.hbar = 1
# initialize simulation engine
eng = sf.Engine('gaussian')
# initialize qumode(s)
prog = sf.Program(1)
with prog.context as q:
# create circuit
ops.Vac | (q[0])
ops.Pgate(-6.419) | (q[0])
ops.Zgate(12.84) | (q[0])
ops.Rgate(-0.0866) | (q[0])
results = eng.run(prog)
### print <x> and <p>
means = results.state.means()
print("<x>,<p>={0}".format(means))
# plot 1
# to the naked eye, this plot appears correct
# <x> and <p> both appear correct here
xvec = np.arange(0, 25, 0.01)
pvec = np.arange(0, 25, 0.1)
sf.plot_quad(state=results.state, modes=[0],
xvec=xvec, pvec=pvec, renderer="notebook")
# plot 2
# this plot appears incorrect
# <x> appears incorrect here
# also curves appear not to integrate to 1
xvec = np.arange(1, 2.5, 0.01)
pvec = np.arange(1, 2.5, 0.01)
sf.plot_quad(state=results.state, modes=[0],
xvec=xvec, pvec=pvec, renderer="notebook")
```