docs: def get_spectra(data, filt_method=dict(lb=0.1, filt_order=256),
                  spect_method=dict(NFFT=1024, n_overlap=1023, BW=2),
                  phase_zero=None, line_broadening=None, zerofill=None):
      """
      Derive the spectra from MRS data
  
      Parameters
      ----------
      data : nitime TimeSeries class instance or array
          Time-series object with data of shape (echos, transients, time-points),
          containing the FID data. If an array is provided, we will assume that a
          sampling rate of 5000.0 Hz was used
  
      filt_method : dict
          Details for the filtering method. A FIR zero phase-delay method is used
          with parameters set according to these parameters
          
      spect_method : dict
          Details for the spectral analysis. Per default, we use 
  
      line_broadening : float
          Linewidth for apodization (in Hz).
  
      zerofill : int
          Number of bins to zero fill with.
          
      Returns
      -------
      f : 
           the center frequency of the frequencies represented in the
          spectra
  
       spectrum_water, spectrum_water_suppressed: 
          The first spectrum is for the data with water not suppressed and
          the second spectrum is for the water-suppressed data.
  
      Notes
      -----
      This function performs the following operations:
  
      1. Filtering.
      2. Apodizing/windowing. Optionally, this is done with line-broadening (see
      page 92 of Keeler2005_.
      3. Spectral analysis.
          
      .. [Keeler2005] Keeler, J (2005). Understanding NMR spectroscopy, second
         edition. Wiley (West Sussex, UK).
  
      """
      if not isinstance(data, nt.TimeSeries):
         data = nt.TimeSeries(data, sampling_rate=5000.0)  
      if filt_method is not None:
          filtered = nta.FilterAnalyzer(data, **filt_method).fir
      else:
          filtered = data
      if line_broadening is not None: 
         lbr_time = line_broadening * np.pi  # Conversion from Hz to
                                             # time-constant, see Keeler page 94 
      else:
         lbr_time = 0
  
      apodized = ut.line_broadening(filtered, lbr_time)
     
      if zerofill is not None:
           new_apodized = np.concatenate([apodized.data,
                      np.zeros(apodized.shape[:-1] + (zerofill,))], -1)
  
           apodized = nt.TimeSeries(new_apodized,
                                    sampling_rate=apodized.sampling_rate)
  
      S = nta.SpectralAnalyzer(apodized,
                               method=dict(NFFT=spect_method['NFFT'],
                                           n_overlap=spect_method['n_overlap']),
                               BW=spect_method['BW'])
      
      f, c = S.spectrum_fourier
  
      return f, c