code: def plot(self, modelx, dlines=False, xmin=None, xmax=None,
               ymin=None, ymax=None, **kwargs):
          
          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
docs: 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.