
def _twopairs(b, orbit3, solve_for=None, **kwargs):

    hier = b.hierarchy

    orbit3_ps = _get_system_ps(b, orbit3)

    q3 = orbit3_ps.get_parameter(qualifier='q', **_skip_filter_checks)

    m = []
    for star in hier.get_stars():
        star_ps = _get_system_ps(b, star)
        m.append(star_ps.get_parameter(qualifier='mass', **_skip_filter_checks))

    if solve_for in [None]:
        lhs = q3
        rhs = (m[2]+m[3])/(m[0]+m[1])
    else:
        raise NotImplementedError

    return lhs, rhs, [q3] + m, {'orbit3': orbit3}


