docs: def estimate_tt_to_rectify(self, order, slitlet2d=None):
          """Estimate the polynomial transformation to rectify the image.
  
          Parameters
          ----------
          order : int
              Order of the polynomial transformation.
          slitlet2d : numpy array
              Slitlet image to be displayed with the computed boundaries
              and intersecting points overplotted. This argument is
              optional.
  
          """
  
          # protections
          if self.x_inter_orig is None \
                  or self.y_inter_orig is None \
                  or self.x_inter_rect is None \
                  or self.y_inter_rect is None:
              raise ValueError('Intersection points not computed')
  
          npoints = len(self.x_inter_orig)
          if len(self.y_inter_orig) != npoints \
                  or len(self.x_inter_rect) != npoints \
                  or len(self.y_inter_rect) != npoints:
              raise ValueError('Unexpected different number of points')
  
          # IMPORTANT: correct coordinates from origin in order to manipulate
          # coordinates corresponding to image indices
          x_inter_orig_shifted = self.x_inter_orig - self.bb_nc1_orig
          y_inter_orig_shifted = self.y_inter_orig - self.bb_ns1_orig
          x_inter_rect_shifted = self.x_inter_rect - self.bb_nc1_orig
          y_inter_rect_shifted = self.y_inter_rect - self.bb_ns1_orig
  
          # compute 2D transformation
          self.ttd_order = order
          self.ttd_aij, self.ttd_bij = compute_distortion(
              x_inter_orig_shifted, y_inter_orig_shifted,
              x_inter_rect_shifted, y_inter_rect_shifted,
              order,
              self.debugplot
          )
          self.tti_aij, self.tti_bij = compute_distortion(
              x_inter_rect_shifted, y_inter_rect_shifted,
              x_inter_orig_shifted, y_inter_orig_shifted,
              order,
              self.debugplot
          )
  
          # display slitlet with intersection points and grid indicating
          # the fitted transformation
          if abs(self.debugplot % 10) != 0 and slitlet2d is not None:
              # display image with zscale cuts
              title = "Slitlet#" + str(self.islitlet) + \
                      " (estimate_tt_to_rectify)"
              ax = ximshow(slitlet2d, title=title,
                           first_pixel=(self.bb_nc1_orig, self.bb_ns1_orig),
                           show=False)
              # intersection points
              ax.plot(self.x_inter_orig, self.y_inter_orig, 'co')
              ax.plot(self.x_inter_rect, self.y_inter_rect, 'bo')
              # grid with fitted transformation: spectrum trails
              xx = np.arange(0, self.bb_nc2_orig - self.bb_nc1_orig + 1,
                             dtype=np.float)
              for spectrail in self.list_spectrails:
                  yy0 = self.corr_yrect_a + \
                        self.corr_yrect_b * spectrail.y_rectified
                  yy = np.tile([yy0 - self.bb_ns1_orig], xx.size)
                  ax.plot(xx + self.bb_nc1_orig, yy + self.bb_ns1_orig, "b")
                  xxx, yyy = fmap(self.ttd_order, self.ttd_aij, self.ttd_bij,
                                  xx, yy)
                  ax.plot(xxx + self.bb_nc1_orig, yyy + self.bb_ns1_orig, "g")
              # grid with fitted transformation: arc lines
              ylower_line = \
                  self.list_spectrails[self.i_lower_spectrail].y_rectified
              ylower_line = self.corr_yrect_a + self.corr_yrect_b * ylower_line
              yupper_line = \
                  self.list_spectrails[self.i_upper_spectrail].y_rectified
              yupper_line = self.corr_yrect_a + self.corr_yrect_b * yupper_line
              n_points = int(yupper_line - ylower_line + 0.5) + 1
              yy = np.linspace(ylower_line - self.bb_ns1_orig,
                               yupper_line - self.bb_ns1_orig,
                               num=n_points,
                               dtype=np.float)
              for arc_line in self.list_arc_lines:
                  xline = arc_line.x_rectified - self.bb_nc1_orig
                  xx = np.array([xline] * n_points)
                  ax.plot(xx + self.bb_nc1_orig,yy + self.bb_ns1_orig, "b" )
                  xxx, yyy = fmap(self.ttd_order, self.ttd_aij, self.ttd_bij,
                                  xx, yy)
                  ax.plot(xxx + self.bb_nc1_orig, yyy + self.bb_ns1_orig, "c")
              # show plot
              pause_debugplot(self.debugplot, pltshow=True)