code: def vmomentsurfacemass(self,*args,**kwargs):
          
          use_physical= kwargs.pop('use_physical',True)
          ro= kwargs.pop('ro',None)
          if ro is None and hasattr(self,'_roSet') and self._roSet:
              ro= self._ro
          if _APY_LOADED and isinstance(ro,units.Quantity):
              ro= ro.to(units.kpc).value
          vo= kwargs.pop('vo',None)
          if vo is None and hasattr(self,'_voSet') and self._voSet:
              vo= self._vo
          if _APY_LOADED and isinstance(vo,units.Quantity):
              vo= vo.to(units.km/units.s).value
          if use_physical and not vo is None and not ro is None:
              fac= surfdens_in_msolpc2(vo,ro)*vo**(args[1]+args[2])
              if _APY_UNITS:
                  u= units.Msun/units.pc**2*(units.km/units.s)**(args[1]+args[2])
              out= self._vmomentsurfacemass(*args,**kwargs)
              if _APY_UNITS:
                  return units.Quantity(out*fac,unit=u)
              else:
                  return out*fac
          else:
              return self._vmomentsurfacemass(*args,**kwargs)
docs: NAME:
  
             vmomentsurfacemass
             
          PURPOSE:
  
             calculate the an arbitrary moment of the velocity distribution 
             at R times the surfacmass
  
          INPUT:
  
             R - radius at which to calculate the moment (in natural units)
  
             n - vR^n
  
             m - vT^m
  
          OPTIONAL INPUT:
  
             nsigma - number of sigma to integrate the velocities over
  
          KEYWORDS:
  
             romberg - if True, use a romberg integrator (default: False)
  
             deriv= None, 'R', or 'phi': calculates derivative of the moment wrt R or phi
  
          OUTPUT:
  
             <vR^n vT^m  x surface-mass> at R (no support for units)
  
          HISTORY:
  
             2011-03-30 - Written - Bovy (NYU)