Modular processing reintroducing artifacts

I think [dcan_signal_processing.m](https://github.com/DCAN-Labs/dcan_bold_processing/blob/main/matlab_code/dcan_signal_processing.m) suffers from the issues described by [Lindquist et al. 2019](https://doi.org/10.1002/hbm.24528). Motion (and other confounds) are regressed from the data, and then the data are subjected to band-pass filtering. If the filters are not orthogonal to the confounds, they will re-introduce confounding. So, for example, the data after band-pass filtering may still have confounds with motion.

```matlab

rng(0)

% number of time points
N = 1000;
bp_order = 2;
lp_hz = 0.009;
hp_hz = 0.080;
F = 0.7;

% true signal, with autocorrelation
a = 0; b = .4; sigma = 1;
y0 = a / (1 - b); %eg initialize to unconditional mean of stationary time series

x = zeros(N,1);
x(1) = a + b * y0 + randn() * sigma;
for t = 2:N
x(t) = a + b * x(t-1) + randn() * sigma;
end

% can also experiment with signal that has no autocorrelation
% x = randn(N, 1);

% confounds correlated with true signal
r = x + randn(N, 1);

% observed signal
y = x + r + randn(N, 1);

% y is correlated with r (by design)
orig_cor = corr(y,r);
fprintf('original correlation: %.5f\n', orig_cor);

b = regress(y, r);
y_clean = y - r*b;

% regression has removed correlation
clean_cor = corr(y_clean, r);
fprintf('clean correlation: %.5f\n', clean_cor);

% now, apply band-pass filter
Ny = F / 2;
BW_Hz = [lp_hz hp_hz];
BW_N = BW_Hz / Ny;
[b, a] = butter(ceil(bp_order/2), BW_N);
% zero-pad the data for filtering by concatenating rows of zeros on either side of the data
padding = zeros(size(y_clean));
pad_amt = size(padding, 1);
temp = cat(1, padding, y_clean, padding);
y_clean_filtered = filtfilt(b, a, temp);
y_clean_filtered = y_clean_filtered((pad_amt + 1):(end-pad_amt), :);

% but now it's back!
filtered_cor = corr(y_clean_filtered, r);
fprintf('correlation after filtering: %.5f\n', filtered_cor);
```

```matlab
% original correlation: 0.88504
% clean correlation: -0.00005
% correlation after filtering: 0.13719
```

Fixes are straightforward. One approach would be to ensure that the band-pass filter is applied to the confounds before regression (though, be careful with the polynomial detrending).