function [features] = extract_features(experiment, N, options)

  switch(experiment)

    case 'mfcc_740ms'
      
      filters = filter_bank(N, options.filt1_opt);
      features = {@(x)(feval(@(x)(...
        [x;...
        circshift(x, [0 +1]); ...
        circshift(x, [0 -1])]),...
        format_scat(log_scat(...
          spec_freq_average(x, filters, options.sc1_opt)))))};
            
    case {'time_scat_l1', 'time_scat_l2', 'time_scat_l3'}
      
      Wop = wavelet_factory_1d(N, options.filt1_opt, options.sc1_opt);
      features = {@(x)(format_scat(log_scat(...
        renorm_scat(scat(x, Wop)))))};

    case {'time_freq_scat_l2', 'time_freq_scat_l2_adap_q1'}

      Wop1 = wavelet_factory_1d(N, options.filt1_opt, options.sc1_opt);      
      fWop1 = wavelet_factory_1d(...
        128, options.ffilt1_opt, options.fsc1_opt);

      scatt_fun1 = @(x)(log_scat(renorm_scat(scat(x, Wop1))));
      fscatt_fun1 = @(x)(func_output(@scat_freq, 2, scatt_fun1(x), fWop1));
      feature_fun1 = @(x)(format_scat(fscatt_fun1(x)));

      switch(experiment)

        case 'time_freq_scat_l2'

          features = {feature_fun1}

        case 'time_freq_scat_l2_adap_q1'

          Wop2 = wavelet_factory_1d(N, options.filt2_opt, options.sc2_opt);
          fWop2 = wavelet_factory_1d(...
            32, options.ffilt2_opt, options.fsc2_opt);
          
          scatt_fun2 = @(x)(log_scat(renorm_scat(scat(x, Wop2))));
          fscatt_fun2 = @(x)(func_output(@scat_freq, 2, scatt_fun2(x), fWop2));
          feature_fun2 = @(x)(format_scat(fscatt_fun2(x)));

          features = {feature_fun1, feature_fun2};

      end  

    for k = 1:length(features)
      fprintf('testing feature #%d...', k);
      tic;
      sz = size(features{k}(randn(N, 1)));
      aa = toc;
      fprintf('Ok (%.2fs) (size [%d, %d])\n', aa, sz(1), sz(2));
    end

    otherwise
        error(['Unknown experiment ', experiment, '. Aborting']);
  end

end
