function out = process(trans)
% goes through experiment folder and processes the raw hdf5 file, and each
% experiment's hdf5 file. The last file, which looks like an experiment, is
% actually a throw-away, so it is not processed.
  
% The argument TRANS should be false if the experiment was run using the
% galene branch, which implements the start/stop transient behavior that we
% describe in the paper (which shaves 0.5 second at the start
% and end of each experiment); it should be true if the experiment was run
% using the galene-whiplash-detection branch, which does not shave off this
% time but instead produces an independent output that attempts to identify
% transients.

  fname = dir('daq*.hdf5');
  fname = fname.name;
  
  thresholds;
  model_config;
  
  out = {};
  out.v = [];
  out.c = [];
  out.s = [];
  out.e = [];
  out.t = [];
  out.p = [];
  out.r0 = [];
  out.r1 = [];
  out.r2 = [];
  out.r3 = [];
  out.dT = [];
  out.detect_trans = [];
  out.time_trans = [];
  
  hf = load_daq(fname);
  out.T = hf.t;
  [out.V, out.C, out.S, out.E] = process_units(hf);

  files = dir('fdi*.hdf5');
  for file = files(1:end-1)'
    fname = file.name;
    res = load_fdi_raw(fname, trans);

    [voltage, current, str_pot, encoder] = process_units(res);

    x_motor = encoder;                  % motor position, rad
    x_battery = str_pot;                % battery position, mm
    i = current;                        % motor current, A
    omega = (voltage - Rm*current)/Km;  % motor speed, rad/s
    pitch = res.p;
    
    %% calculate r0
    r0 = alpha_param*x_motor - x_battery;     % encoder - string pot, rad
    
    %% calculate r1 where defined
    delta_x_motor = x_motor(2:end)-x_motor(1:end-1);
    dT = res.t(2:end)-res.t(1:end-1);

    omega_mech = delta_x_motor./dT;
    r1 = alpha_param * (omega(2:end) - omega_mech);
    r3 = alpha_param * (omega_mech - sign(omega_mech)*omega_nominal);
				% fill in the first element
    r1 = [0, r1];
    r3 = [0, r3];
				% zero out at experiment edges
    idx = find(isnan(res.xmp));
    r1(idx) = 0;
    r3(idx) = 0;
    dT(idx) = 0;

    %% calculate transient mask as 0.5 second before and after idx
    if trans
      start_stop_transient = 0.5; % seconds
      idx_equivalent = 5; % look at 5 samples before and after boundary
      trans_mask = [];
      ends = [idx(2:end) - ones(1,size(idx,2)-1), size(res.t, 2)];
      boundaries = [idx(1:end); ends];
      for k = boundaries
	start = k(1);
	finish = k(2);
	for j = start:start+idx_equivalent
	  time_delta = res.t(j) - res.t(start);
	  if time_delta < start_stop_transient
	    trans_mask = [trans_mask, j];
	  end
	end
	for j = finish-idx_equivalent:finish
	  if j < 1
	    continue
	  end
	  time_delta = res.t(finish) - res.t(j);
	  if time_delta < start_stop_transient
	    trans_mask = [trans_mask, j];
	  end
	end

      end
      trans_mask_out = zeros(size(res.t));
      trans_mask_out(trans_mask) = 1;
    end
    
    %% calculate r2 %% TODO: vectorize
    r2 = nan(size(pitch));
    for j = 1:size(pitch,2)

      if omega(j) >= 0.0
	i_nominal = max(i_cutoff, i_0 + Kpitch*pitch(j));   % current - nominal, A
	r2(j) = i(j) - i_nominal;
      else
	i_nominal = min(-i_cutoff, -i_0 + Kpitch*pitch(j));
	r2(j) = i(j) - i_nominal;
      end
    end

    out.v = [out.v, voltage];
    out.c = [out.c, current];
    out.s = [out.s, str_pot];
    out.e = [out.e, encoder];
    out.t = [out.t, res.t];
    out.p = [out.p, res.p];
    out.r0 = [out.r0, r0];
    out.r1 = [out.r1, r1];
    out.r2 = [out.r2, r2];
    out.r3 = [out.r3, r3];
    out.dT = [out.dT, 0, dT];
    if trans
% if this experiment was run with transient detection, we include two
% ways to mask the output: time_trans is the usual way, used in the paper,
% and detect_trans is the other way, based on the galene-whiplash-detection
% branch.
      out.detect_trans = [out.detect_trans, res.trans];
      out.time_trans = [out.time_trans, trans_mask_out];
    end
  end
end
