def plot(self, modelx, dlines=False, xmin=None, xmax=None,
             ymin=None, ymax=None, **kwargs):
        """Plot the data and model (requires `omega`).

        This assumes that `data` is 1D and that `mfunc` takes one argument
        that should be treated as the X variable.

        """
        import omega as om

        modelx = np.asarray(modelx)
        if modelx.shape != self.data.shape:
            raise ValueError('modelx and data arrays must have same shape')

        modely = self.mfunc(modelx)
        sigmas = self.invsigma**-1 # TODO: handle invsigma = 0

        vb = om.layout.VBox(2)
        vb.pData = om.quickXYErr(modelx, self.data, sigmas,
                                 'Data', lines=dlines, **kwargs)

        vb[0] = vb.pData
        vb[0].addXY(modelx, modely, 'Model')
        vb[0].setYLabel('Y')
        vb[0].rebound(False, True)
        vb[0].setBounds(xmin, xmax, ymin, ymax)

        vb[1] = vb.pResid = om.RectPlot()
        vb[1].defaultField.xaxis = vb[1].defaultField.xaxis
        vb[1].addXYErr(modelx, self.resids, sigmas, None, lines=False)
        vb[1].setLabels('X', 'Residuals')
        vb[1].rebound(False, True)
        # ignore Y values since residuals are on different scale:
        vb[1].setBounds(xmin, xmax)

        vb.setWeight(0, 3)
        return vb