`astropy.table.join` performance compared to pandas

### Description
Comparing `astropy.table.join` to the same operation done in pandas, the difference is a startling factor of 12 for my real-world test case and a factor of 17 for my constructed minimal example:

```python
import numpy as np
import timeit
from astropy.table import Table, join
import pandas as pd

rng = np.random.default_rng(0)

n_obs = 100
n_events_per_obs = rng.integers(5000, 15000, n_obs)
n_events = n_events_per_obs.sum()

table1 = Table({
"obs_id": np.repeat(np.arange(n_obs), n_events_per_obs),
"event_id": np.concatenate([np.arange(n) for n in n_events_per_obs]),
"value1": np.random.normal(size=n_events)
})

table2 = Table({
"obs_id": table1["obs_id"],
"event_id": table1["event_id"],
"value2": np.random.normal(size=n_events)
})

table1.sort("value1")
table2.sort("value2")

print("tables", timeit.timeit("join(table1, table2, keys=('obs_id', 'event_id'))", number=10, globals=globals()) / 10)

df1 = table1.to_pandas()
df2 = table2.to_pandas()

print("dfs", timeit.timeit("pd.merge(df1, df2, on=('obs_id', 'event_id'))", number=10, globals=globals()) / 10)
```

```
tables 3.3859
dfs 0.2128
```

Start by reproducing the supplied timeit comparison between astropy.table.join and pandas.merge, using the minimal example in the issue. Then trace the astropy.table.join entry point to identify where the runtime is spent; done means demonstrating a meaningful improvement on the reported benchmark.

Written by the indexing model from the issue text.