This is a harder version of the problem. In this version, n ≤ 50 000.

There are n distinct points in three-dimensional space numbered from 1 to n. The i-th point has coordinates (x_i, y_i, z_i). The number of points n is even.

You'd like to remove all n points using a sequence of n/2 snaps. In one snap, you can remove any two points a and b that have not been removed yet and form a perfectly balanced pair. A pair of points a and b is perfectly balanced if no other point c (that has not been removed yet) lies within the axis-aligned minimum bounding box of points a and b.

Formally, point c lies within the axis-aligned minimum bounding box of points a and b if and only if min(x_a, x_b) ≤ x_c ≤ max(x_a, x_b), min(y_a, y_b) ≤ y_c ≤ max(y_a, y_b), and min(z_a, z_b) ≤ z_c ≤ max(z_a, z_b). Note that the bounding box might be degenerate. 

Find a way to remove all points in n/2 snaps.

Input

The first line contains a single integer n (2 ≤ n ≤ 50 000; n is even), denoting the number of points.

Each of the next n lines contains three integers x_i, y_i, z_i (-10^8 ≤ x_i, y_i, z_i ≤ 10^8), denoting the coordinates of the i-th point.

No two points coincide.

Output

Output n/2 pairs of integers a_i, b_i (1 ≤ a_i, b_i ≤ n), denoting the indices of points removed on snap i. Every integer between 1 and n, inclusive, must appear in your output exactly once.

We can show that it is always possible to remove all points. If there are many solutions, output any of them.

Examples

Input


6
3 1 0
0 3 0
2 2 0
1 0 0
1 3 0
0 1 0


Output


3 6
5 1
2 4


Input


8
0 1 1
1 0 1
1 1 0
1 1 1
2 2 2
3 2 2
2 3 2
2 2 3


Output


4 5
1 6
2 7
3 8

Note

In the first example, here is what points and their corresponding bounding boxes look like (drawn in two dimensions for simplicity, as all points lie on z = 0 plane). Note that order of removing matters: for example, points 5 and 1 don't form a perfectly balanced pair initially, but they do after point 3 is removed. 

<image>
#include <bits/stdc++.h>
using namespace std;
void out(int a, bool ln) { printf("%d%c", a, ln ? '\n' : ' '); }
void out(long long a, bool ln) { printf("%lld%c", a, ln ? '\n' : ' '); }
void out(double a, int digit, bool ln) {
  printf("%.*f%c", digit, a, ln ? '\n' : ' ');
}
void out(long double a, int digit, bool ln) {
  printf("%.*Lf%c", digit, a, ln ? '\n' : ' ');
}
void out0(int a[], int n) {
  for (int i = 0; i < n; i++) out(a[i], i == n - 1);
}
void out1(int a[], int n) {
  for (int i = 1; i <= n; i++) out(a[i], i == n);
}
void out0(long long a[], int n) {
  for (int i = 0; i < n; i++) out(a[i], i == n - 1);
}
void out1(long long a[], int n) {
  for (int i = 1; i <= n; i++) out(a[i], i == n);
}
int in(int &a, int &b, int &c, int &d) {
  return scanf("%d%d%d%d", &a, &b, &c, &d);
}
int in(int &a, int &b, int &c) { return scanf("%d%d%d", &a, &b, &c); }
int in(int &a, int &b) { return scanf("%d%d", &a, &b); }
int in(long long &a, long long &b, long long &c, long long &d) {
  return scanf("%lld%lld%lld%lld", &a, &b, &c, &d);
}
int in(long long &a, long long &b, long long &c) {
  return scanf("%lld%lld%lld", &a, &b, &c);
}
int in(long long &a, long long &b) { return scanf("%lld%lld", &a, &b); }
int in(long long &a) { return scanf("%lld", &a); }
int in(int &a) { return scanf("%d", &a); }
int in(char *s) { return scanf("%s", s); }
int in(char &c) { return scanf("%c", &c); }
int in(double &a) { return scanf("%lf", &a); }
int in(long double &a) { return scanf("%Lf", &a); }
void in0(int a[], int n) {
  for (int i = 0; i < n; i++) in(a[i]);
}
void in1(int a[], int n) {
  for (int i = 1; i <= n; i++) in(a[i]);
}
void in0(long long a[], int n) {
  for (int i = 0; i < n; i++) in(a[i]);
}
void in1(long long a[], int n) {
  for (int i = 1; i <= n; i++) in(a[i]);
}
const double pi = acos(-1);
const double eps = 1e-8;
int sign(double a) { return a < -eps ? -1 : a > eps; }
int db_cmp(double a, double b) { return sign(a - b); }
int inmid(double k1, double k2, double k3) {
  return sign(k1 - k3) * sign(k2 - k3) <= 0;
}
int inf = 0x3f3f3f3f;
long long llinf = 0x3f3f3f3f3f3f3f3f;
long long mod = 1e9 + 7;
const int M = 2.1e5;
const int N = 5.1e4;
struct da {
  int x, y, z, id;
  bool operator<(const da Y) const {
    if (x != Y.x)
      return x < Y.x;
    else if (y != Y.y)
      return y < Y.y;
    else
      return z < Y.z;
  }
} a[N];
int b[N];
struct dda {
  int x, y, id;
  bool operator<(const dda Y) const {
    if (x != Y.x)
      return x < Y.x;
    else
      return y < Y.y;
  }
};
vector<dda> v[N];
vector<pair<int, int> > vv[N];
int main() {
  int n;
  in(n);
  for (int i = 0; i < n; i++) {
    in(a[i].x, a[i].y, a[i].z);
    b[i] = a[i].x;
    a[i].id = i + 1;
  }
  sort(b, b + n);
  int len = unique(b, b + n) - b;
  for (int i = 0; i < n; i++) {
    int xx = lower_bound(b, b + len, a[i].x) - b;
    v[xx].push_back({a[i].y, a[i].z, a[i].id});
  }
  for (int i = 0; i < len; i++) {
    sort(v[i].begin(), v[i].end());
    n = v[i].size();
    for (int j = 0; j < n; j++) {
      b[j] = v[i][j].x;
    }
    sort(b, b + n);
    int llen = unique(b, b + n) - b;
    for (int j = 0; j < llen; j++) vv[j].clear();
    for (int j = 0; j < n; j++) {
      int xx = lower_bound(b, b + llen, v[i][j].x) - b;
      vv[xx].push_back(pair<int, int>(v[i][j].y, v[i][j].id));
    }
    vector<dda> vvv;
    for (int j = 0; j < llen; j++) {
      sort(vv[j].begin(), vv[j].end());
      for (int k = 1; k < vv[j].size(); k += 2) {
        printf("%d %d\n", vv[j][k].second, vv[j][k - 1].second);
      }
      if (vv[j].size() % 2) {
        vvv.push_back({b[j], vv[j].back().first, vv[j].back().second});
      }
    }
    for (int k = 1; k < vvv.size(); k += 2) {
      printf("%d %d\n", vvv[k].id, vvv[k - 1].id);
    }
    if (vvv.size() % 2) {
      bool is = 0;
      for (int j = 0; j < v[i + 1].size(); j++) {
        if (v[i + 1][j].x == vvv.back().x && v[i + 1][j].y == vvv.back().y) {
          printf("%d %d\n", vvv.back().id, v[i + 1][j].id);
          v[i + 1].erase(v[i + 1].begin() + j);
          is = 1;
          break;
        }
      }
      if (!is) v[i + 1].push_back(vvv.back());
    }
  }
  return 0;
}
