function [data, yes, data_name] = dataread_freesurfer(data_folder,measure,smooth,atlas,vertexwise)
%% read freesurfer data
% Input
%   data_folder   - data folder. e.g. F:/xxx/
%   measure       - measurement, e.g. thickness
%   smooth        - smooth kernel, e.g. 0, 5, 10
%   atlas         - the folder containing all labels of an atlas. e.g. F:/xx/aparc
%   vertexwise    - 1 indicates vertex level, 0 otherwise
% Output
%   data          - subjid x rois
%   yes           - 1 indicates reading data successfully, 0 otherwise.

if ~vertexwise
    labels = dir([atlas,filesep,'lh',filesep,'*.label']);
    for i=1:length(labels)
        labels_lh_Ind{i} = fs_read_label([atlas,filesep,'lh',filesep,labels(i).name]);
    end
    clear labels i
    labels = dir([atlas,filesep,'rh',filesep,'*.label']);
    for i=1:length(labels)
        labels_rh_Ind{i} = fs_read_label([atlas,filesep,'rh',filesep,labels(i).name]);
    end
    clear labels i
end


subs = dir(data_folder);
delet_list = [];
for i=1:length(subs)
    if isequal(subs(i).name,'.');
        delet_list = [delet_list; i];
    elseif isequal(subs(i).name,'..');
        delet_list = [delet_list; i];
    elseif isequal(subs(i).name,'fsaverage')
        delet_list = [delet_list; i];
    elseif ~isdir([data_folder,filesep,subs(i).name])
        delet_list = [delet_list; i];
    end
end
subs(delet_list) = [];

if vertexwise
    Mgh_lh = load_mgh([data_folder,filesep,subs(1).name,filesep,'surf',filesep,['lh.',measure,'.fwhm',smooth,'.fsaverage.mgh']]);
    Mgh_rh = load_mgh([data_folder,filesep,subs(1).name,filesep,'surf',filesep,['rh.',measure,'.fwhm',smooth,'.fsaverage.mgh']]);
    data = zeros(length(subs),length(Mgh_lh)+length(Mgh_rh));
    clear Mgh_lh Mgh_rh
else
    data = zeros(length(subs),length(labels_lh_Ind)+length(labels_rh_Ind));
end

data_name = cell(length(subs),1);
for i=1:length(subs)
    Name = subs(i).name;
    data_name{i} = Name;
    fprintf('Reading subjid : %s \n',Name);
    Mgh_lh = [data_folder,filesep,Name,filesep,'surf',filesep,['lh.',measure,'.fwhm',smooth,'.fsaverage.mgh']];
    Mgh_rh = [data_folder,filesep,Name,filesep,'surf',filesep,['rh.',measure,'.fwhm',smooth,'.fsaverage.mgh']];
    if ~exist(Mgh_lh,'file')
        error('File does not exist. %s \n',Mgh_lh);
    end
    if ~exist(Mgh_rh,'file')
        error('File does not exist. %s \n',Mgh_rh);
    end
    
    if vertexwise
        dl = load_mgh(Mgh_lh)';
        dr = load_mgh(Mgh_rh)';
        data(i,:) = [dl,dr];
        clear dl dr
    else
        d1 = load_mgh(Mgh_lh)';
        for iapar = 1:length(labels_lh_Ind)
            d1_data(1,iapar) = nanmean(d1(labels_lh_Ind{iapar}));
        end
        clear iapar d1
        dr = load_mgh(Mgh_rh)';
        for iapar = 1:length(labels_rh_Ind)
            dr_data(1,iapar) = nanmean(dr(labels_rh_Ind{iapar}));
        end
        clear iapar dr
        data(i,:) = [d1_data,dr_data];
        clear d1_data dr_data
    end
    clear Mgh_lh Mgh_rh
end
yes = 1;
end