docs: def quantile_1D(data, weights, quantile):
      """
      Compute the weighted quantile of a 1D numpy array.
  
      Parameters
      ----------
      data : ndarray
          Input array (one dimension).
      weights : ndarray
          Array with the weights of the same size of `data`.
      quantile : float
          Quantile to compute. It must have a value between 0 and 1.
  
      Returns
      -------
      quantile_1D : float
          The output value.
      """
      # Check the data
      if not isinstance(data, np.matrix):
          data = np.asarray(data)
      if not isinstance(weights, np.matrix):
          weights = np.asarray(weights)
      nd = data.ndim
      if nd != 1:
          raise TypeError("data must be a one dimensional array")
      ndw = weights.ndim
      if ndw != 1:
          raise TypeError("weights must be a one dimensional array")
      if data.shape != weights.shape:
          raise TypeError("the length of data and weights must be the same")
      if ((quantile > 1.) or (quantile < 0.)):
          raise ValueError("quantile must have a value between 0. and 1.")
      # Sort the data
      ind_sorted = np.argsort(data)
      sorted_data = data[ind_sorted]
      sorted_weights = weights[ind_sorted]
      # Compute the auxiliary arrays
      Sn = np.cumsum(sorted_weights)
      # TODO: Check that the weights do not sum zero
      #assert Sn != 0, "The sum of the weights must not be zero"
      Pn = (Sn-0.5*sorted_weights)/Sn[-1]
      # Get the value of the weighted median
      return np.interp(quantile, Pn, sorted_data)