Interpolates data feature

        :param data: Array in a shape of t x nobs, where nobs = h x w x n
        :type data: numpy.ndarray
        :param times: Array of reference times relative to the first timestamp
        :type times: numpy.array
        :param resampled_times: Array of reference times relative to the first timestamp in initial timestamp array.
        :type resampled_times: numpy.array
        :return: Array of interpolated values
        :rtype: numpy.ndarray

def interpolate_data(self, data, times, resampled_times):
        """ Interpolates data feature

        :param data: Array in a shape of t x nobs, where nobs = h x w x n
        :type data: numpy.ndarray
        :param times: Array of reference times relative to the first timestamp
        :type times: numpy.array
        :param resampled_times: Array of reference times relative to the first timestamp in initial timestamp array.
        :type resampled_times: numpy.array
        :return: Array of interpolated values
        :rtype: numpy.ndarray
        """
        # pylint: disable=too-many-locals
        # get size of 2d array t x nobs
        nobs = data.shape[-1]
        if self.interpolate_pixel_wise:
            # initialise array of interpolated values
            new_data = data if self.resample_range is None else np.full(
                (len(resampled_times), nobs),
                np.nan, dtype=data.dtype)

            # Interpolate for each pixel, could be easily parallelized
            for obs in range(nobs):
                valid = ~np.isnan(data[:, obs])

                obs_interpolating_func = self.get_interpolation_function(times[valid], data[valid, obs])

                new_data[:, obs] = obs_interpolating_func(resampled_times[:, np.newaxis])

            # return interpolated values
            return new_data

        # mask representing overlap between reference and resampled times
        time_mask = (resampled_times >= np.min(times)) & (resampled_times <= np.max(times))

        # define time values as linear monotonically increasing over the observations
        const = int(self.filling_factor * (np.max(times) - np.min(times)))
        temp_values = (times[:, np.newaxis] +
                       const * np.arange(nobs)[np.newaxis, :].astype(np.float64))
        res_temp_values = (resampled_times[:, np.newaxis] +
                           const * np.arange(nobs)[np.newaxis, :].astype(np.float64))

        # initialise array of interpolated values
        new_data = np.full((len(resampled_times), nobs), np.nan, dtype=data.dtype)

        # array defining index correspondence between reference times and resampled times
        ori2res = np.array([np.abs(resampled_times - o).argmin()
                            if np.min(resampled_times) <= o <= np.max(resampled_times) else None for o in times])

        # find NaNs that start or end a time-series
        row_nans, col_nans = np.where(self._get_start_end_nans(data))
        nan_row_res_indices = np.array([index for index in ori2res[row_nans] if index is not None], dtype=np.int32)
        nan_col_res_indices = np.array([index is not None for index in ori2res[row_nans]],
                                       dtype=np.bool)
        if nan_row_res_indices.size:
            # mask out from output values the starting/ending NaNs
            res_temp_values[nan_row_res_indices, col_nans[nan_col_res_indices]] = np.nan
        # if temporal values outside the reference dates are required (extrapolation) masked them to NaN
        res_temp_values[~time_mask, :] = np.nan

        # build 1d array for interpolation. Spline functions require monotonically increasing values of x,
        # so .T is used
        input_x = temp_values.T[~np.isnan(data).T]
        input_y = data.T[~np.isnan(data).T]

        # build interpolation function
        if len(input_x) > 1:
            interp_func = self.get_interpolation_function(input_x, input_y)

            # interpolate non-NaN values in resampled time values
            new_data[~np.isnan(res_temp_values)] = interp_func(res_temp_values[~np.isnan(res_temp_values)])

        # return interpolated values
        return new_data