def cov(y=nil, ddof:1, fweights:nil, aweights:nil)
      if y
        m = NArray.vstack([self,y])
      else
        m = self
      end
      w = nil
      if fweights
        f = fweights
        w = f
      end
      if aweights
        a = aweights
        w = w ? w*a : a
      end
      if w
        w_sum = w.sum(axis:-1, keepdims:true)
        if ddof == 0
          fact = w_sum
        elsif aweights.nil?
          fact = w_sum - ddof
        else
          wa_sum = (w*a).sum(axis:-1, keepdims:true)
          fact = w_sum - ddof * wa_sum / w_sum
        end
        if (fact <= 0).any?
          raise StandardError,"Degrees of freedom <= 0 for slice"
        end
      else
        fact = m.shape[-1] - ddof
      end
      if w
        m -= (m*w).sum(axis:-1, keepdims:true) / w_sum
        mw = m*w
      else
        m -= m.mean(axis:-1, keepdims:true)
        mw = m
      end
      mt = (m.ndim < 2) ? m : m.swapaxes(-2,-1)
      mw.dot(mt.conj) / fact
    end
def cov(y=nil, ddof:1, fweights:nil, aweights:nil)
      if y
        m = NArray.vstack([self,y])
      else
        m = self
      end
      w = nil
      if fweights
        f = fweights
        w = f
      end
      if aweights
        a = aweights
        w = w ? w*a : a
      end
      if w
        w_sum = w.sum(axis:-1, keepdims:true)
        if ddof == 0
          fact = w_sum
        elsif aweights.nil?
          fact = w_sum - ddof
        else
          wa_sum = (w*a).sum(axis:-1, keepdims:true)
          fact = w_sum - ddof * wa_sum / w_sum
        end
        if (fact <= 0).any?
          raise StandardError,"Degrees of freedom <= 0 for slice"
        end
      else
        fact = m.shape[-1] - ddof
      end
      if w
        m -= (m*w).sum(axis:-1, keepdims:true) / w_sum
        mw = m*w
      else
        m -= m.mean(axis:-1, keepdims:true)
        mw = m
      end
      mt = (m.ndim < 2) ? m : m.swapaxes(-2,-1)
      mw.dot(mt.conj) / fact
    end
under construction