code: def Er(self,R,z,vR,vz,E,Lz,sinh2u0,u0):
                                     
          u,v= bovy_coords.Rz_to_uv(R,z,self._delta)
          pu= (vR*numpy.cosh(u)*numpy.sin(v)
               +vz*numpy.sinh(u)*numpy.cos(v)) #no delta, bc we will divide it out
          out= (pu**2./2.+Lz**2./2./self._delta**2.*(1./numpy.sinh(u)**2.-1./sinh2u0)
                -E*(numpy.sinh(u)**2.-sinh2u0)
                +(numpy.sinh(u)**2.+1.)*actionAngleStaeckel.potentialStaeckel(u,numpy.pi/2.,self._pot,self._delta)
                -(sinh2u0+1.)*actionAngleStaeckel.potentialStaeckel(u0,numpy.pi/2.,self._pot,self._delta))
  #              +(numpy.sinh(u)**2.+numpy.sin(v)**2.)*actionAngleStaeckel.potentialStaeckel(u,v,self._pot,self._delta)
  #              -(sinh2u0+numpy.sin(v)**2.)*actionAngleStaeckel.potentialStaeckel(u0,v,self._pot,self._delta))
          return out
docs: NAME:
             Er
          PURPOSE:
             calculate the 'radial energy'
          INPUT:
             R, z, vR, vz - coordinates
             E - energy
             Lz - angular momentum
             sinh2u0, u0 - sinh^2 and u0
          OUTPUT:
             Er
          HISTORY:
             2012-11-29 - Written - Bovy (IAS)