docs: def lhoodplot(self, trsig=None, fig=None,
                    piechart=True, figsize=None, logscale=True,
                    constraints='all', suptitle=None, Ltot=None,
                    maxdur=None, maxslope=None, inverse=False,
                    colordict=None, cachefile=None, nbins=20,
                    dur_range=None, slope_range=None, depth_range=None,
                    recalc=False,**kwargs):
          """
          Makes plot of likelihood density function, optionally with transit signal
  
          If ``trsig`` not passed, then just density plot of the likelidhoo
          will be made; if it is passed, then it will be plotted
          over the density plot.
  
          :param trsig: (optional)
              :class:`vespa.TransitSignal` object.
  
          :param fig: (optional)
              Argument for :func:`plotutils.setfig`.
  
          :param piechart: (optional)
              Whether to include a plot of the piechart that describes
              the effect of the constraints on the population.
  
          :param figsize: (optional)
              Passed to :func:`plotutils.setfig`.
  
          :param logscale: (optional)
              If ``True``, then shading will be based on the log-histogram
              (thus showing more detail at low density).  Passed to
              :func:`vespa.stars.StarPopulation.prophist2d`.
  
          :param constraints: (``'all', 'none'`` or ``list``; optional)
              Which constraints to apply in making plot.  Picking
              specific constraints allows you to visualize in more
              detail what the effect of a constraint is.
  
          :param suptitle: (optional)
              Title for the figure.
  
          :param Ltot: (optional)
              Total of ``prior * likelihood`` for all models.  If this is
              passed, then "Probability of scenario" gets a text box
              in the middle.
  
          :param inverse: (optional)
              Intended to allow showing only the instances that are
              ruled out, rather than those that remain.  Not sure if this
              works anymore.
  
          :param colordict: (optional)
              Dictionary to define colors of constraints to be used
              in pie chart.  Intended to unify constraint colors among
              different models.
  
          :param cachefile: (optional)
              Likelihood calculation cache file.
  
          :param nbins: (optional)
              Number of bins with which to make the 2D histogram plot;
              passed to :func:`vespa.stars.StarPopulation.prophist2d`.
  
          :param dur_range, slope_range, depth_range: (optional)
              Define ranges of plots.
  
          :param **kwargs:
              Additional keyword arguments passed to
              :func:`vespa.stars.StarPopulation.prophist2d`.
  
          """
  
          setfig(fig, figsize=figsize)
  
          if trsig is not None:
              dep,ddep = trsig.logdepthfit
              dur,ddur = trsig.durfit
              slope,dslope = trsig.slopefit
  
              ddep = ddep.reshape((2,1))
              ddur = ddur.reshape((2,1))
              dslope = dslope.reshape((2,1))
  
              if dur_range is None:
                  dur_range = (0,dur*2)
              if slope_range is None:
                  slope_range = (2,slope*2)
  
          if constraints == 'all':
              mask = self.distok
          elif constraints == 'none':
              mask = np.ones(len(self.stars)).astype(bool)
          else:
              mask = np.ones(len(self.stars)).astype(bool)
              for c in constraints:
                  if c not in self.distribution_skip:
                      mask &= self.constraints[c].ok
  
          if inverse:
              mask = ~mask
  
          if dur_range is None:
              dur_range = (self.stars[mask]['duration'].min(),
                           self.stars[mask]['duration'].max())
          if slope_range is None:
              slope_range = (2,self.stars[mask]['slope'].max())
          if depth_range is None:
              depth_range = (-5,-0.1)
  
          #This may mess with intended "inverse" behavior, probably?
          mask &= ((self.stars['duration'] > dur_range[0]) &
                   (self.stars['duration'] < dur_range[1]))
          mask &= ((self.stars['duration'] > dur_range[0]) &
                   (self.stars['duration'] < dur_range[1]))
  
          mask &= ((self.stars['slope'] > slope_range[0]) &
                   (self.stars['slope'] < slope_range[1]))
          mask &= ((self.stars['slope'] > slope_range[0]) &
                   (self.stars['slope'] < slope_range[1]))
  
          mask &= ((np.log10(self.depth) > depth_range[0]) &
                   (np.log10(self.depth) < depth_range[1]))
          mask &= ((np.log10(self.depth) > depth_range[0]) &
                   (np.log10(self.depth) < depth_range[1]))
  
  
  
  
          if piechart:
              a_pie = plt.axes([0.07, 0.5, 0.4, 0.5])
              self.constraint_piechart(fig=0, colordict=colordict)
  
          ax1 = plt.subplot(222)
          if not self.is_ruled_out:
              self.prophist2d('duration', 'depth', logy=True, fig=0,
                              mask=mask, interpolation='bicubic',
                              logscale=logscale, nbins=nbins, **kwargs)
          if trsig is not None:
              plt.errorbar(dur,dep,xerr=ddur,yerr=ddep,color='w',marker='x',
                           ms=12,mew=3,lw=3,capsize=3,mec='w')
              plt.errorbar(dur,dep,xerr=ddur,yerr=ddep,color='r',marker='x',
                           ms=10,mew=1.5)
          plt.ylabel(r'log($\delta$)')
          plt.xlabel('')
          plt.xlim(dur_range)
          plt.ylim(depth_range)
          yt = ax1.get_yticks()
          plt.yticks(yt[1:])
          xt = ax1.get_xticks()
          plt.xticks(xt[2:-1:2])
  
          ax3 = plt.subplot(223)
          if not self.is_ruled_out:
              self.prophist2d('depth', 'slope', logx=True, fig=0,
                              mask=mask, interpolation='bicubic',
                              logscale=logscale, nbins=nbins, **kwargs)
          if trsig is not None:
              plt.errorbar(dep,slope,xerr=ddep,yerr=dslope,color='w',marker='x',
                           ms=12,mew=3,lw=3,capsize=3,mec='w')
              plt.errorbar(dep,slope,xerr=ddep,yerr=dslope,color='r',marker='x',
                           ms=10,mew=1.5)
          plt.ylabel(r'$T/\tau$')
          plt.xlabel(r'log($\delta$)')
          plt.ylim(slope_range)
          plt.xlim(depth_range)
          yt = ax3.get_yticks()
          plt.yticks(yt[1:])
  
          ax4 = plt.subplot(224)
          if not self.is_ruled_out:
              self.prophist2d('duration', 'slope', fig=0,
                              mask=mask, interpolation='bicubic',
                              logscale=logscale, nbins=nbins, **kwargs)
          if trsig is not None:
              plt.errorbar(dur,slope,xerr=ddur,yerr=dslope,color='w',marker='x',
                           ms=12,mew=3,lw=3,capsize=3,mec='w')
              plt.errorbar(dur,slope,xerr=ddur,yerr=dslope,color='r',marker='x',
                           ms=10,mew=1.5)
          plt.ylabel('')
          plt.xlabel(r'$T$ [days]')
          plt.ylim(slope_range)
          plt.xlim(dur_range)
          plt.xticks(xt[2:-1:2])
          plt.yticks(ax3.get_yticks())
  
          ticklabels = ax1.get_xticklabels() + ax4.get_yticklabels()
          plt.setp(ticklabels,visible=False)
  
          plt.subplots_adjust(hspace=0.001,wspace=0.001)
  
          if suptitle is None:
              suptitle = self.model
          plt.suptitle(suptitle,fontsize=20)
  
          if Ltot is not None:
              lhood = self.lhood(trsig, recalc=recalc)
              plt.annotate('%s:\nProbability\nof scenario: %.3f' % (trsig.name,
                                                                    self.prior*lhood/Ltot),
                           xy=(0.5,0.5),ha='center',va='center',
                           bbox=dict(boxstyle='round',fc='w'),
                           xycoords='figure fraction',fontsize=15)