function [isActive, p_vol] = brain_LM(input_files, z, isbrain,X,keep,contrast, q,useAdaptiveFDR)
%BRAIN_LM fits a GLM with AR(1) errors to the average signal of every cluster. 
%Y = hrf*X b + e, e_t= rho*e_(t-1) + white noise_t. 
%The data are created by concatenating the frames of the input files.
%
%ISACTIVE = BRAIN_LM(INPUT_FILES, Z, ISBRAIN, X, KEEP, CONTRAST, Q, USEADAPTIVEFDR)
%
%INPUT_FILES: is a cell array of names of 4D fMRI image files with multiple  frames
%in ANALYZE (.img)  format.  
%
%Z = output of BRAIN_CLUSTER
%
%ISBRAIN = output of BRAIN_THRESH
%
%X = output of FMRILM (convolved with hrf only). 
%
%KEEP = the time points (frames) to keep in the analysis from concatenated time series of input
%files.
%
%CONTRAST = a row vector of contrasts to be tested. 
%
%Q = the cluster-wise FDR level.
%
%USEADAPTIVEFDR =a flage with value 1 if adaptive procedure used, 0 if BH procedure used. 




n = length(input_files);

dm = fmris_read_image(char(input_files(1)));
rbrain = dm.data;


if n>1
    for i=2:n
        dm = fmris_read_image(char(input_files(i)));
        rbrain = cat(4,rbrain,dm.data);
        clear dm
    end
end
rbrain = rbrain(:,:,:,keep); %the time points to analyze.

arbrain = cba_sig2avg(rbrain, z,z>0);
clear  rbrain
%correlation matrix  per cluster (estimation based on TS of cluster
%average)
cor_vol = cba_cor( z,arbrain, 1:size(arbrain,4),X,z>0);
[p_vol,tstat_vol, effect_vol, sdeffect_vol] = cba_linmodel(  arbrain,X,cor_vol,z>0,contrast) ;    
isActive = cba_threshold(p_vol, z, isbrain,q, useAdaptiveFDR);
