def JR(self,**kwargs):
        """
        NAME:
           JR
        PURPOSE:
           Calculate the radial action
        INPUT:
           fixed_quad= (False) if True, use n=10 fixed_quad
           +scipy.integrate.quad keywords
        OUTPUT:
           J_R(R,vT,vT)/ro/vc + estimate of the error (nan for fixed_quad)
        HISTORY:
           2012-11-27 - Written - Bovy (IAS)
        """
        if hasattr(self,'_JR'): #pragma: no cover
            return self._JR
        umin, umax= self.calcUminUmax()
        #print self._ux, self._pux, (umax-umin)/umax
        if (umax-umin)/umax < 10.**-6: return nu.array([0.])
        order= kwargs.pop('order',10)
        if kwargs.pop('fixed_quad',False):
            # factor in next line bc integrand=/2delta^2
            self._JR= 1./nu.pi*nu.sqrt(2.)*self._delta\
                *integrate.fixed_quad(_JRStaeckelIntegrand,
                                      umin,umax,
                                      args=(self._E,self._Lz,self._I3U,
                                            self._delta,
                                            self._u0,self._sinhu0**2.,
                                            self._vx,self._sinvx**2.,
                                            self._potu0v0,self._pot),
                                      n=order,
                                      **kwargs)[0]
        else:
            self._JR= 1./nu.pi*nu.sqrt(2.)*self._delta\
                *integrate.quad(_JRStaeckelIntegrand,
                                umin,umax,
                                args=(self._E,self._Lz,self._I3U,
                                      self._delta,
                                      self._u0,self._sinhu0**2.,
                                      self._vx,self._sinvx**2.,
                                      self._potu0v0,self._pot),
                                **kwargs)[0]
        return self._JR
NAME:
           JR
        PURPOSE:
           Calculate the radial action
        INPUT:
           fixed_quad= (False) if True, use n=10 fixed_quad
           +scipy.integrate.quad keywords
        OUTPUT:
           J_R(R,vT,vT)/ro/vc + estimate of the error (nan for fixed_quad)
        HISTORY:
           2012-11-27 - Written - Bovy (IAS)