def box_actions(results, times, N_matrix, ifprint):
    """
        Finds actions, angles and frequencies for box orbit.
        Takes a series of phase-space points from an orbit integration at times t and returns
        L = (act,ang,n_vec,toy_aa, pars) -- explained in find_actions() below.
    """
    if(ifprint):
        print("\n=====\nUsing triaxial harmonic toy potential")

    t = time.time()
    # Find best toy parameters
    omega = toy.findbestparams_ho(results)
    if(ifprint):
        print("Best omega "+str(omega)+" found in "+str(time.time()-t)+" seconds")

    # Now find toy actions and angles
    AA = np.array([toy.angact_ho(i,omega) for i in results])
    AA = AA[~np.isnan(AA).any(1)]
    if(len(AA)==0):
        return

    t = time.time()
    act = solver.solver(AA, N_matrix)
    if act==None:
        return

    if(ifprint):
        print("Action solution found for N_max = "+str(N_matrix)+", size "+str(len(act[0]))+" symmetric matrix in "+str(time.time()-t)+" seconds")

    np.savetxt("GF.Sn_box",np.vstack((act[1].T,act[0][3:])).T)

    ang = solver.angle_solver(AA,times,N_matrix,np.ones(3))
    if(ifprint):
        print("Angle solution found for N_max = "+str(N_matrix)+", size "+str(len(ang))+" symmetric matrix in "+str(time.time()-t)+" seconds")

    # Just some checks
    if(len(ang)>len(AA)):
        print("More unknowns than equations")

    return act[0], ang, act[1], AA, omega
Finds actions, angles and frequencies for box orbit.
        Takes a series of phase-space points from an orbit integration at times t and returns
        L = (act,ang,n_vec,toy_aa, pars) -- explained in find_actions() below.