diff --git a/+bids/+internal/add_missing_field.m b/+bids/+internal/add_missing_field.m
new file mode 100644
index 00000000..9e7ae6f1
--- /dev/null
+++ b/+bids/+internal/add_missing_field.m
@@ -0,0 +1,5 @@
+function structure = add_missing_field(structure, field)
+  if ~isfield(structure, field)
+    structure(1).(field) = '';
+  end
+end
diff --git a/+bids/+internal/append_to_layout.m b/+bids/+internal/append_to_layout.m
new file mode 100644
index 00000000..11c3a76e
--- /dev/null
+++ b/+bids/+internal/append_to_layout.m
@@ -0,0 +1,91 @@
+function subject = append_to_layout(file, subject, modality, schema)
+  %
+  % appends a file to the BIDS layout by parsing it according to the provided schema
+  %
+  % USAGE::
+  %
+  %   subject = append_to_layout(file, subject, modality, schema == [])
+  %
+  % :param file:
+  % :type  file: string
+  % :param subject: subject sub-structure from the BIDS layout
+  % :type  subject: strcture
+  % :param modality:
+  % :type  modality: string
+  % :param schema:
+  % :type  schema: strcture
+  %
+  %
+  % Copyright (C) 2021--, BIDS-MATLAB developers
+
+  if ~exist('schema', 'var')
+    schema = [];
+  end
+
+  % Parse file fist to identify the suffix group in the template.
+  % Then reparse the file using the entity-label pairs defined in the schema.
+  p = bids.internal.parse_filename(file);
+
+  idx = find_suffix_group(modality, p.suffix, schema);
+
+  if ~isempty(schema)
+
+    if isempty(idx)
+      warning('append_to_structure:noMatchingSuffix', ...
+              'Skipping file with no valid suffix in schema: %s', file);
+      return
+    end
+
+    entities = bids.schema.return_modality_entities(schema.datatypes.(modality)(idx), schema);
+    p = bids.internal.parse_filename(file, entities);
+
+  end
+
+  % Check any new entity field that needs to be added into the layout or the output
+  % of the parsing to make sure the 2 structures can be concatenated
+  if ~isempty(subject.(modality))
+
+    [subject.(modality), p] = bids.internal.match_structure_fields(subject.(modality), p);
+
+  end
+
+  if isempty(subject.(modality))
+    subject.(modality) = p;
+  else
+    subject.(modality)(end + 1, 1) = p;
+  end
+
+end
+
+function idx = find_suffix_group(modality, suffix, schema)
+
+  idx = [];
+
+  if isempty(schema)
+    return
+  end
+
+  % the following loop could probably be improved with some cellfun magic
+  %   cellfun(@(x, y) any(strcmp(x,y)), {p.type}, suffix_groups)
+  for i = 1:size(schema.datatypes.(modality), 1)
+
+    this_suffix_group = schema.datatypes.(modality)(i);
+
+    % for CI
+    if iscell(this_suffix_group)
+      this_suffix_group = this_suffix_group{1};
+    end
+
+    if any(strcmp(suffix, this_suffix_group.suffixes))
+      idx = i;
+      break
+    end
+
+  end
+
+  if isempty(idx)
+    warning('findSuffix:noMatchingSuffix', ...
+            'No corresponding suffix in schema for %s for datatype %s', suffix, modality);
+  end
+
+end
diff --git a/+bids/+internal/file_utils.m b/+bids/+internal/file_utils.m
index 70b8083b..c565803f 100644
--- a/+bids/+internal/file_utils.m
+++ b/+bids/+internal/file_utils.m
@@ -282,7 +282,6 @@
     files = dirs;
 
   else
-
     t = regexp(files, expr);
 
     if numel(files) == 1 && ~iscell(t)
diff --git a/+bids/+internal/get_metadata.m b/+bids/+internal/get_metadata.m
index 612658f1..106ddf39 100644
--- a/+bids/+internal/get_metadata.m
+++ b/+bids/+internal/get_metadata.m
@@ -21,7 +21,7 @@
   N = 3;
 
   % -There is a session level in the hierarchy
-  if isfield(p, 'ses') && ~isempty(p.ses)
+  if isfield(p.entities, 'ses') && ~isempty(p.entities.ses)
     N = N + 1;
   end
 
@@ -31,7 +31,7 @@
 
     % -List the potential metadata files associated with this file suffix type
     % Default is to assume it is a JSON file
-    metafile = bids.internal.file_utils('FPList', pth, sprintf(pattern, p.type));
+    metafile = bids.internal.file_utils('FPList', pth, sprintf(pattern, p.suffix));
 
     if isempty(metafile)
       metafile = {};
@@ -44,13 +44,17 @@
     for i = 1:numel(metafile)
 
       p2 = bids.internal.parse_filename(metafile{i});
-      fn = setdiff(fieldnames(p2), {'filename', 'ext', 'type'});
+      entities = {};
+      if isfield(p2, 'entities')
+        entities = fieldnames(p2.entities);
+      end
 
       % -Check if this metadata file contains the same entity-label pairs as its
       % data file counterpart
       ismeta = true;
-      for j = 1:numel(fn)
-        if ~isfield(p, fn{j}) || ~strcmp(p.(fn{j}), p2.(fn{j}))
+      for j = 1:numel(entities)
+        if ~isfield(p.entities, entities{j}) || ...
+                ~strcmp(p.entities.(entities{j}), p2.entities.(entities{j}))
           ismeta = false;
           break
         end
@@ -73,9 +77,15 @@
 
   end
 
-  % ==========================================================================
-  % -Inheritance principle
-  % ==========================================================================
+  if isempty(meta)
+    warning('No metadata for %s', filename);
+  end
+
+end
+
+% ==========================================================================
+% -Inheritance principle
+% ==========================================================================
 function s1 = update_metadata(s1, s2, file)
   if isempty(s2)
     return
@@ -88,3 +98,4 @@
       s1.(fn{i}) = s2.(fn{i});
     end
   end
+end
diff --git a/+bids/+internal/keep_file_for_query.m b/+bids/+internal/keep_file_for_query.m
new file mode 100644
index 00000000..391116ef
--- /dev/null
+++ b/+bids/+internal/keep_file_for_query.m
@@ -0,0 +1,32 @@
+function status = keep_file(file_struct, options)
+
+  status = true;
+
+  % suffix is treated separately as it is not one of the entities
+  for l = 1:size(options, 1)
+    if strcmp(options{l, 1}, 'suffix') && ~ismember(file_struct.suffix, options{l, 2})
+      status = false;
+      return
+    end
+  end
+
+  for l = 1:size(options, 1)
+
+    if ~strcmp(options{l, 1}, 'suffix')
+
+      if ~ismember(options{l, 1}, fieldnames(file_struct.entities))
+        status = false;
+        break
+      end
+
+      if isfield(file_struct.entities, options{l, 1}) && ...
+              ~ismember(file_struct.entities.(options{l, 1}), options{l, 2})
+        status = false;
+        break
+      end
+
+    end
+
+  end
+
+end
diff --git a/+bids/+internal/match_structure_fields.m b/+bids/+internal/match_structure_fields.m
new file mode 100644
index 00000000..1234f84a
--- /dev/null
+++ b/+bids/+internal/match_structure_fields.m
@@ -0,0 +1,14 @@
+function [s1, s2] = match_structure_fields(s1, s2)
+
+  missing_fields = setxor(fieldnames(s1), fieldnames(s2));
+
+  if ~isempty(missing_fields)
+    for iField = 1:numel(missing_fields)
+
+      s1 = bids.internal.add_missing_field(s1, missing_fields{iField});
+      s2 = bids.internal.add_missing_field(s2, missing_fields{iField});
+
+    end
+  end
+
+end
diff --git a/+bids/+internal/parse_filename.m b/+bids/+internal/parse_filename.m
index 681a2c8e..ebaba2bf 100644
--- a/+bids/+internal/parse_filename.m
+++ b/+bids/+internal/parse_filename.m
@@ -1,23 +1,35 @@
 function p = parse_filename(filename, fields)
+  %
   % Split a filename into its building constituents
-  % FORMAT p = bids.internal.parse_filename(filename, fields)
+  %
+  % USAGE::
+  %
+  %   p = bids.internal.parse_filename(filename, fields)
+  %
+  % :param filename: fielname to parse that follows the pattern
+  %                  ``sub-label[_entity-label]*_suffix.extension``
+  % :type  filename: string
+  % :param fields:   cell of strings of the entities to use for parsing
+  % :type  fields:   cell
   %
   % Example:
   %
-  % >> filename = '../sub-16/anat/sub-16_ses-mri_run-1_echo-2_FLASH.nii.gz';
-  % >> bids.internal.parse_filename(filename)
+  %   filename = '../sub-16/anat/sub-16_ses-mri_run-1_acq-hd_T1w.nii.gz';
   %
-  % ans =
+  %   bids.internal.parse_filename(filename)
+  %
+  %   ans =
   %
   %   struct with fields:
   %
-  %     filename: 'sub-16_ses-mri_run-1_echo-2_FLASH.nii.gz'
-  %         type: 'FLASH'
-  %          ext: '.nii.gz'
-  %          sub: '16'
-  %          ses: 'mri'
-  %          run: '1'
-  %         echo: '2'
+  %     'filename', 'sub-16_ses-mri_run-1_acq-hd_T1w.nii.gz', ...
+  %     'suffix', 'T1w', ...
+  %     'ext', '.nii.gz', ...
+  %     'entities', struct('sub', '16', ...
+  %                        'ses', 'mri', ...
+  %                        'run', '1', ...
+  %                        'acq', 'hd');
+  %
   % __________________________________________________________________________
 
   % Copyright (C) 2016-2018, Guillaume Flandin, Wellcome Centre for Human Neuroimaging
@@ -26,31 +38,31 @@
   filename = bids.internal.file_utils(filename, 'filename');
 
   % -Identify all the BIDS entity-label pairs present in the filename (delimited by "_")
-  % https://bids-specification.readthedocs.io/en/stable/99-appendices/04-entity-table.html
   [parts, dummy] = regexp(filename, '(?:_)+', 'split', 'match'); %#ok<ASGLU>
   p.filename = filename;
 
   % -Identify the suffix and extension of this file
-  % https://bids-specification.readthedocs.io/en/stable/02-common-principles.html#file-name-structure
-  [p.type, p.ext] = strtok(parts{end}, '.');
+  [p.suffix, p.ext] = strtok(parts{end}, '.');
 
   % -Separate the entity from the label for each pair identified above
   for i = 1:numel(parts) - 1
     [d, dummy] = regexp(parts{i}, '(?:\-)+', 'split', 'match'); %#ok<ASGLU>
-    p.(d{1}) = d{2};
+    p.entities.(d{1}) = d{2};
   end
 
   % -Extra fields can be added to the structure and ordered specifically.
   if nargin == 2
     for i = 1:numel(fields)
-      if ~isfield(p, fields{i})
-        p.(fields{i}) = '';
-      end
+      p.entities = bids.internal.add_missing_field(p.entities, fields{i});
     end
     try
-      p = orderfields(p, ['filename', 'ext', 'type', fields]);
+      p = orderfields(p, {'filename', 'ext', 'suffix', 'entities'});
+      p.entities = orderfields(p.entities, fields);
     catch
-      warning('Ignoring file ''%s'' not matching template.', filename);
+      warning('bidsMatlab:noMatchingTemplate', ...
+              'Ignoring file %s not matching template.', filename);
       p = struct([]);
     end
   end
+
+end
diff --git a/+bids/+internal/return_modality_extensions.m b/+bids/+internal/return_modality_extensions.m
new file mode 100644
index 00000000..8a1d4eea
--- /dev/null
+++ b/+bids/+internal/return_modality_extensions.m
@@ -0,0 +1,19 @@
+function extensions = return_modality_extensions(modality)
+
+  extensions = '(';
+
+  % for CI
+  if iscell(modality)
+    modality = modality{1};
+  end
+
+  for iExt = 1:numel(modality.extensions)
+    if ~strcmp(modality.extensions{iExt}, '.json')
+      extensions = [extensions,  modality.extensions{iExt}, '|']; %#ok<AGROW>
+    end
+  end
+
+  % Replace final "|" by a "){1}"
+  extensions(end:end + 3) = '){1}';
+
+end
diff --git a/+bids/+internal/return_modality_regular_expression.m b/+bids/+internal/return_modality_regular_expression.m
new file mode 100644
index 00000000..938e8f72
--- /dev/null
+++ b/+bids/+internal/return_modality_regular_expression.m
@@ -0,0 +1,8 @@
+function regular_expression = return_modality_regular_expression(modality)
+
+  suffixes = bids.internal.return_modality_suffixes(modality);
+  extensions = bids.internal.return_modality_extensions(modality);
+
+  regular_expression = ['^%s.*' suffixes extensions '$'];
+
+end
diff --git a/+bids/+internal/return_modality_suffixes.m b/+bids/+internal/return_modality_suffixes.m
new file mode 100644
index 00000000..10d81750
--- /dev/null
+++ b/+bids/+internal/return_modality_suffixes.m
@@ -0,0 +1,17 @@
+function suffixes = return_modality_suffixes(modality)
+
+  suffixes = '_(';
+
+  % For CI
+  if iscell(modality)
+    modality = modality{1};
+  end
+
+  for iExt = 1:numel(modality(:).suffixes)
+    suffixes = [suffixes,  modality.suffixes{iExt}, '|']; %#ok<AGROW>
+  end
+
+  % Replace final "|" by a "){1}"
+  suffixes(end:end + 3) = '){1}';
+
+end
diff --git a/+bids/+schema/load_schema.m b/+bids/+schema/load_schema.m
new file mode 100644
index 00000000..3f656e69
--- /dev/null
+++ b/+bids/+schema/load_schema.m
@@ -0,0 +1,79 @@
+function schema = load_schema(use_schema)
+  % Loads a json schema by recursively looking through a folder structure.
+  %
+  % The nesting of the output structure reflects a combination of the folder structure and
+  % any eventual nesting within each json.
+  %
+  %
+  % Copyright (C) 2021--, BIDS-MATLAB developers
+
+  % TODO:
+  %  - folders that do not contain json files themselves but contain
+  %  subfolders that do, are not reflected in the output structure (they are
+  %  skipped). This can lead to "name conflicts". See "silenced" unit tests
+  %  for more info.
+
+  if nargin < 1
+    use_schema = true();
+  end
+
+  if ~use_schema
+    schema = [];
+    return
+  end
+
+  if ischar(use_schema)
+    schema_dir = use_schema;
+  else
+    schema_dir = fullfile(fileparts(mfilename('fullpath')), '..', '..', 'schema');
+  end
+
+  if ~exist(schema_dir, 'dir')
+    error('The schema directory %s does not exist.', schema_dir);
+  end
+
+  schema = struct();
+
+  [json_file_list, dirs] = bids.internal.file_utils('FPList', schema_dir, '^.*.json$');
+
+  schema = append_json_content_to_structure(schema, json_file_list);
+
+  schema = inspect_subdir(schema, dirs);
+
+end
+
+function structure = append_json_content_to_structure(structure, json_file_list)
+
+  for iFile = 1:size(json_file_list, 1)
+
+    file = deblank(json_file_list(iFile, :));
+
+    field_name = bids.internal.file_utils(file, 'basename');
+
+    structure.(field_name) = bids.util.jsondecode(file);
+  end
+
+end
+
+function structure = inspect_subdir(structure, subdir_list)
+  % recursively inspects subdirectory for json files and reflects folder
+  % hierarchy in the output structure.
+
+  for iDir = 1:size(subdir_list, 1)
+
+    directory = deblank(subdir_list(iDir, :));
+
+    [json_file_list, dirs] = bids.internal.file_utils('FPList', directory, '^.*.json$');
+
+    if ~isempty(json_file_list)
+      field_name = bids.internal.file_utils(directory, 'basename');
+      structure.(field_name) = struct();
+      structure.(field_name) = append_json_content_to_structure(structure.(field_name), ...
+                                                                json_file_list);
+    end
+
+    structure = inspect_subdir(structure, dirs);
+
+  end
+
+end
diff --git a/+bids/+schema/return_modalities.m b/+bids/+schema/return_modalities.m
new file mode 100644
index 00000000..07643c2b
--- /dev/null
+++ b/+bids/+schema/return_modalities.m
@@ -0,0 +1,13 @@
+function modalities = return_modalities(subject, schema, modality_group)
+
+  % if we go schema-less we list directories in the subject/session folder
+  % as proxy of the modalities that we have to parse
+  modalities = cellstr(bids.internal.file_utils('List', ...
+                                                subject.path, ...
+                                                'dir', ...
+                                                '.*'));
+  if ~isempty(schema)
+    modalities = schema.modalities.(modality_group).datatypes;
+  end
+
+end
diff --git a/+bids/+schema/return_modality_entities.m b/+bids/+schema/return_modality_entities.m
new file mode 100644
index 00000000..d90a9720
--- /dev/null
+++ b/+bids/+schema/return_modality_entities.m
@@ -0,0 +1,14 @@
+function entities = return_modality_entities(suffix_group, schema)
+
+  % for CI
+  if iscell(suffix_group)
+    suffix_group = suffix_group{1};
+  end
+
+  entity_names = fieldnames(suffix_group.entities);
+
+  for i = 1:size(entity_names, 1)
+    entities{1, i} = schema.entities.(entity_names{i}).entity; %#ok<*AGROW>
+  end
+
+end
diff --git a/+bids/+schema/return_modality_groups.m b/+bids/+schema/return_modality_groups.m
new file mode 100644
index 00000000..3651f232
--- /dev/null
+++ b/+bids/+schema/return_modality_groups.m
@@ -0,0 +1,9 @@
+function modality_groups = return_modality_groups(schema)
+
+  % dummy variable if we go schema less
+  modality_groups = {nan()};
+  if ~isempty(schema)
+    modality_groups = fieldnames(schema.modalities);
+  end
+
+end
diff --git a/+bids/layout.m b/+bids/layout.m
index 872fde01..4d3f0fcf 100644
--- a/+bids/layout.m
+++ b/+bids/layout.m
@@ -1,9 +1,23 @@
-function BIDS = layout(root, tolerant)
+function BIDS = layout(root, use_schema)
+  %
   % Parse a directory structure formated according to the BIDS standard
-  % FORMAT BIDS = bids.layout(root)
-  % root     - directory formated according to BIDS [Default: pwd]
-  % tolerant - if set to 0 (default) only files g
-  % BIDS     - structure containing the BIDS file layout
+  %
+  % USAGE::
+  %
+  %   BIDS = bids.layout(root = pwd, use_schema = false)
+  %
+  % :param root:       directory of the dataset formated according to BIDS [default: ``pwd``]
+  % :type  root:       string
+  % :param use_schema: If set to ``true`` (default), the parsing of the dataset
+  %                    will follow the bids-schema provided with bids-matlab.
+  %                    If set to ``false`` files just have to be of the form
+  %                    ``sub-label_[entity-label]_suffix.ext`` to be parsed.
+  %                    If a folder path is provided, then the schema contained
+  %                    in that folder willl be used for parsing.
+  % :type  use_schema: boolean
+  %
+  %
+
   % __________________________________________________________________________
   %
   % BIDS (Brain Imaging Data Structure): https://bids.neuroimaging.io/
@@ -11,7 +25,7 @@
   %   describing outputs of neuroimaging experiments.
   %   K. J. Gorgolewski et al, Scientific Data, 2016.
   % __________________________________________________________________________
-
+  %
   % Copyright (C) 2016-2018, Guillaume Flandin, Wellcome Centre for Human Neuroimaging
   % Copyright (C) 2018--, BIDS-MATLAB developers
 
@@ -19,31 +33,40 @@
   % ==========================================================================
   if ~nargin
     root = pwd;
+
   elseif nargin == 1
+
     if ischar(root)
       root = bids.internal.file_utils(root, 'CPath');
+
     elseif isstruct(root)
-      BIDS = root; % or BIDS = bids.layout(root.root);
+      BIDS = root; % for bids.query
       return
+
     else
       error('Invalid syntax.');
+
     end
+
   elseif nargin > 2
     error('Too many input arguments.');
+
+  end
+
+  if ~exist('use_schema', 'var')
+    use_schema = true;
   end
 
-  if ~exist('tolerant', 'var')
-    tolerant = false;
+  if ~exist(root, 'dir')
+    error('BIDS directory does not exist: ''%s''', root);
   end
 
-  % -BIDS structure
+  % BIDS structure
   % ==========================================================================
 
   % BIDS.dir          -- BIDS directory
   % BIDS.description  -- content of dataset_description.json
   % BIDS.sessions     -- cellstr of sessions
-  % BIDS.scans        -- for sub-<participant_label>_scans.tsv (should go within subjects)
-  % BIDS.sess         -- for sub-<participants_label>_sessions.tsv (should go within subjects)
   % BIDS.participants -- for participants.tsv
   % BIDS.subjects'    -- structure array of subjects
 
@@ -51,716 +74,232 @@
                 'dir', root, ...
                 'description', struct([]), ...
                 'sessions', {{}}, ...
-                'scans', struct([]), ...
-                'sess', struct([]), ...
                 'participants', struct([]), ...
                 'subjects', struct([]));
 
-  % -Validation of BIDS root directory
-  % ==========================================================================
-  if ~exist(BIDS.dir, 'dir')
-    error('BIDS directory does not exist: ''%s''', BIDS.dir);
-
-  elseif ~exist(fullfile(BIDS.dir, 'dataset_description.json'), 'file')
-
-    msg = sprintf('BIDS directory not valid: missing dataset_description.json: ''%s''', ...
-                  BIDS.dir);
-
-    tolerant_message(tolerant, msg);
-
-  end
-
-  % -Dataset description
-  % ==========================================================================
-  try
-    BIDS.description = bids.util.jsondecode(fullfile(BIDS.dir, 'dataset_description.json'));
-  catch err
-    msg = sprintf('BIDS dataset description could not be read: %s', err.message);
-    tolerant_message(tolerant, msg);
-  end
-
-  fields_to_check = {'BIDSVersion', 'Name'};
-  for iField = 1:numel(fields_to_check)
-
-    if ~isfield(BIDS.description, fields_to_check{iField})
-      msg = sprintf( ...
-                    'BIDS dataset description not valid: missing %s field.', ...
-                    fields_to_check{iField});
-      tolerant_message(tolerant, msg);
-    end
-
-  end
+  BIDS = validate_description(BIDS, use_schema);
 
-  % -Optional directories
+  % Optional directories
   % ==========================================================================
-  % [code/]
+  % [code/] - ignore
   % [derivatives/]
-  % [stimuli/]
-  % [sourcedata/]
+  % [stimuli/] - ingore
+  % [sourcedata/] - ignore
   % [phenotype/]
 
-  % -Scans key file
-  % ==========================================================================
-
-  % sub-<participant_label>/[ses-<session_label>/]
-  %     sub-<participant_label>_scans.tsv
+  BIDS.participants = [];
+  BIDS.participants = manage_tsv(BIDS.participants, BIDS.dir, 'participants.tsv');
 
-  % See also optional README and CHANGES files
-
-  % -Participant key file
+  % Subjects
   % ==========================================================================
-  p = bids.internal.file_utils('FPList', BIDS.dir, '^participants\.tsv$');
-  if ~isempty(p)
-    try
-      BIDS.participants = bids.util.tsvread(p);
-    catch
-      msg = ['unable to read ' p];
-      tolerant_message(tolerant, msg);
-    end
-  end
-  p = bids.internal.file_utils('FPList', BIDS.dir, '^participants\.json$');
-  if ~isempty(p)
-    BIDS.participants.meta = bids.util.jsondecode(p);
+  subjects = cellstr(bids.internal.file_utils('List', BIDS.dir, 'dir', '^sub-.*$'));
+  if isequal(subjects, {''})
+    error('No subjects found in BIDS directory.');
   end
 
-  % -Sessions file
-  % ==========================================================================
-
-  % sub-<participant_label>/[ses-<session_label>/]
-  %      sub-<participant_label>[_ses-<session_label>]_sessions.tsv
+  schema = bids.schema.load_schema(use_schema);
 
-  % -Tasks: JSON files are accessed through metadata
-  % ==========================================================================
-  % t = bids.internal.file_utils('FPList',BIDS.dir,...
-  %    '^task-.*_(beh|bold|events|channels|physio|stim|meg)\.(json|tsv)$');
-
-  % -Subjects
-  % ==========================================================================
-  sub = cellstr(bids.internal.file_utils('List', BIDS.dir, 'dir', '^sub-.*$'));
-  if isequal(sub, {''})
-    error('No subjects found in BIDS directory.');
-  end
+  for iSub = 1:numel(subjects)
+    sessions = cellstr(bids.internal.file_utils('List', ...
+                                                fullfile(BIDS.dir, subjects{iSub}), ...
+                                                'dir', ...
+                                                '^ses-.*$'));
 
-  for iSub = 1:numel(sub)
-    sess = cellstr(bids.internal.file_utils('List', ...
-                                            fullfile(BIDS.dir, sub{iSub}), ...
-                                            'dir', ...
-                                            '^ses-.*$'));
-                                        
-    for iSess = 1:numel(sess)
+    for iSess = 1:numel(sessions)
       if isempty(BIDS.subjects)
-        BIDS.subjects = parse_subject(BIDS.dir, sub{iSub}, sess{iSess});
+        BIDS.subjects = parse_subject(BIDS.dir, subjects{iSub}, sessions{iSess}, schema);
       else
-        BIDS.subjects(end + 1) = parse_subject(BIDS.dir, sub{iSub}, sess{iSess});
+        new_subject = parse_subject(BIDS.dir, subjects{iSub}, sessions{iSess}, schema);
+
+        [BIDS.subjects, new_subject] = bids.internal.match_structure_fields(BIDS.subjects, ...
+                                                                            new_subject);
+        % TODO: this can be added to "match_structure_fields"
+        BIDS.subjects(end + 1) = new_subject;
+
       end
     end
-    
-  end
 
-end
-
-function tolerant_message(tolerant, msg)
-  if tolerant
-    warning(msg);
-  else
-    error(msg);
   end
-end
 
-% ==========================================================================
-% -Parse a subject's directory
-% ==========================================================================
-function subject = parse_subject(pth, subjname, sesname)
+end
 
+function subject = parse_subject(pth, subjname, sesname, schema)
+  %
+  % Parse a subject's directory
+  %
   % For each modality (anat, func, eeg...) all the files from the
-  % corresponding directory are listed and their filenames parsed with extra
-  % BIDS valid entities listed (e.g. 'acq','ce','rec','fa'...).
+  % corresponding directory are listed and their filenames parsed with
+  % BIDS valid entities as listed in the schema (if the schema is not empty).
 
   subject.name    = subjname;   % subject name ('sub-<participant_label>')
   subject.path    = fullfile(pth, subjname, sesname); % full path to subject directory
-  subject.session = sesname; % session name ('' or 'ses-<label>')
-  subject.anat    = struct([]); % anatomy imaging data
-  subject.func    = struct([]); % task imaging data
-  subject.fmap    = struct([]); % fieldmap data
-  subject.beh     = struct([]); % behavioral experiment data
-  subject.dwi     = struct([]); % diffusion imaging data
-  subject.perf    = struct([]); % ASL perfusion imaging data
-  subject.eeg     = struct([]); % EEG data
-  subject.meg     = struct([]); % MEG data
-  subject.ieeg    = struct([]); % iEEG data
-  subject.pet     = struct([]); % PET imaging data
-
-  subject = parse_anat(subject);
-  subject = parse_func(subject);
-  subject = parse_fmap(subject);
-  subject = parse_eeg(subject);
-  subject = parse_meg(subject);
-  subject = parse_beh(subject);
-  subject = parse_dwi(subject);
-  subject = parse_perf(subject);
-  subject = parse_pet(subject);
-  subject = parse_ieeg(subject);
-
-end
-
-function subject = parse_anat(subject)
-
-  % --------------------------------------------------------------------------
-  % -Anatomy imaging data
-  % --------------------------------------------------------------------------
-  pth = fullfile(subject.path, 'anat');
-
-  if exist(pth, 'dir')
-
-    entities = return_entities('anat');
-
-    file_list = return_file_list('anat', subject);
-
-    for i = 1:numel(file_list)
-
-      subject = append_to_structure(file_list{i}, entities, subject, 'anat');
-
+  subject.session = sesname;    % session name ('' or 'ses-<label>')
+  subject.scans   = struct([]); % for sub-<participant_label>_scans.tsv
+  subject.sess    = struct([]); % for sub-<participants_label>_sessions.tsv
+
+  modality_groups = bids.schema.return_modality_groups(schema);
+
+  for iGroup = 1:numel(modality_groups)
+
+    modalities = bids.schema.return_modalities(subject, schema, modality_groups{iGroup});
+
+    % if we go schema-less, we pass an empty schema to all the parsing functions
+    % so the parsing is unconstrained
+    for iModality = 1:numel(modalities)
+      switch modalities{iModality}
+        case {'anat', 'func', 'beh', 'meg', 'eeg', 'ieeg', 'pet'}
+          subject = parse_using_schema(subject, modalities{iModality}, schema);
+        case 'dwi'
+          subject = parse_dwi(subject, schema);
+        case 'fmap'
+          subject = parse_fmap(subject, schema);
+        case 'perf'
+          subject = parse_perf(subject, schema);
+        otherwise
+          % in case we are going schemaless
+          % and the modality is not one of the usual suspect
+          subject.(modalities{iModality}) = struct([]);
+          subject = parse_using_schema(subject, modalities{iModality}, []);
+      end
     end
+
   end
 
 end
 
-function subject = parse_func(subject)
+function subject = parse_using_schema(subject, modality, schema)
 
-  % --------------------------------------------------------------------------
-  % -Task imaging data
-  % --------------------------------------------------------------------------
-  pth = fullfile(subject.path, 'func');
+  pth = fullfile(subject.path, modality);
 
   if exist(pth, 'dir')
 
-    entities = return_entities('func');
-
-    file_list = return_file_list('func', subject);
-
-    for i = 1:numel(file_list)
-
-      subject = append_to_structure(file_list{i}, entities, subject, 'func');
-
-      subject.func(end).meta = struct([]); % ?
-
-    end
-
-    file_list = return_event_file_list('func', subject);
-
-    for i = 1:numel(file_list)
-
-      subject = append_to_structure(file_list{i}, entities, subject, 'func');
+    subject = bids.internal.add_missing_field(subject, modality);
 
-      subject.func(end).meta = bids.util.tsvread(fullfile(pth, file_list{i})); % ?
-
-    end
-
-    file_list = return_physio_stim_file_list('func', subject);
+    file_list = return_file_list(modality, subject);
 
     for i = 1:numel(file_list)
 
-      subject = append_to_structure(file_list{i}, entities, subject, 'func');
-
-      subject.func(end).meta = struct([]); % ?
-
-    end
-  end
-end
+      subject = bids.internal.append_to_layout(file_list{i}, subject, modality, schema);
 
-function subject = parse_perf(subject)
+      if ~isempty(subject.(modality)) && strcmp(subject.(modality)(end).ext, '.tsv')
+        % events
+        % stim
+        % channels
+        % electrodes
+        %
+        % does not cover physio.tsv.gz or stim.tsv.gz
 
-  % --------------------------------------------------------------------------
-  % -ASL perfusion imaging data
-  % --------------------------------------------------------------------------
-  pth = fullfile(subject.path, 'perf');
+        subject.(modality)(end).content = [];
+        subject.(modality)(end).meta = [];
 
-  if exist(pth, 'dir')
+        subject.(modality)(end) = manage_tsv( ...
+                                             subject.(modality)(end), ...
+                                             pth, ...
+                                             subject.(modality)(end).filename);
 
-    entities = return_entities('perf');
-
-    file_list = return_file_list('perf', subject);
-
-    file_list = convert_to_cell(file_list);
-
-    j = 1;
-
-    % ASL timeseries NIfTI file
-    % ----------------------------------------------------------------------
-    labels = regexp(file_list, [ ...
-                                '^sub-[a-zA-Z0-9]+' ...              % sub-<participant_label>
-                                '_asl\.nii(\.gz)?$'], 'names'); % NIfTI file suffix/extension
-
-    if any(~cellfun(@isempty, labels))
-
-      idx = find(~cellfun(@isempty, labels));
-
-      for i = 1:numel(idx)
-
-        % Parse filename
-        % ---------------------------
-        fb = bids.internal.file_utils(bids.internal.file_utils(file_list{idx(i)}, ...
-                                                               'basename'), ...
-                                      'basename');
-
-        p = bids.internal.parse_filename(file_list{idx(i)}, entities);
-
-        subject = append_to_perf(subject, p, j);
-
-        % add type
-        subject.perf(j).type = 'asl';
-
-        % Manage JSON-sidecar metadata (REQUIRED)
-        % ---------------------------
-        metafile = fullfile(pth, bids.internal.file_utils(fb, 'ext', 'json'));
-
-        if exist(metafile, 'file')
-          [~, Ffile] = fileparts(metafile);
-          subject.perf(j).json_sidecar_filename = [Ffile '.json'];
-          subject.perf(j).meta = bids.util.jsondecode(metafile);
-        else
-          warning(['Missing: ' metafile]);
-
-        end
-
-        % Manage ASLCONTEXT-sidecar metadata (REQUIRED)
-        % ---------------------------
-        metafile = fullfile(pth, bids.internal.file_utils([fb(1:end - 4) '_aslcontext'], 'ext', 'tsv'));
-
-        if exist(metafile, 'file')
-          [~, Ffile] = fileparts(metafile);
-          subject.perf(j).context_sidecar_filename = [Ffile '.tsv'];
-          subject.perf(j).context = bids.util.tsvread(metafile);
-
-        else
-          warning(['Missing: ' metafile]);
-
-        end
-
-        % Manage M0 (REQUIRED)
-        % ---------------------------
-        % M0 field is flexible:
-
-        if ~isfield(subject.perf(j).meta, 'M0Type')
-          warning(['M0Type field missing in ' subject.perf(j).json_sidecar_filename]);
-
-        else
-
-          m0_type = [];
-          m0_explanation = [];
-          m0_volume_index = [];
-          m0_value = [];
-
-          switch subject.perf(j).meta.M0Type
-
-            case 'Separate'
-              % the M0 was obtained as a separate scan
-              m0_type = 'separate_scan';
-              m0_explanation = 'M0 was obtained as a separate scan';
-
-              % M0scan.nii filename
-              % assuming the (.nii|.nii.gz) extension choice is the same throughout
-              m0_filename = [subject.perf(j).filename(1:end - 10) 'm0scan' subject.perf(j).ext];
-              if ~exist(fullfile(pth, m0_filename), 'file')
-                warning(['Missing: ' m0_filename]);
-              else
-                % subject.perf(j).m0_filename = m0_filename;
-                % -> this is included in the same structure for the m0scan.nii
-              end
-
-              % M0 sidecar filename
-              m0_json_sidecar_filename = [subject.perf(j).filename(1:end - 10) 'm0scan.json'];
-              if ~exist(fullfile(pth, m0_json_sidecar_filename), 'file')
-                warning(['Missing: ' m0_json_sidecar_filename]);
-              else
-                % subject.perf(j).m0_json_sidecar_filename = m0_json_sidecar_filename;
-                % -> this is included in the same structure for the m0scan.nii
-              end
-
-            case 'Included'
-              % M0 is one or more image(s) in the *asl.nii[.gz] timeseries
-              if ~isfield(subject.perf(j), 'context') || ...
-                      ~isfield(subject.perf(j).context, 'volume_type')
-                warning('Cannot find M0 volume in aslcontext, context-information missing');
-
-              else
-                m0indices = find(cellfun(@(x) strcmp(x, 'm0scan'), ...
-                                         subject.perf(j).context.volume_type) == true);
-                if isempty(m0indices)
-                  warning('No M0 volume found in aslcontext');
-
-                else
-                  m0_type = 'within_timeseries';
-                  m0_explanation = 'M0 is one or more image(s) in the *asl.nii[.gz] timeseries';
-                  m0_volume_index = m0indices;
-
-                end
-              end
-
-            case 'Estimate'
-              m0_type = 'single_value';
-              m0_explanation = [ ...
-                                'this is a single estimated M0 value, ', ...
-                                'e.g. when the M0 is obtained from an external scan and/or study'];
-              m0_value = subject.perf(j).meta.M0;
-
-            case 'Absent'
-              m0_type = 'use_control_as_m0';
-              m0_explanation = [ ...
-                                'M0 is absent, so we can use the (average) control volume ', ...
-                                'as pseudo-M0 (if no background suppression was used)'];
-
-              if subject.perf(j).meta.BackgroundSuppression == true
-                warning('Caution when using control as M0, background suppression was applied');
-              end
-
-            otherwise
-              warning(['Unknown M0Type:', ...
-                       subject.perf(j).meta.M0Type, ...
-                       ' in ', ...
-                       subject.perf(j).json_sidecar_filename]);
-          end
-
-          if ~isempty(m0_type)
-            subject.perf(j).m0_type = m0_type;
-            subject.perf(j).m0_explanation = m0_explanation;
-          end
-
-          if ~isempty(m0_volume_index)
-            subject.perf(j).m0_volume_index = m0_volume_index;
-          end
-          if ~isempty(m0_value)
-            subject.perf(j).m0_value = m0_value;
-          end
-
-        end
-
-        % Manage labeling image metadata (OPTIONAL)
-        % ---------------------------
-        metafile = fullfile(pth, bids.internal.file_utils([fb(1:end - 4) '_labeling'], 'ext', 'jpg'));
-
-        if exist(metafile, 'file')
-          [~, Ffile] = fileparts(metafile);
-          subject.perf(j).labeling_image_filename = [Ffile '.jpg'];
-        end
-
-        j = j + 1;
-      end % for i = 1:numel(idx)
-    end % if any(~cellfun(@isempty, labels))
-
-    % -M0scan NIfTI file
-    % ---------------------------------------------------------------------
-    labels = regexp(file_list, [ ...
-                                '^sub-[a-zA-Z0-9]+' ...              % sub-<participant_label>
-                                '_m0scan\.nii(\.gz)?$'], 'names'); % NIfTI file suffix/extension
-
-    if any(~cellfun(@isempty, labels))
-      idx = find(~cellfun(@isempty, labels));
-      for i = 1:numel(idx)
-        % Parse filename
-        % ---------------------------
-        fb = bids.internal.file_utils(bids.internal.file_utils(file_list{idx(i)}, 'basename'), 'basename');
-        p = bids.internal.parse_filename(file_list{idx(i)}, entities);
-
-        subject = append_to_perf(subject, p, j);
-
-        % Manage JSON-sidecar metadata (REQUIRED)
-        % ---------------------------
-        metafile = fullfile(pth, bids.internal.file_utils(fb, 'ext', 'json'));
-
-        if ~exist(metafile, 'file')
-          warning(['Missing: ' metafile]);
-
-        else
-          [~, Ffile] = fileparts(metafile);
-          subject.perf(j).json_sidecar_filename = [Ffile '.json'];
-          subject.perf(j).meta = bids.util.jsondecode(metafile);
-
-          % Manage intended-for (REQUIRED)
-          % ---------------------------
-
-          % Get all NIfTIs that this m0scan is intended for
-          path_intended_for = {};
-          if ~isfield(subject.perf(j).meta, 'IntendedFor')
-            warning(['Missing field IntendedFor in ' metafile]);
-
-          elseif ischar(subject.perf(j).meta.IntendedFor)
-            path_intended_for{1} = subject.perf(j).meta.IntendedFor;
-
-          elseif isstruct(subject.perf(j).meta.IntendedFor)
-            for iPath = 1:length(subject.perf(j).meta.IntendedFor)
-              path_intended_for{iPath} = subject.perf(j).meta.IntendedFor(iPath); %#ok<*AGROW>
-            end
-
-          end
-
-          for iPath = 1:length(path_intended_for)
-            % check if this NIfTI is not missing
-            if ~exist(fullfile(fileparts(pth), path_intended_for{iPath}), 'file')
-              warning(['Missing: ' fullfile(fileparts(pth), path_intended_for{iPath})]);
-
-            else
-              % also check that this NIfTI aims to the same m0scan
-              [~, path2check, ext2check] = fileparts(path_intended_for{iPath});
-              filename_found = max(arrayfun(@(x) strcmp(x.filename, [path2check ext2check]), subject.perf));
-              if ~filename_found
-                warning(['Did not find NIfTI for which is intended: ' subject.perf(j).filename]);
-
-              else
-                subject.perf(j).intended_for = path_intended_for{iPath};
-
-              end
-            end
-          end
-        end
-      end % for i = 1:numel(idx)
-      j = j + 1;
-    end  % if any(~cellfun(@isempty, labels))
-
-  end % if exist(pth, 'dir')
-end % function subject = parse_perf(subject)
+      end
 
-function subject = append_to_perf(subject, p, j)
+      % case {'photo', 'coordsystem'}
 
-  if j == 1
-    subject.perf = p;
-  else
-    fields_p = fieldnames(p);
-    for iField = 1:length(fields_p)
-      subject.perf(j).(fields_p{iField}) = p.(fields_p{iField});
     end
-  end
 
-  % default to run 1 ((!) TODO: but could be that we need to check this
-  % at the end!)
-  if isempty(subject.perf(j).run)
-    subject.perf(j).run = '1';
   end
 
 end
 
-function subject = parse_fmap(subject)
-  %
-  % TODO:
-  %
-  % 20210114 - From Remi:
-  % For other modalities, metadata are fetched upon query.
-  % It is unclear why we do it differently for fmaps
+function subject = parse_dwi(subject, schema)
 
   % --------------------------------------------------------------------------
-  % -Fieldmap data
+  %  Diffusion imaging data
   % --------------------------------------------------------------------------
-  pth = fullfile(subject.path, 'fmap');
-
-  if exist(pth, 'dir')
-
-    file_list = return_file_list('fmap', subject);
-
-    j = 1;
-
-    % -Phase difference image and at least one magnitude image
-    % ----------------------------------------------------------------------
-    labels = return_labels_fieldmap(file_list, 'phase_difference_image');
-
-    if any(~cellfun(@isempty, labels))
 
-      idx = find(~cellfun(@isempty, labels));
-
-      for i = 1:numel(idx)
-
-        subject.fmap(j).type = 'phasediff';
-        subject.fmap(j).filename = file_list{idx(i)};
-        subject.fmap(j).magnitude = { ...
-                                     strrep(file_list{idx(i)}, ...
-                                            '_phasediff.nii', ...
-                                            '_magnitude1.nii'), ...
-                                     strrep(file_list{idx(i)}, ...
-                                            '_phasediff.nii', ...
-                                            '_magnitude2.nii')}; % optional
-
-        subject = append_common_fmap_fields_to_structure(subject, labels{idx(i)}, j);
+  modality = 'dwi';
+  pth = fullfile(subject.path, modality);
 
-        metafile = return_fmap_metadata_file(subject, file_list{idx(i)});
-        subject.fmap(j).meta = struct([]);
-        % (!) TODO: file can also be stored at higher levels (inheritance principle)
-        if ~isempty(metafile)
-          subject.fmap(j).meta = bids.util.jsondecode(metafile);
-        end
+  if exist(pth, 'dir')
 
-        j = j + 1;
+    subject = bids.internal.add_missing_field(subject, modality);
 
-      end
-    end
+    file_list = return_file_list(modality, subject);
 
-    % -Two phase images and two magnitude images
-    % ----------------------------------------------------------------------
-    labels = return_labels_fieldmap(file_list, 'two_phase_image');
-
-    if any(~cellfun(@isempty, labels))
-
-      idx = find(~cellfun(@isempty, labels));
-
-      for i = 1:numel(idx)
-
-        subject.fmap(j).type = 'phase12';
-        subject.fmap(j).filename = { ...
-                                    file_list{idx(i)}, ...
-                                    strrep(file_list{idx(i)}, ...
-                                           '_phase1.nii', ...
-                                           '_phase2.nii')};
-        subject.fmap(j).magnitude = { ...
-                                     strrep(file_list{idx(i)}, ...
-                                            '_phase1.nii', ...
-                                            '_magnitude1.nii'), ...
-                                     strrep(file_list{idx(i)}, ...
-                                            '_phase1.nii', ...
-                                            '_magnitude2.nii')};
-
-        subject = append_common_fmap_fields_to_structure(subject, labels{idx(i)}, j);
-
-        metafile = return_fmap_metadata_file(subject, file_list{idx(i)});
-        subject.fmap(j).meta = struct([]);
-        % (!) TODO: file can also be stored at higher levels (inheritance principle)
-        if ~isempty(metafile)
-          subject.fmap(j).meta = { ...
-                                  bids.util.jsondecode(metafile), ...
-                                  bids.util.jsondecode(strrep(metafile, ...
-                                                              '_phase1.json', ...
-                                                              '_phase2.json'))};
-        end
+    for i = 1:numel(file_list)
 
-        j = j + 1;
+      subject = bids.internal.append_to_layout(file_list{i}, subject, modality, schema);
 
+      % bval & bvec file
+      % ------------------------------------------------------------------
+      % TODO: they can also be stored at higher levels (inheritance principle)
+      bvalfile = bids.internal.get_metadata(file_list{i}, '^.*%s\\.bval$');
+      if isfield(bvalfile, 'filename')
+        subject.dwi(end).bval = bids.util.tsvread(bvalfile.filename);
       end
 
-    end
-
-    % -A single, real fieldmap image
-    % ----------------------------------------------------------------------
-    labels = return_labels_fieldmap(file_list, 'fieldmap_image');
-
-    if any(~cellfun(@isempty, labels))
-
-      idx = find(~cellfun(@isempty, labels));
-
-      for i = 1:numel(idx)
-
-        subject.fmap(j).type = 'fieldmap';
-        subject.fmap(j).filename = file_list{idx(i)};
-        subject.fmap(j).magnitude = strrep(file_list{idx(i)}, ...
-                                           '_fieldmap.nii', ...
-                                           '_magnitude.nii');
-
-        subject = append_common_fmap_fields_to_structure(subject, labels{idx(i)}, j);
-
-        metafile = return_fmap_metadata_file(subject, file_list{idx(i)});
-        subject.fmap(j).meta = struct([]);
-        % (!) TODO: file can also be stored at higher levels (inheritance principle)
-        if ~isempty(metafile)
-          subject.fmap(j).meta = bids.util.jsondecode(metafile);
-        end
-
-        j = j + 1;
-
+      bvecfile = bids.internal.get_metadata(file_list{i}, '^.*%s\\.bvec$');
+      if isfield(bvalfile, 'filename')
+        subject.dwi(end).bvec = bids.util.tsvread(bvecfile.filename);
       end
-    end
-
-    % -Multiple phase encoded directions (topup)
-    % ----------------------------------------------------------------------
-    labels = return_labels_fieldmap(file_list, 'phase_encoded_direction_image');
 
-    if any(~cellfun(@isempty, labels))
-
-      idx = find(~cellfun(@isempty, labels));
-
-      for i = 1:numel(idx)
-
-        subject.fmap(j).filename = file_list{idx(i)};
-        if ~isempty(regexp(subject.fmap(j).filename, 'm0scan', 'ONCE'))
-          subject.fmap(j).type = 'm0scan';
-        else
-          subject.fmap(j).type = 'epi';
-        end
-        subject.fmap(j).dir = labels{idx(i)}.dir;
-
-        subject = append_common_fmap_fields_to_structure(subject, labels{idx(i)}, j);
-
-        metafile = return_fmap_metadata_file(subject, file_list{idx(i)});
-        subject.fmap(j).meta = struct([]);
-        % (!) TODO: file can also be stored at higher levels (inheritance principle)
-        if ~isempty(metafile)
-          subject.fmap(j).meta = bids.util.jsondecode(metafile);
-        end
-
-        j = j + 1;
-
-      end
     end
   end
-
 end
 
-function subject = parse_eeg(subject)
+function subject = parse_perf(subject, schema)
+
   % --------------------------------------------------------------------------
-  % -EEG data
+  % ASL perfusion imaging data
   % --------------------------------------------------------------------------
-  pth = fullfile(subject.path, 'eeg');
+
+  modality = 'perf';
+  pth = fullfile(subject.path, 'perf');
 
   if exist(pth, 'dir')
 
-    entities = return_entities('eeg');
+    subject = bids.internal.add_missing_field(subject, modality);
 
-    file_list = return_file_list('eeg', subject);
+    file_list = return_file_list(modality, subject);
 
     for i = 1:numel(file_list)
 
-      % European data format (.edf)
-      % BrainVision Core Data Format (.vhdr, .vmrk, .eeg) by Brain Products GmbH
-      % The format used by the MATLAB toolbox EEGLAB (.set and .fdt files)
-      % Biosemi data format (.bdf)
-
-      p = bids.internal.parse_filename(file_list{i}, entities);
-      switch p.ext
-        case {'.edf', '.vhdr', '.set', '.bdf'}
-          % each recording is described with a single file,
-          % even though the data can consist of multiple
-          subject.eeg = [subject.eeg p];
-          subject.eeg(end).meta = struct([]); % ?
-        case {'.vmrk', '.eeg', '.fdt'}
-          % skip the additional files that come with certain data formats
-        otherwise
-          % skip unknown files
-      end
-
-    end
-
-    file_list = return_event_file_list('eeg', subject);
+      subject = bids.internal.append_to_layout(file_list{i}, subject, modality, schema);
 
-    for i = 1:numel(file_list)
+      switch subject.perf(i).suffix
 
-      subject = append_to_structure(file_list{i}, entities, subject, 'eeg');
+        case 'asl'
 
-      subject.eeg(end).meta = bids.util.tsvread(fullfile(pth, file_list{i})); % ?
+          subject.perf(i).meta = [];
+          subject.perf(i).dependencies = [];
 
-    end
+          subject.perf(i).meta = bids.internal.get_metadata( ...
+                                                            fullfile( ...
+                                                                     subject.path, ...
+                                                                     modality, ...
+                                                                     file_list{i}));
 
-    file_list = return_channel_description_file_list('eeg', subject);
+          aslcontext_file = strrep(subject.perf(i).filename, ...
+                                   ['_asl' subject.perf(i).ext], ...
+                                   '_aslcontext.tsv');
+          subject.perf(i).dependencies.context = manage_tsv( ...
+                                                            struct('content', [], 'meta', []), ...
+                                                            pth, ...
+                                                            aslcontext_file);
 
-    for i = 1:numel(file_list)
+          subject.perf(i) = manage_asllabeling(subject.perf(i), pth);
 
-      subject = append_to_structure(file_list{i}, entities, subject, 'eeg');
+          subject.perf(i) = manage_M0(subject.perf(i), pth);
 
-      subject.eeg(end).meta = bids.util.tsvread(fullfile(pth, file_list{i})); % ?
+        case 'm0scan'
 
-    end
+          subject.perf(i).intended_for = [];
 
-    file_list = return_session_specific_file_list('eeg', subject);
+          subject.perf(i).meta = bids.internal.get_metadata( ...
+                                                            fullfile( ...
+                                                                     subject.path, ...
+                                                                     modality, ...
+                                                                     file_list{i}));
 
-    for i = 1:numel(file_list)
-
-      subject = append_to_structure(file_list{i}, entities, subject, 'eeg');
+          subject.perf(i) = manage_intended_for(subject.perf(i), subject, pth);
 
-      subject.eeg(end).meta = struct([]); % ?
+      end
 
     end
 
@@ -768,53 +307,58 @@ function tolerant_message(tolerant, msg)
 
 end
 
-function subject = parse_meg(subject)
-  % --------------------------------------------------------------------------
-  % -MEG data
-  % --------------------------------------------------------------------------
-  pth = fullfile(subject.path, 'meg');
-
-  if exist(pth, 'dir')
-
-    entities = return_entities('meg');
-
-    file_list = return_file_list('meg', subject);
-
-    for i = 1:numel(file_list)
-
-      subject = append_to_structure(file_list{i}, entities, subject, 'meg');
-
-      subject.meg(end).meta = struct([]); % ?
-
-    end
-
-    file_list = return_event_file_list('meg', subject);
-
-    for i = 1:numel(file_list)
-
-      subject = append_to_structure(file_list{i}, entities, subject, 'meg');
-
-      subject.meg(end).meta = bids.util.tsvread(fullfile(pth, file_list{i})); % ?
-
-    end
-
-    file_list = return_channel_description_file_list('meg', subject);
-
-    for i = 1:numel(file_list)
+function subject = parse_fmap(subject, schema)
 
-      subject = append_to_structure(file_list{i}, entities, subject, 'meg');
+  modality = 'fmap';
+  pth = fullfile(subject.path, modality);
 
-      subject.meg(end).meta = bids.util.tsvread(fullfile(pth, file_list{i})); % ?
+  if exist(pth, 'dir')
 
-    end
+    subject = bids.internal.add_missing_field(subject, modality);
 
-    file_list = return_session_specific_file_list('meg', subject);
+    file_list = return_file_list(modality, subject);
 
     for i = 1:numel(file_list)
 
-      subject = append_to_structure(file_list{i}, entities, subject, 'meg');
+      subject = bids.internal.append_to_layout(file_list{i}, subject, modality, schema);
+
+      subject.fmap(i).meta = bids.internal.get_metadata( ...
+                                                        fullfile( ...
+                                                                 subject.path, ...
+                                                                 modality, ...
+                                                                 file_list{i}));
+      %       subject.perf(i).intended_for = [];
+      %       subject.fmap(i) = manage_intended_for(subject.fmap(i), subject, pth);
+
+      switch subject.fmap(i).suffix
+
+        % -A single, real fieldmap image
+        case {'fieldmap', 'magnitude'}
+          subject.fmap(i).dependencies.magnitude = strrep(file_list{idx(i)}, ...
+                                                          '_fieldmap.nii', ...
+                                                          '_magnitude.nii');
+
+          % Phase difference image and at least one magnitude image
+        case {'phasediff'}
+          subject.fmap(i).dependencies.magnitude = { ...
+                                                    strrep(file_list{i}, ...
+                                                           '_phasediff.nii', ...
+                                                           '_magnitude1.nii'), ...
+                                                    strrep(file_list{i}, ...
+                                                           '_phasediff.nii', ...
+                                                           '_magnitude2.nii')}; % optional
+
+          % Two phase images and two magnitude images
+        case {'phase1', 'phase2'}
+          subject.fmap(i).dependencies.magnitude = { ...
+                                                    strrep(file_list{i}, ...
+                                                           '_phase1.nii', ...
+                                                           '_magnitude1.nii'), ...
+                                                    strrep(file_list{i}, ...
+                                                           '_phase1.nii', ...
+                                                           '_magnitude2.nii')};
 
-      subject.meg(end).meta = struct([]); % ?
+      end
 
     end
 
@@ -822,140 +366,50 @@ function tolerant_message(tolerant, msg)
 
 end
 
-function subject = parse_beh(subject)
-  % --------------------------------------------------------------------------
-  % -Behavioral experiments data
-  %
-  % - Event timing, metadata, physiological and other continuous recordings
-  % --------------------------------------------------------------------------
-  pth = fullfile(subject.path, 'beh');
-
-  if exist(pth, 'dir')
+% --------------------------------------------------------------------------
+%                            HELPER FUNCTIONS
+% --------------------------------------------------------------------------
 
-    entities = return_entities('beh');
+function BIDS = validate_description(BIDS, use_schema)
 
-    file_list = return_file_list('beh', subject);
+  if ~exist(fullfile(BIDS.dir, 'dataset_description.json'), 'file')
 
-    for i = 1:numel(file_list)
+    msg = sprintf('BIDS directory not valid: missing dataset_description.json: ''%s''', ...
+                  BIDS.dir);
 
-      subject = append_to_structure(file_list{i}, entities, subject, 'beh');
+    tolerant_message(use_schema, msg);
 
-    end
   end
-end
-
-function subject = parse_dwi(subject)
-  % --------------------------------------------------------------------------
-  % -Diffusion imaging data
-  % --------------------------------------------------------------------------
-  pth = fullfile(subject.path, 'dwi');
-
-  if exist(pth, 'dir')
-
-    entities = return_entities('dwi');
-
-    file_list = return_file_list('dwi', subject);
-
-    for i = 1:numel(file_list)
-
-      subject = append_to_structure(file_list{i}, entities, subject, 'dwi');
-
-      % -bval file
-      % ------------------------------------------------------------------
-      % bval file can also be stored at higher levels (inheritance principle)
-      bvalfile = bids.internal.get_metadata(file_list{i}, '^.*%s\\.bval$');
-      if isfield(bvalfile, 'filename')
-        subject.dwi(end).bval = bids.util.tsvread(bvalfile.filename); % ?
-      end
-
-      % -bvec file
-      % ------------------------------------------------------------------
-      % bvec file can also be stored at higher levels (inheritance principle)
-      bvecfile = bids.internal.get_metadata(file_list{i}, '^.*%s\\.bvec$');
-      if isfield(bvalfile, 'filename')
-        subject.dwi(end).bvec = bids.util.tsvread(bvecfile.filename); % ?
-      end
-
-    end
+  try
+    BIDS.description = bids.util.jsondecode(fullfile(BIDS.dir, 'dataset_description.json'));
+  catch err
+    msg = sprintf('BIDS dataset description could not be read: %s', err.message);
+    tolerant_message(use_schema, msg);
   end
-end
-
-function subject = parse_pet(subject)
-  % --------------------------------------------------------------------------
-  % -Positron Emission Tomography imaging data
-  % --------------------------------------------------------------------------
-  pth = fullfile(subject.path, 'pet');
-
-  if exist(pth, 'dir')
-
-    entities = return_entities('pet');
 
-    file_list = return_file_list('pet', subject);
-
-    for i = 1:numel(file_list)
-
-      subject = append_to_structure(file_list{i}, entities, subject, 'pet');
+  fields_to_check = {'BIDSVersion', 'Name'};
+  for iField = 1:numel(fields_to_check)
 
+    if ~isfield(BIDS.description, fields_to_check{iField})
+      msg = sprintf( ...
+                    'BIDS dataset description not valid: missing %s field.', ...
+                    fields_to_check{iField});
+      tolerant_message(use_schema, msg);
     end
-  end
-end
-
-function subject = parse_ieeg(subject)
-  % --------------------------------------------------------------------------
-  % -Human intracranial electrophysiology
-  % --------------------------------------------------------------------------
-  pth = fullfile(subject.path, 'ieeg');
-
-  if exist(pth, 'dir')
-
-    entities = return_entities('ieeg');
-
-    file_list = return_file_list('ieeg', subject);
-
-    for i = 1:numel(file_list)
 
-      % European Data Format (.edf)
-      % BrainVision Core Data Format (.vhdr, .eeg, .vmrk) by Brain Products GmbH
-      % The format used by the MATLAB toolbox EEGLAB (.set and .fdt files)
-      % Neurodata Without Borders (.nwb)
-      % MEF3 (.mef)
-
-      p = bids.internal.parse_filename(file_list{i}, entities);
-      switch p.ext
-        case {'.edf', '.vhdr', '.set', '.nwb', '.mef'}
-          % each recording is described with a single file,
-          % even though the data can consist of multiple
-          subject.ieeg = [subject.ieeg p];
-          subject.ieeg(end).meta = struct([]); % ?
-        case {'.vmrk', '.eeg', '.fdt'}
-          % skip the additional files that come with certain data formats
-        otherwise
-          % skip unknown files
-      end
-
-    end
+    % TODO
+    % Add warning if bids version does not match schema version
 
   end
 
 end
 
-% --------------------------------------------------------------------------
-%                            HELPER FUNCTIONS
-% --------------------------------------------------------------------------
-
-function subject = append_to_structure(file, entities, subject, modality)
-
-  p = bids.internal.parse_filename(file, entities);
-  subject.(modality) = [subject.(modality) p];
-
-end
-
-function subject = append_common_fmap_fields_to_structure(subject, labels, idx)
-
-  subject.fmap(idx).ses = regexprep(labels.ses, '^_[a-zA-Z0-9]+-', '');
-  subject.fmap(idx).acq = regexprep(labels.acq, '^_[a-zA-Z0-9]+-', '');
-  subject.fmap(idx).run = regexprep(labels.run, '^_[a-zA-Z0-9]+-', '');
-
+function tolerant_message(use_schema, msg)
+  if use_schema
+    error(msg);
+  else
+    warning(msg);
+  end
 end
 
 function f = convert_to_cell(f)
@@ -966,243 +420,231 @@ function tolerant_message(tolerant, msg)
   end
 end
 
-function entities = return_entities(modality)
-
-  switch modality
-
-    case 'anat'
-      entities = {'sub', 'ses', 'acq', 'ce', 'rec', 'fa', 'echo', 'inv', 'run'};
+function file_list = return_file_list(modality, subject)
 
-    case 'func'
-      entities = {'sub', ...
-                  'ses', ...
-                  'task', ...
-                  'acq', ...
-                  'rec', ...
-                  'fa', ...
-                  'echo', ...
-                  'dir', ...
-                  'inv', ...
-                  'run', ...
-                  'recording', ...
-                  'meta'};
+  % We list anything but json files
 
-    case {'eeg', 'ieeg'}
-      entities = {'sub', 'ses', 'task', 'acq', 'run', 'meta'};
+  % TODO
+  % it should be possible to create some of those patterns for the regexp
+  % based on some of the required entities written down in the schema
 
-    case 'meg'
-      entities = {'sub', 'ses', 'task', 'acq', 'run', 'proc', 'meta'};
+  % TODO
+  % this does not cover coordsystem.json
 
-    case 'beh'
-      entities = {'sub', 'ses', 'task'};
+  % jn to omit json but not .pos file for headshape.pos
+  pattern = '_([a-zA-Z0-9]+){1}\\..*[^jn]';
 
-    case 'dwi'
-      entities = {'sub', 'ses', 'acq', 'run', 'bval', 'bvec'};
+  pth = fullfile(subject.path, modality);
 
-    case 'pet'
-      entities = {'sub', 'ses', 'task', 'acq', 'rec', 'run'};
+  [file_list, d] = bids.internal.file_utils('List', ...
+                                            pth, ...
+                                            sprintf(['^%s.*' pattern '$'], ...
+                                                    subject.name));
 
-    case 'perf'
-      entities = {'sub', 'ses', 'acq', 'dir', 'rec', 'run'};
+  file_list = convert_to_cell(file_list);
 
+  if strcmp(modality, 'meg') && ~isempty(d)
+    for i = 1:size(d, 1)
+      file_list{end + 1, 1} = d(i, :);
+    end
   end
-end
-
-function labels = return_labels_fieldmap(file_list, fiefmap_type)
-
-  direction_pattern = '';
 
-  switch fiefmap_type
+end
 
-    case 'phase_difference_image'
-      suffix = 'phasediff';
+function structure = manage_tsv(structure, pth, filename)
 
-    case 'two_phase_image'
-      suffix = 'phase1';
+  p = bids.internal.file_utils('FPList', pth,  ['^' strrep(filename, '.tsv', '\.tsv') '$']);
 
-    case 'fieldmap_image'
-      suffix = 'fieldmap';
+  if isempty(p)
+    warning('Missing: %s', fullfile(pth, filename));
 
-    case 'phase_encoded_direction_image'
-      suffix = '(epi|m0scan)';
+  else
+    structure.content = bids.util.tsvread(p);
 
-      direction_pattern = '_dir-(?<dir>[a-zA-Z0-9]+)?';
+    p = bids.internal.file_utils('FPList', pth,  ['^' strrep(filename, '.tsv', '\.json') '$']);
+    if ~isempty(p)
+      structure.meta = bids.util.jsondecode(p);
+    end
 
   end
 
-  labels = regexp(file_list, [ ...
-                              '^sub-[a-zA-Z0-9]+', ...          % sub-<participant_label>
-                              '(?<ses>_ses-[a-zA-Z0-9]+)?', ... % ses-<label>
-                              '(?<acq>_acq-[a-zA-Z0-9]+)?', ... % acq-<label>
-                              direction_pattern, ...            % dir-<index>
-                              '(?<run>_run-[a-zA-Z0-9]+)?', ... % run-<index>
-                              '_' suffix '\.nii(\.gz)?$'], 'names');   % NIfTI file extension
-
 end
 
-% TODO
-%
-% more refactoring can be done across the several 'return_X_file_list' functions
-%
-
-function file_list = return_file_list(modality, subject)
-
-  switch modality
-
-    case 'anat'
-      pattern = '_([a-zA-Z0-9]+){1}\\.nii(\\.gz)?';
-
-    case 'func'
-      pattern = '_task-.*_bold\\.nii(\\.gz)?';
-
-    case 'fmap'
-      pattern = '\\.nii(\\.gz)?';
+function structure = manage_intended_for(structure, subject, pth)
 
-    case 'eeg'
-      pattern = '_task-.*_eeg\\..*[^json]';
+  if isempty(structure.meta)
+    return
 
-    case 'meg'
-      pattern = '_task-.*_meg\\..*[^json]';
+  else
 
-    case 'beh'
-      pattern = '_(events\\.tsv|beh\\.json|physio\\.tsv\\.gz|stim\\.tsv\\.gz)';
+    % Get all NIfTIs that this m0scan is intended for
+    path_intended_for = {};
+    if ~isfield(structure.meta, 'IntendedFor')
+      warning('Missing field IntendedFor for %s', structure.filename);
 
-    case 'dwi'
-      pattern = '_([a-zA-Z0-9]+){1}\\.nii(\\.gz)?';
+    elseif ischar(structure.meta.IntendedFor)
+      path_intended_for{1} = structure.meta.IntendedFor;
 
-    case 'pet'
-      pattern = '_task-.*_pet\\.nii(\\.gz)?';
+    elseif isstruct(structure.meta.IntendedFor)
+      for iPath = 1:length(structure.meta.IntendedFor)
+        path_intended_for{iPath} = structure.meta.IntendedFor(iPath); %#ok<*AGROW>
+      end
 
-    case 'ieeg'
-      pattern = '_task-.*_ieeg\\..*[^json]';
+    end
 
-    case 'perf'
-      pattern = '_(asl|m0scan)\\.nii(\\.gz)?';
+    for iPath = 1:length(path_intended_for)
+      % check if this NIfTI is not missing
+      if ~exist(fullfile(fileparts(pth), path_intended_for{iPath}), 'file')
+        warning(['Missing: ' fullfile(fileparts(pth), path_intended_for{iPath})]);
 
-  end
+      else
+        % also check that this NIfTI aims to the same m0scan
+        [~, path2check, ext2check] = fileparts(path_intended_for{iPath});
+        filename_found = max(arrayfun(@(x) strcmp(x.filename, ...
+                                                  [path2check ext2check]), ...
+                                      subject.perf));
+        if ~filename_found
+          warning(['Did not find NIfTI for which is intended: ' structure.filename]);
 
-  pth = fullfile(subject.path, modality);
+        else
+          structure.intended_for = path_intended_for{iPath};
 
-  [file_list, d] = bids.internal.file_utils('List', ...
-                                            pth, ...
-                                            sprintf(['^%s.*' pattern '$'], ...
-                                                    subject.name));
+        end
+      end
+    end
 
-  if strcmp(modality, 'meg') && isempty(file_list)
-    file_list = d;
   end
 
-  file_list = convert_to_cell(file_list);
-
 end
 
-function file_list = return_event_file_list(modality, subject)
-  %
-  % TODO: events file can also be stored at higher levels (inheritance principle)
-  %
+function perf = manage_asllabeling(perf, pth)
+  % labeling image metadata (OPTIONAL)
+  % ---------------------------
+  metafile = fullfile(pth, strrep(perf.filename, ...
+                                  ['_asl' perf.ext], ...
+                                  '_asllabeling.jpg'));
 
-  switch modality
-
-    case {'func', 'eeg', 'meg'}
-      pattern = '_task-.*_events\\.tsv';
+  if exist(metafile, 'file')
+    [~, Ffile] = fileparts(metafile);
+    perf.dependencies.labeling_image.filename = [Ffile '.jpg'];
 
   end
 
-  pth = fullfile(subject.path, modality);
-
-  [file_list, d] = bids.internal.file_utils('List', ...
-                                            pth, ...
-                                            sprintf(['^%s.*' pattern '$'], ...
-                                                    subject.name));
-
-  file_list = convert_to_cell(file_list);
-
 end
 
-function file_list = return_physio_stim_file_list(modality, subject)
-  %
-  % Physiological and other continuous recordings file
-  %
-  % TODO: stim files can also be stored at higher levels (inheritance principle)
-  %
-
-  switch modality
+function perf = manage_M0(perf, pth)
 
-    case {'func'}
-      pattern = '_task-.*_(physio|stim)\\.tsv\\.gz';
+  % M0 field is flexible:
 
-  end
-
-  pth = fullfile(subject.path, modality);
+  if ~isfield(perf.meta, 'M0Type')
+    warning('M0Type field missing for %s', perf.filename);
 
-  [file_list, d] = bids.internal.file_utils('List', ...
-                                            pth, ...
-                                            sprintf(['^%s.*' pattern '$'], ...
-                                                    subject.name));
+  else
 
-  file_list = convert_to_cell(file_list);
+    m0_type = [];
+    m0_explanation = [];
+    m0_volume_index = [];
+    m0_value = [];
+    m0_filename = [];
+    m0_sidecar = [];
 
-end
+    switch perf.meta.M0Type
 
-function metafile = return_fmap_metadata_file(subject, fmap_file)
+      case 'Separate'
+        % the M0 was obtained as a separate scan
+        m0_type = 'separate_scan';
+        m0_explanation = 'M0 was obtained as a separate scan';
 
-  pth = fullfile(subject.path, 'fmap');
+        % M0scan.nii filename
+        % assuming the (.nii|.nii.gz) extension choice is the same throughout
+        m0_filename = strrep(perf.filename, ...
+                             ['_asl' perf.ext], ...
+                             ['_m0scan' perf.ext]);
 
-  fb = bids.internal.file_utils(bids.internal.file_utils( ...
-                                                         fmap_file, ...
-                                                         'basename'), ...
-                                'basename');
-  metafile = fullfile(pth, bids.internal.file_utils(fb, 'ext', 'json'));
+        if ~exist(fullfile(pth, m0_filename), 'file')
+          warning(['Missing: ' m0_filename]);
+        else
+          % subject.perf(j).m0_filename = m0_filename;
+          % -> this is included in the same structure for the m0scan.nii
+        end
 
-  if ~exist(metafile, 'file')
-    metafile = [];
-  end
+        % M0 sidecar filename
+        m0_sidecar = strrep(perf.filename, ...
+                            ['_asl' perf.ext], ...
+                            '_m0scan.json');
 
-end
+        if ~exist(fullfile(pth, m0_sidecar), 'file')
+          warning(['Missing: ' m0_sidecar]);
 
-function file_list = return_channel_description_file_list(modality, subject)
-  %
-  % Channel description table
-  %
-  % TODO: those files can also be stored at higher levels (inheritance principle)
-  %
+        else
+          % subject.perf(j).m0_json_sidecar_filename = m0_json_sidecar_filename;
+          % -> this is included in the same structure for the m0scan.nii
+        end
 
-  switch modality
+      case 'Included'
+        % M0 is one or more image(s) in the *asl.nii[.gz] timeseries
+        if ~isfield(perf.dependencies, 'context') || ...
+                ~isfield(perf.dependencies.context.content, 'volume_type')
+          warning('Cannot find M0 volume in aslcontext, context-information missing');
 
-    case {'eeg', 'meg'}
-      pattern = '_task-.*_channels\\.tsv';
+        else
+          m0indices = find(cellfun(@(x) strcmp(x, 'm0scan'), ...
+                                   perf.dependencies.context.content.volume_type) == true);
 
-  end
+          if isempty(m0indices)
+            warning('No M0 volume found in aslcontext');
 
-  pth = fullfile(subject.path, modality);
+          else
+            m0_type = 'within_timeseries';
+            m0_explanation = 'M0 is one or more image(s) in the *asl.nii[.gz] timeseries';
+            m0_volume_index = m0indices;
 
-  [file_list, d] = bids.internal.file_utils('List', ...
-                                            pth, ...
-                                            sprintf(['^%s.*' pattern '$'], ...
-                                                    subject.name));
+          end
+        end
 
-  file_list = convert_to_cell(file_list);
+      case 'Estimate'
+        m0_type = 'single_value';
+        m0_explanation = [ ...
+                          'this is a single estimated M0 value, ', ...
+                          'e.g. when the M0 is obtained from an external scan and/or study'];
+        m0_value = perf.meta.M0;
+
+      case 'Absent'
+        m0_type = 'use_control_as_m0';
+        m0_explanation = [ ...
+                          'M0 is absent, so we can use the (average) control volume ', ...
+                          'as pseudo-M0 (if no background suppression was used)'];
+
+        if perf.meta.BackgroundSuppression == true
+          warning('Caution when using control as M0, background suppression was applied');
+        end
 
-end
+      otherwise
+        warning(['Unknown M0Type:', perf.meta.M0Type]);
 
-function file_list = return_session_specific_file_list(modality, subject)
+    end
 
-  switch modality
+    if ~isempty(m0_type)
+      perf.dependencies.m0.type = m0_type;
+      perf.dependencies.m0.explanation = m0_explanation;
+    end
 
-    case {'eeg', 'meg'}
-      pattern = [ ...
-                 '(_ses-[a-zA-Z0-9]+)?.*_', ...
-                 '(electrodes\\.tsv|photo\\.jpg|coordsystem\\.json|headshape\\..*)'];
+    if ~isempty(m0_volume_index)
+      perf.dependencies.m0.volume_index = m0_volume_index;
+    end
 
-  end
+    if ~isempty(m0_value)
+      perf.dependencies.m0.value = m0_value;
+    end
 
-  pth = fullfile(subject.path, modality);
+    if ~isempty(m0_filename)
+      perf.dependencies.m0.filename = m0_filename;
+    end
 
-  [file_list, d] = bids.internal.file_utils('List', ...
-                                            pth, ...
-                                            sprintf(['^%s' pattern '$'], ...
-                                                    subject.name));
+    if ~isempty(m0_sidecar)
+      perf.dependencies.m0.sidecar = m0_sidecar;
+    end
 
-  file_list = convert_to_cell(file_list);
+  end
 
 end
diff --git a/+bids/query.m b/+bids/query.m
index 80bbe10f..cbed700d 100644
--- a/+bids/query.m
+++ b/+bids/query.m
@@ -1,13 +1,25 @@
 function result = query(BIDS, query, varargin)
-  % QUERY Query a directory structure formatted according to the BIDS standard
-  % FORMAT result = bids.query(BIDS,query,...)
-  % BIDS   - BIDS directory name or BIDS structure (from bids.layout)
-  % query  - type of query: {'data', 'metadata', 'sessions', 'subjects',
-  %          'runs', 'tasks', 'runs', 'types', 'modalities'}
-  % result - outcome of query
   %
-  % See also:
-  % bids
+  % Queries a directory structure formatted according to the BIDS standard
+  %
+  % USAGE::
+  %
+  %   result = bids.query(BIDS, query, ...)
+  %
+  % :param BIDS: BIDS directory name or BIDS structure (from bids.layout)
+  % :type  BIDS: (strcuture or string)
+  % :param query: type of query:
+  %                          - 'data',
+  %                          - 'metadata',
+  %                          - 'sessions',
+  %                          - 'subjects',
+  %                          - 'runs',
+  %                          - 'tasks',
+  %                          - 'suffixes',
+  %                          - 'modalities'
+  % :type  query: string
+  %
+  %
 
   % __________________________________________________________________________
   %
@@ -23,161 +35,225 @@
   %#ok<*AGROW>
   narginchk(2, Inf);
 
+  VALID_QUERIES = { ...
+                   'sessions', ...
+                   'subjects', ...
+                   'modalities', ...
+                   'tasks', ...
+                   'runs', ...
+                   'suffixes', ...
+                   'data', ...
+                   'metadata', ...
+                   'dependencies'};
+
+  if ~any(strcmp(query, VALID_QUERIES))
+    error('Invalid query input: ''%s''', query);
+  end
+
   BIDS = bids.layout(BIDS);
 
-  opts = parse_query(varargin);
+  options = parse_query(varargin);
+
+  % For subjects and modality we pass only the subjects/modalities asked for
+  % otherwise we pass all of them
+
+  [subjects, options] = get_subjects(BIDS, options);
+
+  [modalities, options] = get_modalities(BIDS, options);
+
+  % Get optional target option for metadata query
+  [target, options] = get_target(query, options);
+
+  result = perform_query(BIDS, query, options, subjects, modalities, target);
 
+  %% Postprocessing output variable
   switch query
-    case {'sessions', 'subjects', 'modalities', 'tasks', 'runs', 'types', 'data', 'metadata'}
-      % -Initialise output variable
-      result = {};
-
-      % -For subjects and modality we pass only the subjects/modalities asked for
-      % otherwise we pass all of them
-
-      % -Filter according to subjects
-      if any(ismember(opts(:, 1), 'sub'))
-        subs = opts{ismember(opts(:, 1), 'sub'), 2};
-        opts(ismember(opts(:, 1), 'sub'), :) = [];
-      else
-        subs = unique({BIDS.subjects.name});
-        subs = regexprep(subs, '^[a-zA-Z0-9]+-', '');
-      end
-      % -Filter according to modality
-      if any(ismember(opts(:, 1), 'modality'))
-        mods = opts{ismember(opts(:, 1), 'modality'), 2};
-        opts(ismember(opts(:, 1), 'modality'), :) = [];
-      else
-        hasmod = arrayfun(@(y) structfun(@(x) isstruct(x) & ~isempty(x), y), ...
-                          BIDS.subjects, 'UniformOutput', false);
-        hasmod = any([hasmod{:}], 2);
-        mods   = fieldnames(BIDS.subjects)';
-        mods   = mods(hasmod);
-      end
 
-      % -Get optional target option for metadata query
-      if strcmp(query, 'metadata') && any(ismember(opts(:, 1), 'target'))
-        target = opts{ismember(opts(:, 1), 'target'), 2};
-        opts(ismember(opts(:, 1), 'target'), :) = [];
-        if iscellstr(target)
-          target = substruct('.', target{1});
-        end
-      else
-        target = [];
-      end
+    case 'subjects'
+      result = unique(result);
+      result = regexprep(result, '^[a-zA-Z0-9]+-', '');
 
-      % -Perform query
-      % Loop through all the subjects and modalities filtered previously
-      for i = 1:numel(BIDS.subjects)
-        % -Only continue if this subject is one of those filtered
-        if ~ismember(BIDS.subjects(i).name(5:end), subs)
-          continue
-        end
-        for j = 1:numel(mods)
-          d = BIDS.subjects(i).(mods{j});
-          for k = 1:numel(d)
-            % -sts is kept true only if this modality is one of those filtered
-            sts = true;
-            for l = 1:size(opts, 1)
-              if ~isfield(d(k), opts{l, 1}) || ~ismember(d(k).(opts{l, 1}), opts{l, 2})
-                sts = false;
-              end
-            end
-            switch query
-              case 'subjects'
-                if sts
-                  result{end + 1} = BIDS.subjects(i).name;
-                end
-              case 'sessions'
-                if sts
-                  result{end + 1} = BIDS.subjects(i).session;
-                end
-              case 'modalities'
-                if sts
-                  hasmod = structfun(@(x) isstruct(x) & ~isempty(x), ...
-                                     BIDS.subjects(i));
-                  allmods = fieldnames(BIDS.subjects(i))';
-                  result = union(result, allmods(hasmod));
-                end
-              case 'data'
-                if sts && isfield(d(k), 'filename')
-                  result{end + 1} = fullfile(BIDS.subjects(i).path, mods{j}, d(k).filename);
-                end
-              case 'metadata'
-                if sts && isfield(d(k), 'filename')
-                  f = fullfile(BIDS.subjects(i).path, mods{j}, d(k).filename);
-                  result{end + 1} = bids.internal.get_metadata(f);
-                  if ~isempty(target)
-                    try
-                      result{end} = subsref(result{end}, target);
-                    catch
-                      warning('Non-existent field for metadata.');
-                      result{end} = [];
-                    end
-                  end
-                end
-                %                             if sts && isfield(d(k),'meta')
-                %                                 result{end+1} = d(k).meta;
-                %                             end
-              case 'runs'
-                if sts && isfield(d(k), 'run')
-                  result{end + 1} = d(k).run;
-                end
-              case 'tasks'
-                if sts && isfield(d(k), 'task')
-                  result{end + 1} = d(k).task;
-                end
-              case 'types'
-                if sts && isfield(d(k), 'type')
-                  result{end + 1} = d(k).type;
-                end
-            end
-          end
-        end
-      end
+    case 'sessions'
+      result = unique(result);
+      result = regexprep(result, '^[a-zA-Z0-9]+-', '');
+      result(cellfun('isempty', result)) = [];
 
-      % -Postprocessing output variable
-      switch query
-        case 'subjects'
-          result = unique(result);
-          result = regexprep(result, '^[a-zA-Z0-9]+-', '');
-        case 'sessions'
-          result = unique(result);
-          result = regexprep(result, '^[a-zA-Z0-9]+-', '');
-          result(cellfun('isempty', result)) = [];
-        case {'modalities', 'data'}
-          result = result';
-        case 'metadata'
-          if numel(result) == 1
-            result = result{1};
-          end
-        case {'tasks', 'runs', 'types'}
-          result = unique(result);
-          result(cellfun('isempty', result)) = [];
+    case {'modalities', 'data'}
+      result = result';
+
+    case 'metadata'
+      if numel(result) == 1
+        result = result{1};
       end
-    otherwise
-      error('Invalid query input: ''%s''', query);
+
+    case {'tasks', 'runs', 'suffixes'}
+      result = unique(result);
+      result(cellfun('isempty', result)) = [];
   end
 
-  % ==========================================================================
-  % -Parse BIDS query
-  % ==========================================================================
-function query = parse_query(query)
-  if numel(query) == 1 && isstruct(query{1})
-    query = [fieldnames(query{1}), struct2cell(query{1})];
+end
+
+function options = parse_query(options)
+
+  if numel(options) == 1 && isstruct(options{1})
+    options = [fieldnames(options{1}), struct2cell(options{1})];
+
   else
-    if mod(numel(query), 2)
+    if mod(numel(options), 2)
       error('Invalid input syntax: each BIDS entity requires an associated label');
     end
-    query = reshape(query, 2, [])';
+    options = reshape(options, 2, [])';
+
+  end
+
+  for i = 1:size(options, 1)
+
+    if ischar(options{i, 2})
+      options{i, 2} = cellstr(options{i, 2});
+    end
+
+    for j = 1:numel(options{i, 2})
+      if iscellstr(options{i, 2})
+        options{i, 2}{j} = regexprep(options{i, 2}{j}, sprintf('^%s-', options{i, 1}), '');
+      end
+    end
+
+  end
+
+end
+
+function [subjects, options] = get_subjects(BIDS, options)
+
+  if any(ismember(options(:, 1), 'sub'))
+    subjects = options{ismember(options(:, 1), 'sub'), 2};
+    options(ismember(options(:, 1), 'sub'), :) = [];
+  else
+    subjects = unique({BIDS.subjects.name});
+    subjects = regexprep(subjects, '^[a-zA-Z0-9]+-', '');
+  end
+
+end
+
+function [modalities, options] = get_modalities(BIDS, options)
+
+  if any(ismember(options(:, 1), 'modality'))
+    modalities = options{ismember(options(:, 1), 'modality'), 2};
+    options(ismember(options(:, 1), 'modality'), :) = [];
+  else
+    hasmod = arrayfun(@(y) structfun(@(x) isstruct(x) & ~isempty(x), y), ...
+                      BIDS.subjects, 'UniformOutput', false);
+    hasmod = any([hasmod{:}], 2);
+    modalities   = fieldnames(BIDS.subjects)';
+    modalities   = modalities(hasmod);
+  end
+
+end
+
+function [target, options] = get_target(query, options)
+
+  target = [];
+
+  if strcmp(query, 'metadata') && any(ismember(options(:, 1), 'target'))
+
+    target = options{ismember(options(:, 1), 'target'), 2};
+    options(ismember(options(:, 1), 'target'), :) = [];
+
+    if iscellstr(target)
+      target = substruct('.', target{1});
+    end
+
   end
-  for i = 1:size(query, 1)
-    if ischar(query{i, 2})
-      query{i, 2} = cellstr(query{i, 2});
+
+end
+
+function result = perform_query(BIDS, query, options, subjects, modalities, target)
+
+  % Initialise output variable
+  result = {};
+
+  % Loop through all the subjects and modalities filtered previously
+  for i = 1:numel(BIDS.subjects)
+
+    % Only continue if this subject is one of those filtered
+    if ~ismember(BIDS.subjects(i).name(5:end), subjects)
+      continue
     end
-    for j = 1:numel(query{i, 2})
-      if iscellstr(query{i, 2})
-        query{i, 2}{j} = regexprep(query{i, 2}{j}, sprintf('^%s-', query{i, 1}), '');
+
+    for j = 1:numel(modalities)
+
+      d = BIDS.subjects(i).(modalities{j});
+
+      for k = 1:numel(d)
+
+        % status is kept true only if this file matches
+        % the options of the query
+        status = bids.internal.keep_file_for_query(d(k), options);
+
+        if status
+
+          switch query
+
+            case 'subjects'
+              result{end + 1} = BIDS.subjects(i).name;
+
+            case 'sessions'
+              result{end + 1} = BIDS.subjects(i).session;
+
+            case 'modalities'
+              hasmod = structfun(@(x) isstruct(x) & ~isempty(x), ...
+                                 BIDS.subjects(i));
+              allmods = fieldnames(BIDS.subjects(i))';
+              result = union(result, allmods(hasmod));
+
+            case 'data'
+              if isfield(d(k), 'filename')
+                result{end + 1} = fullfile(BIDS.subjects(i).path, modalities{j}, d(k).filename);
+              end
+
+            case 'metadata'
+              if isfield(d(k), 'filename')
+
+                f = fullfile(BIDS.subjects(i).path, modalities{j}, d(k).filename);
+                result{end + 1} = bids.internal.get_metadata(f);
+                if ~isempty(target)
+                  try
+                    result{end} = subsref(result{end}, target);
+                  catch
+                    warning('Non-existent field for metadata.');
+                    result{end} = [];
+                  end
+                end
+
+              end
+              % if status && isfield(d(k),'meta')
+              %   result{end+1} = d(k).meta;
+              % end
+
+            case 'runs'
+              if isfield(d(k).entities, 'run')
+                result{end + 1} = d(k).entities.run;
+              end
+
+            case 'tasks'
+              if isfield(d(k).entities, 'task')
+                result{end + 1} = d(k).entities.task;
+              end
+
+            case 'suffixes'
+              if isfield(d(k), 'suffix')
+                result{end + 1} = d(k).suffix;
+              end
+
+            case 'dependencies'
+              if isfield(d(k), 'dependencies')
+                result{end + 1} = d(k).dependencies;
+              end
+
+          end
+
+        end
       end
     end
   end
+
+end
diff --git a/.github/workflows/run_tests.yml b/.github/workflows/run_tests.yml
index c64684d5..88e50fee 100644
--- a/.github/workflows/run_tests.yml
+++ b/.github/workflows/run_tests.yml
@@ -29,7 +29,7 @@ jobs:
       uses: actions/checkout@v2
       with:
         submodules: true
-        fetch-depth: 1
+        fetch-depth: 2
 
     - name: Install BIDS validator
       run: |        
@@ -58,6 +58,12 @@ jobs:
         with_coverage: true
         cover_xml_file: coverage.xml
 
+    - name: Upload coverage
+      uses: actions/upload-artifact@v1 
+      with: 
+        name: coverage_file 
+        path: coverage.xml         
+
     - name: Code coverage
       uses: codecov/codecov-action@v1
       with:
diff --git a/tests/schema/base.json b/tests/schema/base.json
new file mode 100644
index 00000000..73cd85a8
--- /dev/null
+++ b/tests/schema/base.json
@@ -0,0 +1,3 @@
+{
+    "test": 1
+}
\ No newline at end of file
diff --git a/tests/schema/subfolder_1/sub.json b/tests/schema/subfolder_1/sub.json
new file mode 100644
index 00000000..73cd85a8
--- /dev/null
+++ b/tests/schema/subfolder_1/sub.json
@@ -0,0 +1,3 @@
+{
+    "test": 1
+}
\ No newline at end of file
diff --git a/tests/schema/subfolder_2/sub.txt b/tests/schema/subfolder_2/sub.txt
new file mode 100644
index 00000000..e69de29b
diff --git a/tests/schema/subfolder_2/subfolder_3/sub.json b/tests/schema/subfolder_2/subfolder_3/sub.json
new file mode 100644
index 00000000..73cd85a8
--- /dev/null
+++ b/tests/schema/subfolder_2/subfolder_3/sub.json
@@ -0,0 +1,3 @@
+{
+    "test": 1
+}
\ No newline at end of file
diff --git a/tests/schema/subfolder_4/sub.txt b/tests/schema/subfolder_4/sub.txt
new file mode 100644
index 00000000..e69de29b
diff --git a/tests/test_append_to_layout.m b/tests/test_append_to_layout.m
new file mode 100644
index 00000000..492b81e3
--- /dev/null
+++ b/tests/test_append_to_layout.m
@@ -0,0 +1,95 @@
+function test_suite = test_append_to_layout %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_append_to_layout_basic()
+
+  schema = bids.schema.load_schema();
+
+  subject = struct('anat', struct([]));
+
+  modality = 'anat';
+
+  file = '../sub-16/anat/sub-16_ses-mri_run-1_acq-hd_T1w.nii.gz';
+  subject = bids.internal.append_to_layout(file, subject, modality, schema);
+
+  expected.anat = struct( ...
+                         'filename', 'sub-16_ses-mri_run-1_acq-hd_T1w.nii.gz', ...
+                         'suffix', 'T1w', ...
+                         'ext', '.nii.gz', ...
+                         'entities', struct('sub', '16', ...
+                                            'ses', 'mri', ...
+                                            'run', '1', ...
+                                            'acq', 'hd', ...
+                                            'ce', '', ...
+                                            'rec', '', ...
+                                            'part', ''));
+
+  assertEqual(subject, expected);
+
+end
+
+function test_append_to_structure_basic_test()
+
+  schema = bids.schema.load_schema();
+
+  subject = struct('anat', struct([]));
+  modality = 'anat';
+
+  file = '../sub-16/anat/sub-16_ses-mri_run-1_acq-hd_T1w.nii.gz';
+  subject = bids.internal.append_to_layout(file, subject, modality, schema);
+
+  file = '../sub-16/anat/sub-16_ses-mri_run-1_T1map.nii.gz';
+  subject = bids.internal.append_to_layout(file, subject, modality, schema);
+
+  expected.anat(1, 1) = struct( ...
+                               'filename', 'sub-16_ses-mri_run-1_acq-hd_T1w.nii.gz', ...
+                               'suffix', 'T1w', ...
+                               'ext', '.nii.gz', ...
+                               'entities', struct('sub', '16', ...
+                                                  'ses', 'mri', ...
+                                                  'run', '1', ...
+                                                  'acq', 'hd', ...
+                                                  'ce', '', ...
+                                                  'rec', '', ...
+                                                  'part', ''));
+
+  expected.anat(2, 1) = struct( ...
+                               'filename', 'sub-16_ses-mri_run-1_T1map.nii.gz', ...
+                               'suffix', 'T1map', ...
+                               'ext', '.nii.gz', ...
+                               'entities', struct('sub', '16', ...
+                                                  'ses', 'mri', ...
+                                                  'run', '1', ...
+                                                  'acq', '', ...
+                                                  'ce', '', ...
+                                                  'rec', ''));     %#ok<*STRNU>
+
+  assertEqual(subject, expected);
+
+end
+
+function test_append_to_layout_schemaless()
+
+  subject = struct('newmod', struct([]));
+
+  modality = 'newmod';
+
+  file = '../sub-16/newmod/sub-16_schema-less_anything-goes_newsuffix.EXT';
+  subject = bids.internal.append_to_layout(file, subject, modality);
+
+  expected.newmod = struct( ...
+                           'filename', 'sub-16_schema-less_anything-goes_newsuffix.EXT', ...
+                           'suffix', 'newsuffix', ...
+                           'ext', '.EXT', ...
+                           'entities', struct('sub', '16', ...
+                                              'schema', 'less', ...
+                                              'anything', 'goes'));
+
+  assertEqual(subject.newmod, expected.newmod);
+
+end
diff --git a/tests/test_bids_examples.m b/tests/test_bids_examples.m
index 1cd30282..8db00981 100644
--- a/tests/test_bids_examples.m
+++ b/tests/test_bids_examples.m
@@ -30,17 +30,17 @@ function test_bids_examples_basic()
 
   % -Try to run bids.layout on each dataset directory and keep track of any
   % failure with a try/catch
-  sts = false(1, numel(d));
+  status = false(1, numel(d));
   msg = cell(1, numel(d));
   for i = 1:numel(d)
     if exist(fullfile(pth_bids_example, d(i).name, '.SKIP_VALIDATION'), 'file')
-      sts(i) = true;
+      status(i) = true;
       fprintf('-');
       continue
     end
     try
       BIDS = bids.layout(fullfile(pth_bids_example, d(i).name));
-      sts(i) = true;
+      status(i) = true;
       fprintf('.');
     catch err
       fprintf('X');
@@ -50,8 +50,8 @@ function test_bids_examples_basic()
   fprintf('\n');
 
   % lists all the folder for which bids.layout failed
-  if ~all(sts)
-    for i = find(~sts)
+  if ~all(status)
+    for i = find(~status)
       fprintf('* %s: %s\n', d(i).name, msg{i});
     end
     error('Parsing of BIDS-compatible datasets failed.');
diff --git a/tests/test_bids_query.m b/tests/test_bids_query.m
index 211c04d2..07e9f287 100644
--- a/tests/test_bids_query.m
+++ b/tests/test_bids_query.m
@@ -4,10 +4,8 @@
   catch % no problem; early Matlab versions can use initTestSuite fine
   end
   initTestSuite;
-end
 
-function test_bids_query_basic()
-  % Test BIDS queries on ds007
+  % Test BIDS queries
   % This dataset comes from https://github.com/bids-standard/bids-examples
   % and is downloaded automatically by the continuous integration framework
   % and is required for the tests to be run.
@@ -23,91 +21,66 @@ function test_bids_query_basic()
   % Copyright (C) 2019, Guillaume Flandin, Wellcome Centre for Human Neuroimaging
   % Copyright (C) 2019--, BIDS-MATLAB developers
 
-  pth_bids_example = get_test_data_dir();
-
-  BIDS = bids.layout(fullfile(pth_bids_example, 'ds007'));
-
-  subjs = arrayfun(@(x) sprintf('%02d', x), 1:20, 'UniformOutput', false);
-  assert(isequal(bids.query(BIDS, 'subjects'), subjs));
-
-  assert(isempty(bids.query(BIDS, 'sessions')));
+end
 
-  assert(isequal(bids.query(BIDS, 'runs'), {'01', '02'}));
+function test_bids_query_data()
 
-  tasks = { ...
-           'stopsignalwithletternaming', ...
-           'stopsignalwithmanualresponse', ...
-           'stopsignalwithpseudowordnaming'};
-  assert(isequal(bids.query(BIDS, 'tasks'), tasks));
-
-  types = {'T1w', 'bold', 'events', 'inplaneT2'};
-  assert(isequal(bids.query(BIDS, 'types'), types));
+  pth_bids_example = get_test_data_dir();
 
-  mods = {'anat', 'func'};
-  assert(isequal(bids.query(BIDS, 'modalities'), mods));
-  assert(isequal(bids.query(BIDS, 'modalities', 'sub', '01'), mods));
+  BIDS = bids.layout(fullfile(pth_bids_example, 'ds007'));
 
-  assert(isempty(bids.query(BIDS, 'runs', 'type', 'T1w')));
+  t1 = bids.query(BIDS, 'data', 'suffix', 'T1w');
+  assert(iscellstr(t1));
+  assert(numel(t1) == numel(bids.query(BIDS, 'subjects')));
 
-  runs = {'01', '02'};
-  assert(isequal(bids.query(BIDS, 'runs', 'type', 'bold'), runs));
+  data = bids.query(BIDS, 'data', 'sub', '01', 'task', 'stopsignalwithpseudowordnaming');
+  assertEqual(size(data, 1), 4);
 
   bold = bids.query(BIDS, 'data', ...
                     'sub', '05', ...
                     'run', '02', ...
                     'task', 'stopsignalwithmanualresponse', ...
-                    'type', 'bold');
+                    'suffix', 'bold');
   assert(iscellstr(bold));
   assert(numel(bold) == 1);
 
+end
+
+function test_bids_query_metadata()
+
+  pth_bids_example = get_test_data_dir();
+
+  BIDS = bids.layout(fullfile(pth_bids_example, 'ds007'));
+
   md = bids.query(BIDS, 'metadata', ...
                   'sub', '05', ...
                   'run', '02', ...
                   'task', 'stopsignalwithmanualresponse', ...
-                  'type', 'bold');
+                  'suffix', 'bold');
   assert(isstruct(md) & isfield(md, 'RepetitionTime') & isfield(md, 'TaskName'));
   assert(md.RepetitionTime == 2);
   assert(strcmp(md.TaskName, 'stop signal with manual response'));
 
-  t1 = bids.query(BIDS, 'data', 'type', 'T1w');
-  assert(iscellstr(t1));
-  assert(numel(t1) == numel(bids.query(BIDS, 'subjects')));
-
 end
 
-function test_bids_query_sessions()
-  %
-  %   parse a folder with sessions
-  %
+function test_bids_query_modalities()
 
   pth_bids_example = get_test_data_dir();
 
-  BIDS = bids.layout(fullfile(pth_bids_example, 'synthetic'));
-
-  %   test
-  sessions = {'01', '02'};
-  assert(isequal(bids.query(BIDS, 'sessions'), sessions));
-  assert(isequal(bids.query(BIDS, 'sessions', 'sub', '02'), sessions));
-
-end
-
-function test_bids_query_modalities()
-  %
-  %   parse a folder with different modalities per session
-  %
+  BIDS = bids.layout(fullfile(pth_bids_example, 'ds007'));
 
-  pth_bids_example = get_test_data_dir();
+  modalities = {'anat', 'func'};
+  assertEqual(bids.query(BIDS, 'modalities'), modalities);
+  assertEqual(bids.query(BIDS, 'modalities', 'sub', '01'), modalities);
 
   BIDS = bids.layout(fullfile(pth_bids_example, '7t_trt'));
 
-  %   test
-  mods = {'anat', 'fmap', 'func'};
+  modalities = {'anat', 'fmap', 'func'};
 
-  assert(isequal(bids.query(BIDS, 'modalities'), mods));
-  assert(isequal(bids.query(BIDS, 'modalities', 'sub', '01'), mods));
-  assert(isequal(bids.query(BIDS, 'modalities', 'sub', '01', 'ses', '1'), mods));
+  assertEqual(bids.query(BIDS, 'modalities'), modalities);
+  assertEqual(bids.query(BIDS, 'modalities', 'sub', '01'), modalities);
+  assertEqual(bids.query(BIDS, 'modalities', 'sub', '01', 'ses', '1'), modalities);
 
-  %
   % this now fails on octave 4.2.2 but not on Matlab
   %
   % bids.query(BIDS, 'modalities', 'sub', '01', 'ses', '2')
@@ -121,6 +94,63 @@ function test_bids_query_modalities()
   %
   % when it should return
 
-  % assert(isequal(bids.query(BIDS, 'modalities', 'sub', '01', 'ses', '2'), mods(2:3)));
+  % assertEqual(bids.query(BIDS, 'modalities', 'sub', '01', 'ses', '2'), mods(2:3)));
+
+end
+
+function test_bids_query_basic()
+
+  pth_bids_example = get_test_data_dir();
+
+  BIDS = bids.layout(fullfile(pth_bids_example, 'ds007'));
+
+  tasks = { ...
+           'stopsignalwithletternaming', ...
+           'stopsignalwithmanualresponse', ...
+           'stopsignalwithpseudowordnaming'};
+  assertEqual(bids.query(BIDS, 'tasks'), tasks);
+
+  assert(isempty(bids.query(BIDS, 'runs', 'suffix', 'T1w')));
+
+  runs = {'01', '02'};
+  assertEqual(bids.query(BIDS, 'runs'), runs);
+  assertEqual(bids.query(BIDS, 'runs', 'suffix', 'bold'), runs);
+
+end
+
+function test_bids_query_subjects()
+
+  pth_bids_example = get_test_data_dir();
+
+  BIDS = bids.layout(fullfile(pth_bids_example, 'ds007'));
+
+  subjs = arrayfun(@(x) sprintf('%02d', x), 1:20, 'UniformOutput', false);
+  assertEqual(bids.query(BIDS, 'subjects'), subjs);
+
+end
+
+function test_bids_query_sessions()
+
+  pth_bids_example = get_test_data_dir();
+
+  BIDS = bids.layout(fullfile(pth_bids_example, 'synthetic'));
+  sessions = {'01', '02'};
+  assertEqual(bids.query(BIDS, 'sessions'), sessions);
+  assertEqual(bids.query(BIDS, 'sessions', 'sub', '02'), sessions);
+
+  BIDS = bids.layout(fullfile(pth_bids_example, 'ds007'));
+
+  assert(isempty(bids.query(BIDS, 'sessions')));
+
+end
+
+function test_bids_query_suffixes()
+
+  pth_bids_example = get_test_data_dir();
+
+  BIDS = bids.layout(fullfile(pth_bids_example, 'ds007'));
+
+  suffixes = {'T1w', 'bold', 'events', 'inplaneT2'};
+  assertEqual(bids.query(BIDS, 'suffixes'), suffixes);
 
 end
diff --git a/tests/test_bids_query_asl.m b/tests/test_bids_query_asl.m
index 4211f0ca..a8740b77 100644
--- a/tests/test_bids_query_asl.m
+++ b/tests/test_bids_query_asl.m
@@ -8,32 +8,45 @@
 
 function test_bids_query_asl_basic()
   %
-  %   asl specific queries
+  %   asl queries
   %
 
   pth_bids_example = get_test_data_dir();
 
   %% 'asl001'
-  BIDS = bids.layout(fullfile(pth_bids_example, 'asl001 '));
+  BIDS = bids.layout(fullfile(pth_bids_example, 'asl001'));
 
   modalities = {'anat', 'perf'};
   assertEqual(bids.query(BIDS, 'modalities'), modalities);
 
-  types = {'T1w', 'asl'};
-  %   types = {'T1w', 'asl', 'aslcontext', 'asllabelling'};
-  assertEqual(bids.query(BIDS, 'types'), types);
+  suffixes = {'T1w', 'asl', 'aslcontext', 'asllabeling'};
+  assertEqual(bids.query(BIDS, 'suffixes'), suffixes);
+
+  BIDS.subjects(1).perf(1);
+
+  filename = bids.query(BIDS, 'data', 'sub', 'Sub103', 'suffix', 'asl');
+  basename = bids.internal.file_utils(filename, 'basename');
+  assertEqual(basename, {'sub-Sub103_asl.nii'});
+
+  meta = bids.query(BIDS, 'metadata', 'sub', 'Sub103', 'suffix', 'asl');
+
+  dependencies = bids.query(BIDS, 'dependencies', 'sub', 'Sub103', 'suffix', 'asl');
+  dependencies{1}.labeling_image;
+  dependencies{1}.context;
+  dependencies{1}.m0;
 
   %% 'asl002'
-  BIDS = bids.layout(fullfile(pth_bids_example, 'asl002 '));
+  BIDS = bids.layout(fullfile(pth_bids_example, 'asl002'));
 
   modalities = {'anat', 'perf'};
   assertEqual(bids.query(BIDS, 'modalities'), modalities);
 
-  types = {'T1w', 'asl', 'm0scan'};
-  %   types = {'T1w', 'asl', 'aslcontext', 'asllabelling', 'm0scan'};
-  assertEqual(bids.query(BIDS, 'types'), types);
-  assertEqual(bids.internal.file_utils(bids.query(BIDS, 'data', 'type', 'm0scan'), 'basename'), ...
-              {'sub-Sub103_m0scan.nii'});
+  suffixes = {'T1w', 'asl', 'aslcontext', 'asllabeling', 'm0scan'};
+  assertEqual(bids.query(BIDS, 'suffixes'), suffixes);
+
+  filename = bids.query(BIDS, 'data', 'sub', 'Sub103', 'suffix', 'm0scan');
+  basename = bids.internal.file_utils(filename, 'basename');
+  assertEqual(basename, {'sub-Sub103_m0scan.nii'});
 
   %% 'asl003'
   BIDS = bids.layout(fullfile(pth_bids_example, 'asl003'));
@@ -41,9 +54,8 @@ function test_bids_query_asl_basic()
   modalities = {'anat', 'perf'};
   assertEqual(bids.query(BIDS, 'modalities'), modalities);
 
-  types = {'T1w', 'asl', 'm0scan'};
-  %   types = {'T1w', 'asl', 'aslcontext', 'asllabelling', 'm0scan'};
-  assertEqual(bids.query(BIDS, 'types'), types);
+  suffixes = {'T1w', 'asl', 'aslcontext', 'asllabeling', 'm0scan'};
+  assertEqual(bids.query(BIDS, 'suffixes'), suffixes);
 
   %% 'asl004'
   BIDS = bids.layout(fullfile(pth_bids_example, 'asl004'));
@@ -51,8 +63,11 @@ function test_bids_query_asl_basic()
   modalities = {'anat', 'fmap', 'perf'};
   assertEqual(bids.query(BIDS, 'modalities'), modalities);
 
-  types = {'T1w', 'asl', 'm0scan'};
-  %   types = {'T1w', 'asl', 'aslcontext', 'asllabelling', 'm0scan'};
-  assertEqual(bids.query(BIDS, 'types'), types);
+  suffixes = {'T1w', 'asl', 'aslcontext', 'asllabeling', 'm0scan'};
+  assertEqual(bids.query(BIDS, 'suffixes'), suffixes);
+
+  filename = bids.query(BIDS, 'data', 'suffix', 'm0scan', 'dir', 'pa');
+  basename = bids.internal.file_utils(filename, 'basename');
+  assertEqual(basename, {'sub-Sub1_dir-pa_m0scan.nii'});
 
 end
diff --git a/tests/test_bids_query_dwi.m b/tests/test_bids_query_dwi.m
new file mode 100644
index 00000000..59475244
--- /dev/null
+++ b/tests/test_bids_query_dwi.m
@@ -0,0 +1,27 @@
+function test_suite = test_bids_query_dwi %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_bids_query_dwi_basic()
+  %
+  %   eeg queries
+  %
+
+  pth_bids_example = get_test_data_dir();
+
+  %%
+  BIDS = bids.layout(fullfile(pth_bids_example, 'ds000117'));
+
+  modalities = {'anat',    'beh',    'dwi',    'fmap',    'func',    'meg'};
+  assertEqual(bids.query(BIDS, 'modalities'), modalities);
+
+  suffixes = {'T1w', 'bold', 'dwi', 'events', 'headshape', ...
+              'magnitude1', 'magnitude2', 'meg', 'phasediff'};
+  % Missing: 'FLASH'
+  assertEqual(bids.query(BIDS, 'suffixes'), suffixes);
+
+end
diff --git a/tests/test_bids_query_eeg.m b/tests/test_bids_query_eeg.m
new file mode 100644
index 00000000..7421638c
--- /dev/null
+++ b/tests/test_bids_query_eeg.m
@@ -0,0 +1,65 @@
+function test_suite = test_bids_query_eeg %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_bids_query_eeg_basic()
+  %
+  %   eeg queries
+  %
+
+  pth_bids_example = get_test_data_dir();
+
+  %%
+  BIDS = bids.layout(fullfile(pth_bids_example, 'eeg_cbm'));
+
+  modalities = {'eeg'};
+  assertEqual(bids.query(BIDS, 'modalities'), modalities);
+
+  suffixes = {'channels', 'eeg', 'events'};
+  assertEqual(bids.query(BIDS, 'suffixes'), suffixes);
+
+  %%
+  BIDS = bids.layout(fullfile(pth_bids_example, 'eeg_ds000117'));
+
+  modalities = {'anat', 'eeg'};
+  assertEqual(bids.query(BIDS, 'modalities'), modalities);
+
+  suffixes = {'T1w', 'eeg', 'electrodes', 'events'};
+  % Missing: 'coordsystem',
+  % Missing: 'channels' in root folder
+  assertEqual(bids.query(BIDS, 'suffixes'), suffixes);
+
+  %%
+  BIDS = bids.layout(fullfile(pth_bids_example, 'eeg_face13'));
+
+  modalities = {'eeg'};
+  assertEqual(bids.query(BIDS, 'modalities'), modalities);
+
+  suffixes = {'channels', 'eeg', 'events'};
+  % Missing: 'electrodes'
+  % skipped as it contains a task entity and thus does not match the schema
+  assertEqual(bids.query(BIDS, 'suffixes'), suffixes);
+
+  %%
+  BIDS = bids.layout(fullfile(pth_bids_example, 'eeg_matchingpennies'));
+
+  modalities = {'eeg'};
+  assertEqual(bids.query(BIDS, 'modalities'), modalities);
+
+  suffixes = {'channels', 'eeg', 'events'};
+  assertEqual(bids.query(BIDS, 'suffixes'), suffixes);
+
+  %%
+  BIDS = bids.layout(fullfile(pth_bids_example, 'eeg_rishikesh'));
+
+  modalities = {'eeg'};
+  assertEqual(bids.query(BIDS, 'modalities'), modalities);
+
+  suffixes = {'channels', 'eeg', 'events'};
+  assertEqual(bids.query(BIDS, 'suffixes'), suffixes);
+
+end
diff --git a/tests/test_bids_query_ieeg.m b/tests/test_bids_query_ieeg.m
new file mode 100644
index 00000000..89c22396
--- /dev/null
+++ b/tests/test_bids_query_ieeg.m
@@ -0,0 +1,36 @@
+function test_suite = test_bids_query_ieeg %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_bids_query_ieeg_basic()
+  %
+  %   eeg queries
+  %
+
+  pth_bids_example = get_test_data_dir();
+
+  %%
+  BIDS = bids.layout(fullfile(pth_bids_example, 'ieeg_epilepsy'));
+
+  modalities = {'anat', 'ieeg'};
+  assertEqual(bids.query(BIDS, 'modalities'), modalities);
+
+  suffixes = {'T1w', 'channels', 'electrodes', 'events', 'ieeg'};
+  % Missing: 'coordsystem'
+  assertEqual(bids.query(BIDS, 'suffixes'), suffixes);
+
+  %%
+  BIDS = bids.layout(fullfile(pth_bids_example, 'ieeg_epilepsy_ecog'));
+
+  modalities = {'anat', 'ieeg'};
+  assertEqual(bids.query(BIDS, 'modalities'), modalities);
+
+  suffixes = {'T1w', 'channels', 'electrodes', 'events', 'ieeg', 'photo'};
+  % Missing: 'coordsystem'
+  assertEqual(bids.query(BIDS, 'suffixes'), suffixes);
+
+end
diff --git a/tests/test_bids_query_meg.m b/tests/test_bids_query_meg.m
new file mode 100644
index 00000000..55c27b32
--- /dev/null
+++ b/tests/test_bids_query_meg.m
@@ -0,0 +1,28 @@
+function test_suite = test_bids_query_meg %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_bids_query_meg_basic()
+  %
+  %   meg queries
+  %
+
+  pth_bids_example = get_test_data_dir();
+
+  BIDS = bids.layout(fullfile(pth_bids_example, 'ds000246'));
+
+  modalities = {'anat', 'meg'};
+  assertEqual(bids.query(BIDS, 'modalities'), modalities);
+
+  suffixes = {'T1w', 'channels', 'headshape', 'meg', 'photo'};
+  % missing: 'coordsystem'
+  assertEqual(bids.query(BIDS, 'suffixes'), suffixes);
+
+  BIDS = bids.layout(fullfile(pth_bids_example, 'ds000247'));
+  BIDS = bids.layout(fullfile(pth_bids_example, 'ds000248'));
+
+end
diff --git a/tests/test_file_utils.m b/tests/test_file_utils.m
index 6279b2fb..f7607873 100644
--- a/tests/test_file_utils.m
+++ b/tests/test_file_utils.m
@@ -1,4 +1,12 @@
-function test_file_utils()
+function test_suite = test_file_utils %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_file_utils_basic()
 
   %% test to get certain part of a filename
   % {'path', 'basename', 'ext', 'filename', 'cpath', 'fpath'}
@@ -6,21 +14,21 @@ function test_file_utils()
   str = fullfile('folder', 'filename.extension');
 
   p = bids.internal.file_utils(str, 'path');
-  assert(isequal(p, 'folder'));
+  assertEqual(p, 'folder');
 
   p = bids.internal.file_utils(str, 'fpath');
-  assert(isequal(p, fullfile(pwd, 'folder')));
+  assertEqual(p, fullfile(pwd, 'folder'));
 
   filename = bids.internal.file_utils(str, 'basename');
-  assert(isequal(filename, 'filename'));
+  assertEqual(filename, 'filename');
 
   ext = bids.internal.file_utils(str, 'ext');
-  assert(isequal(ext, 'extension'));
+  assertEqual(ext, 'extension');
 
   str = fullfile('folder', 'subfolder', '..', 'filename.extension');
   cpath = bids.internal.file_utils(str, 'cpath');
-  assert(isequal(cpath, ...
-                 fullfile(pwd, 'folder', 'filename.extension')));
+  assertEqual(cpath, ...
+              fullfile(pwd, 'folder', 'filename.extension'));
 
   %% test to set certain part of a filename
   % {'path', 'basename', 'ext', 'filename', 'prefix', 'suffix'}
@@ -28,27 +36,27 @@ function test_file_utils()
   str = fullfile('folder', 'filename.extension');
 
   new_str = bids.internal.file_utils(str, 'ext', 'newext');
-  assert(isequal(new_str, fullfile('folder', 'filename.newext')));
+  assertEqual(new_str, fullfile('folder', 'filename.newext'));
 
   new_str = bids.internal.file_utils(str, 'basename', 'new_name');
-  assert(isequal(new_str, fullfile('folder', 'new_name.extension')));
+  assertEqual(new_str, fullfile('folder', 'new_name.extension'));
 
   new_str = bids.internal.file_utils(str, 'filename', 'new_name.newext');
-  assert(isequal(new_str, fullfile('folder', 'new_name.newext')));
+  assertEqual(new_str, fullfile('folder', 'new_name.newext'));
 
   new_str = bids.internal.file_utils(str, 'prefix', 'pre_');
-  assert(isequal(new_str, fullfile('folder', 'pre_filename.extension')));
+  assertEqual(new_str, fullfile('folder', 'pre_filename.extension'));
 
   new_str = bids.internal.file_utils(str, 'suffix', '_suffix');
-  assert(isequal(new_str, fullfile('folder', 'filename_suffix.extension')));
+  assertEqual(new_str, fullfile('folder', 'filename_suffix.extension'));
 
   new_str = bids.internal.file_utils(str, 'path', fullfile(pwd, 'new_folder'));
-  assert(isequal(new_str, fullfile(pwd, 'new_folder', 'filename.extension')));
+  assertEqual(new_str, fullfile(pwd, 'new_folder', 'filename.extension'));
 
   new_str = bids.internal.file_utils(str, ...
                                      'prefix', 'pre_', ...
                                      'suffix', '_suffix');
-  assert(isequal(new_str, fullfile('folder', 'pre_filename_suffix.extension')));
+  assertEqual(new_str, fullfile('folder', 'pre_filename_suffix.extension'));
 
   %% test to list files
 
@@ -57,24 +65,29 @@ function test_file_utils()
   file = bids.internal.file_utils('List', ...
                                   test_directory, ...
                                   '^test_file_utils.m$');
-  assert(isequal(file, 'test_file_utils.m'));
+  assertEqual(file, 'test_file_utils.m');
+
+  file = bids.internal.file_utils('List', ...
+                                  test_directory, ...
+                                  '^.*.md$');
+  assertEqual(file, 'README.md');
 
   directory = bids.internal.file_utils('List', ...
                                        test_directory, ...
                                        'dir', ...
                                        '^data$');
-  assert(isequal(directory, 'data'));
+  assertEqual(directory, 'data');
 
   fp_file = bids.internal.file_utils('FPList', ...
                                      test_directory, ...
                                      '^test_file_utils.m$');
-  assert(isequal(fp_file, [mfilename('fullpath') '.m']));
+  assertEqual(fp_file, [mfilename('fullpath') '.m']);
 
   fp_directory = bids.internal.file_utils('FPList', ...
                                           test_directory, ...
                                           'dir', ...
                                           '^data$');
-  assert(isequal(fp_directory, ...
-                 fullfile(test_directory, 'data')));
+  assertEqual(fp_directory, ...
+              fullfile(test_directory, 'data'));
 
 end
diff --git a/tests/test_functionSignatures.m b/tests/test_function_signatures.m
similarity index 61%
rename from tests/test_functionSignatures.m
rename to tests/test_function_signatures.m
index 59730c28..f63b4ad8 100644
--- a/tests/test_functionSignatures.m
+++ b/tests/test_function_signatures.m
@@ -1,4 +1,12 @@
-function test_functionSignatures(pth)
+function test_suite = test_function_signatures %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_function_signatures_basic()
   % Test functionSignatures.json file
   % The functionSignatures file is used by Matlab to provide code suggestions
   % and completions for functions.
@@ -13,12 +21,11 @@ function test_functionSignatures(pth)
   %
   % Copyright (C) 2020--, BIDS-MATLAB developers
 
-  if ~nargin
-    % parent directory of the tests directory
-    pth = fullfile(pwd, '..');
-  end
+  root_dir = fullfile(fileparts(mfilename('fullpath')), '..');
 
   % Run a smoke test - see if the file is a readable json
   signatures = bids.util.jsondecode( ...
-                                    fullfile(pth, 'functionSignatures.json'));
+                                    fullfile(root_dir, 'functionSignatures.json'));
   assert(isstruct(signatures));
+
+end
diff --git a/tests/test_get_metadata.m b/tests/test_get_metadata.m
index be0e32be..383c3860 100644
--- a/tests/test_get_metadata.m
+++ b/tests/test_get_metadata.m
@@ -39,20 +39,20 @@ function test_get_metadata_basic()
   BIDS = bids.layout(pth);
 
   %% test func metadata base directory
-  metadata = bids.query(BIDS, 'metadata', 'type', 'bold');
+  metadata = bids.query(BIDS, 'metadata', 'suffix', 'bold');
   % assert(metadata.RepetitionTime == func.RepetitionTime);
 
   %% test func metadata subject 01
-  metadata = bids.query(BIDS, 'metadata', 'sub', '01', 'type', 'bold');
+  metadata = bids.query(BIDS, 'metadata', 'sub', '01', 'suffix', 'bold');
   % assert(metadata.RepetitionTime == func_sub_01.RepetitionTime);
 
   %% test anat metadata base directory
-  metadata = bids.query(BIDS, 'metadata', 'type', 'T1w');
+  metadata = bids.query(BIDS, 'metadata', 'suffix', 'T1w');
   % assert(metadata.FlipAngle == anat.FlipAngle);
 
   %% test anat metadata subject 01
-  metadata = bids.query(BIDS, 'metadata', 'sub', '01', 'type', 'T1w');
-  assert(metadata.FlipAngle == anat_sub_01.FlipAngle);
-  assert(strcmp(metadata.Manufacturer, anat_sub_01.Manufacturer));
+  metadata = bids.query(BIDS, 'metadata', 'sub', '01', 'suffix', 'T1w');
+  assertEqual(metadata.FlipAngle, anat_sub_01.FlipAngle);
+  assertEqual(metadata.Manufacturer, anat_sub_01.Manufacturer);
 
 end
diff --git a/tests/test_keep_file.m b/tests/test_keep_file.m
new file mode 100644
index 00000000..71226937
--- /dev/null
+++ b/tests/test_keep_file.m
@@ -0,0 +1,57 @@
+function test_suite = test_keep_file %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+
+end
+
+function test_keep_file_basic()
+
+  file_struct = struct('suffix', 'bold', ...
+                       'entities', struct('task', 'balloon', ...
+                                          'ses', '01', ...
+                                          'sub', '02'));
+
+  options = {'ses', {'01'}};
+  assertEqual(bids.internal.keep_file_for_query(file_struct, options), true);
+
+  options = {'ses', {'02'}};
+  assertEqual(bids.internal.keep_file_for_query(file_struct, options), false);
+
+  options = {'ses', {'01', '02'}};
+  assertEqual(bids.internal.keep_file_for_query(file_struct, options), true);
+
+  options = {'suffix', {'T1w'}};
+  assertEqual(bids.internal.keep_file_for_query(file_struct, options), false);
+
+  options = {'suffix', {'bold'}};
+  assertEqual(bids.internal.keep_file_for_query(file_struct, options), true);
+
+  options = {'suffix', {'T1w', 'bold'}};
+  assertEqual(bids.internal.keep_file_for_query(file_struct, options), true);
+
+  options = {'sub', {'02'}
+             'ses', {'01'}};
+  assertEqual(bids.internal.keep_file_for_query(file_struct, options), true);
+
+  options = {'sub', {'02'}
+             'ses', {'02'}};
+  assertEqual(bids.internal.keep_file_for_query(file_struct, options), false);
+
+  options = {'sub', {'02'}
+             'ses', {'01'}
+             'suffix', {'T1w'}};
+  assertEqual(bids.internal.keep_file_for_query(file_struct, options), false);
+
+  options = {'sub', {'02'}
+             'task', {'balloon'}};
+  assertEqual(bids.internal.keep_file_for_query(file_struct, options), true);
+
+  file_struct = struct('suffix', 'T1w', ...
+                       'entities', struct('ses', '01', ...
+                                          'sub', '02'));
+  assertEqual(bids.internal.keep_file_for_query(file_struct, options), false);
+
+end
diff --git a/tests/test_layout.m b/tests/test_layout.m
new file mode 100644
index 00000000..30eb4c83
--- /dev/null
+++ b/tests/test_layout.m
@@ -0,0 +1,17 @@
+function test_suite = test_layout %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_layout_smoke_test()
+
+  pth_bids_example = get_test_data_dir();
+
+  BIDS = bids.layout(fullfile(pth_bids_example, 'genetics_ukbb'));
+
+  BIDS = bids.layout(fullfile(pth_bids_example, 'ds210'));
+
+end
diff --git a/tests/test_layout_derivatives.m b/tests/test_layout_derivatives.m
new file mode 100644
index 00000000..ae6923ad
--- /dev/null
+++ b/tests/test_layout_derivatives.m
@@ -0,0 +1,71 @@
+function test_suite = test_layout_derivatives %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_layout_schemaless()
+
+  pth_bids_example = get_test_data_dir();
+
+  use_schema = false();
+
+  BIDS = bids.layout(fullfile(pth_bids_example, ...
+                              'ds000001-fmriprep'), use_schema);
+
+  modalities = {'anat', 'figures', 'func'};
+  assertEqual(bids.query(BIDS, 'modalities'), modalities);
+
+  % Those fail for now
+  data = bids.query(BIDS, 'data', ...
+                    'sub', '10', ...
+                    'modality', 'func', ...
+                    'suffix', 'bold', ...
+                    'run', '1', ...
+                    'res', '2');
+
+  basename = bids.internal.file_utils(data, 'basename');
+  assertEqual(basename, {
+                         ['sub-10_task-balloonanalogrisktask_run-1', ...
+                          '_space-MNI152NLin2009cAsym_res-2_desc-preproc_bold.nii']
+                        });
+
+  data = bids.query(BIDS, 'data', ...
+                    'sub', '10', ...
+                    'modality', 'func', ...
+                    'suffix', 'bold', ...
+                    'run', '1', ...
+                    'space', 'MNI152NLin6Asym');
+
+  basename = bids.internal.file_utils(data, 'basename');
+  assertEqual(basename, {
+                         ['sub-10_task-balloonanalogrisktask_run-1', ...
+                          '_space-MNI152NLin6Asym_desc-smoothAROMAnonaggr_bold.nii']
+                        });
+end
+
+function test_layout_nested_derivatives()
+
+  pth_bids_example = get_test_data_dir();
+
+  use_schema = false();
+
+  BIDS = bids.layout(fullfile(pth_bids_example, ...
+                              'ds000117', ...
+                              'derivatives', ...
+                              'meg_derivatives'), use_schema);
+
+  modalities = {'meg'};
+  assertEqual(bids.query(BIDS, 'modalities'), modalities);
+
+  data = bids.query(BIDS, 'data', ...
+                    'sub', '01', ...
+                    'run', '01', ...
+                    'proc', 'sss', ...
+                    'suffix', 'meg');
+  basename = bids.internal.file_utils(data, 'basename');
+  assertEqual(basename, {'sub-01_ses-meg_task-facerecognition_run-01_proc-sss_meg'});
+
+end
diff --git a/tests/test_load_schema.m b/tests/test_load_schema.m
new file mode 100644
index 00000000..6b267a06
--- /dev/null
+++ b/tests/test_load_schema.m
@@ -0,0 +1,34 @@
+function test_suite = test_load_schema %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_load_schema_path()
+
+  use_schema = fullfile(fileparts(mfilename('fullpath')), 'schema');
+
+  schema = bids.schema.load_schema(use_schema);
+
+  assert(isfield(schema, 'base'));
+  assert(isfield(schema, 'subfolder_1'));
+  assert(isfield(schema.subfolder_1, 'sub'));
+  assert(~isfield(schema, 'subfolder_4'));
+
+  % some recursive aspects are not implemented yet
+  %     assert(isfield(schema.subfolder_2, 'subfolder_3'));
+  %     assert(isfield(schema.subfolder_2.subfolder_3, 'sub'));
+
+end
+
+function test_load_schema_schemaless()
+
+  use_schema = false();
+
+  schema = bids.schema.load_schema(use_schema);
+
+  assertEqual(schema, []);
+
+end
diff --git a/tests/test_parse_filename.m b/tests/test_parse_filename.m
new file mode 100644
index 00000000..36ebdfe2
--- /dev/null
+++ b/tests/test_parse_filename.m
@@ -0,0 +1,54 @@
+function test_suite = test_parse_filename %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_parse_filename_basic()
+
+  filename = '../sub-16/anat/sub-16_ses-mri_run-1_acq-hd_T1w.nii.gz';
+  output = bids.internal.parse_filename(filename);
+
+  expected = struct( ...
+                    'filename', 'sub-16_ses-mri_run-1_acq-hd_T1w.nii.gz', ...
+                    'suffix', 'T1w', ...
+                    'ext', '.nii.gz', ...
+                    'entities', struct('sub', '16', ...
+                                       'ses', 'mri', ...
+                                       'run', '1', ...
+                                       'acq', 'hd'));
+
+  assertEqual(output, expected);
+
+end
+
+function test_parse_filename_fields()
+
+  filename = '../sub-16/anat/sub-16_ses-mri_run-1_acq-hd_T1w.nii.gz';
+  fields = {'sub', 'ses', 'run', 'acq', 'ce', 'rec', 'part'};
+  output = bids.internal.parse_filename(filename);
+
+  expected = struct( ...
+                    'filename', 'sub-16_ses-mri_run-1_acq-hd_T1w.nii.gz', ...
+                    'suffix', 'T1w', ...
+                    'ext', '.nii.gz', ...
+                    'entities', struct('sub', '16', ...
+                                       'ses', 'mri', ...
+                                       'run', '1', ...
+                                       'acq', 'hd'));
+
+  assertEqual(output, expected);
+
+end
+
+function test_parse_filename_wrong_template()
+
+  filename = '../sub-16/anat/sub-16_ses-mri_run-1_acq-hd_T1w.nii.gz';
+
+  assertWarning( ...
+                @()bids.internal.parse_filename(filename, {'echo'}), ...
+                'bidsMatlab:noMatchingTemplate');
+
+end
diff --git a/tests/test_return_modality_entities.m b/tests/test_return_modality_entities.m
new file mode 100644
index 00000000..d9a38cd6
--- /dev/null
+++ b/tests/test_return_modality_entities.m
@@ -0,0 +1,19 @@
+function test_suite = test_return_modality_entities %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_return_modality_entities_basic
+
+  schema = bids.schema.load_schema();
+
+  entities = bids.schema.return_modality_entities(schema.datatypes.func(1), schema);
+
+  expected_output = {'sub', 'ses', 'task', 'acq', 'ce', 'rec', 'dir', 'run', 'echo', 'part'};
+
+  assertEqual(entities, expected_output);
+
+end
diff --git a/tests/test_return_modality_extensions.m b/tests/test_return_modality_extensions.m
new file mode 100644
index 00000000..551adc49
--- /dev/null
+++ b/tests/test_return_modality_extensions.m
@@ -0,0 +1,16 @@
+function test_suite = test_return_modality_extensions %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_return_modality_extensions_basic
+
+  schema = bids.schema.load_schema();
+
+  extensions = bids.internal.return_modality_extensions(schema.datatypes.func(1));
+  assertEqual(extensions, '(.nii.gz|.nii){1}');
+
+end
diff --git a/tests/test_return_modality_regular_expression.m b/tests/test_return_modality_regular_expression.m
new file mode 100644
index 00000000..b3984f5b
--- /dev/null
+++ b/tests/test_return_modality_regular_expression.m
@@ -0,0 +1,29 @@
+function test_suite = test_return_modality_regular_expression %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_return_modality_regular_expression_basic
+
+  schema = bids.schema.load_schema();
+
+  regular_expression = bids.internal.return_modality_regular_expression(schema.datatypes.anat(1));
+
+  expected_expression = ['^%s.*', ...
+                         '_(T1w|T2w|PDw|T2starw|FLAIR|inplaneT1|inplaneT2|PDT2|angio|', ...
+                         'T2star|FLASH|PD', ... % deprecated suffixes
+                         '){1}', ...
+                         '(.nii.gz|.nii){1}$'];
+
+  assertEqual(regular_expression, expected_expression);
+
+  data_dir = fullfile(fileparts(mfilename('fullpath')), 'data', 'MoAEpilot', 'sub-01', 'anat');
+  subject_name = 'sub-01';
+  file = bids.internal.file_utils('List', data_dir, sprintf(expected_expression, subject_name));
+
+  assertEqual(file, 'sub-01_T1w.nii.gz');
+
+end
diff --git a/tests/test_return_modality_suffixes.m b/tests/test_return_modality_suffixes.m
new file mode 100644
index 00000000..19212185
--- /dev/null
+++ b/tests/test_return_modality_suffixes.m
@@ -0,0 +1,16 @@
+function test_suite = test_return_modality_suffixes %#ok<*STOUT>
+  try % assignment of 'localfunctions' is necessary in Matlab >= 2016
+    test_functions = localfunctions(); %#ok<*NASGU>
+  catch % no problem; early Matlab versions can use initTestSuite fine
+  end
+  initTestSuite;
+end
+
+function test_return_modality_suffixes_basic
+
+  schema = bids.schema.load_schema();
+
+  suffixes = bids.internal.return_modality_suffixes(schema.datatypes.func(1));
+  assertEqual(suffixes, '_(bold|cbv|sbref){1}');
+
+end
diff --git a/tests/test_tsvread.m b/tests/test_tsvread.m
index fec8e190..d5ca578a 100644
--- a/tests/test_tsvread.m
+++ b/tests/test_tsvread.m
@@ -30,13 +30,13 @@ function test_tsvread_basic()
   %% test tsvread on tsv file
   tsv_file = fullfile(pth, 'sub-01', 'func', 'sub-01_task-auditory_events.tsv');
   output = bids.util.tsvread(tsv_file);
-  assert(isequal(output.onset', events.onset));
+  assertEqual(output.onset', events.onset);
 
   %% test tsvread on zipped tsv file
   output = bids.util.tsvread(fullfile( ...
                                       fileparts(mfilename('fullpath')), ....
                                       'data', ...
                                       'sub-01_task-auditory_events.tsv.gz'));
-  assert(isequal(output.onset', events.onset));
+  assertEqual(output.onset', events.onset);
 
 end
diff --git a/tests/test_tsvwrite.m b/tests/test_tsvwrite.m
index 3285000b..751685eb 100644
--- a/tests/test_tsvwrite.m
+++ b/tests/test_tsvwrite.m
@@ -51,19 +51,19 @@ function test_tsvwrite_basic()
   C = textscan(FID, '%s%s%s%s%s', 'Delimiter', '\t', 'EndOfLine', '\n');
 
   % check header
-  assert(isequal(C{4}{1}, 'speed'));
+  assertEqual(C{4}{1}, 'speed');
 
   % check that empty values are entered as NaN: logFile(1,1).speed
-  assert(isequal(C{4}{2}, 'n/a'));
+  assertEqual(C{4}{2}, 'n/a');
 
   % check that missing fields are entered as NaN: logFile(2,1).speed
-  assert(isequal(C{4}{3}, 'n/a'));
+  assertEqual(C{4}{3}, 'n/a');
 
   % check that NaN are written as : logFile(2,1).onset
-  assert(isequal(C{1}{3}, 'n/a')); %
+  assertEqual(C{1}{3}, 'n/a'); %
 
   % check values entered properly: logFile(2,1).is_fixation
-  assert(isequal(str2double(C{5}{3}), 3));
+  assertEqual(str2double(C{5}{3}), 3);
 
   %% test tsvread on tsv file using cell input
   % TO DO?
@@ -79,6 +79,6 @@ function test_tsvwrite_basic()
 
   FID = fopen(tsv_file, 'r');
   C = textscan(FID, '%s%s', 'Delimiter', '\t', 'EndOfLine', '\n');
-  assert(isequal(C{1}{2}, 'n/a')); %
+  assertEqual(C{1}{2}, 'n/a'); %
 
 end