code: def _multitaper_cross_spectrum(self, clm, slm, k, convention='power',
                                     unit='per_l', lmax=None, taper_wt=None):
          
          if lmax is None:
              lmax = min(clm.lmax, slm.lmax)
  
          sh1 = clm.to_array(normalization='4pi', csphase=1, lmax=lmax)
          sh2 = slm.to_array(normalization='4pi', csphase=1, lmax=lmax)
  
          if taper_wt is None:
              mtse, sd = _shtools.SHMultiTaperMaskCSE(sh1, sh2, self.tapers,
                                                      lmax=lmax, k=k)
          else:
              mtse, sd = _shtools.SHMultiTaperMaskCSE(sh1, sh2, self.tapers,
                                                      lmax=lmax, k=k,
                                                      taper_wt=taper_wt)
  
          if (unit == 'per_l'):
              pass
          elif (unit == 'per_lm'):
              degree_l = _np.arange(len(mtse))
              mtse /= (2.0 * degree_l + 1.0)
              sd /= (2.0 * degree_l + 1.0)
          else:
              raise ValueError(
                  "unit must be 'per_l' or 'per_lm'." +
                  "Input value was {:s}".format(repr(unit)))
  
          if (convention == 'power'):
              return mtse, sd
          elif (convention == 'energy'):
              return mtse * 4.0 * _np.pi, sd * 4.0 * _np.pi
          else:
              raise ValueError(
                  "convention must be 'power' or 'energy'." +
                  "Input value was {:s}".format(repr(convention)))
docs: Return the multitaper cross-spectrum estimate and standard error for
          two input SHCoeffs class instances.