How to use the galpy.actionAngle.actionAngleStaeckel function in galpy

To help you get started, we’ve selected a few galpy examples, based on popular ways it is used in public projects.

Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.

github GalacticDynamics-Oxford / Agama / pytests / example_galpy.py View on Github external
plt.figure(figsize=(12,8))
    plt.axes([0.04, 0.54, 0.45, 0.45])
    plotCoords(delta, 1.5)
    plt.plot(g_orb[:,0],g_orb[:,3], 'b', label='galpy')  # R,z
    plt.plot(c_orb[:,0],c_orb[:,3], 'g', label='Agama')
    plt.plot(a_orb[:,0],a_orb[:,3], 'r', label='galpy using Agama potential')
    plt.xlabel("R/8kpc")
    plt.ylabel("z/8kpc")
    plt.xlim(0, 1.2)
    plt.ylim(-1,1)
    plt.legend(loc='lower left')

    ### create galpy action/angle finder for the given value of Delta
    ### note: using c=False in the routine below is much slower but apparently more accurate,
    ### comparable to the Agama for the same value of delta
    g_actfinder = galpy.actionAngle.actionAngleStaeckel(pot=g_pot, delta=delta, c=True)

    ### find the actions for each point of the orbit using galpy action finder
    dt = time.time()
    g_act = g_actfinder(g_orb[:,0],g_orb[:,1],g_orb[:,2],g_orb[:,3],g_orb[:,4],fixed_quad=True)
    print 'Time to compute actions in galpy: %.4g s' % (time.time()-dt)
    print 'Jr = %.6g +- %.4g, Jz = %.6g +- %.4g' % \
        (numpy.mean(g_act[0]), numpy.std(g_act[0]), numpy.mean(g_act[2]), numpy.std(g_act[2]))

    ### use the Agama action routine for the same value of Delta as in galpy -
    ### the result is almost identical but computed much faster
    dt = time.time()
    c_act = agama.actions(point=c_orb_car, potential=c_pot, fd=delta)   # explicitly specify focal distance
    print 'Time to compute actions in Agama using Galpy-estimated focal distance: %.4g s' % (time.time()-dt)
    print 'Jr = %.6g +- %.4g, Jz = %.6g +- %.4g' % \
        (numpy.mean(c_act[:,0]), numpy.std(c_act[:,0]), numpy.mean(c_act[:,1]), numpy.std(c_act[:,1]))