dsCheckSpecification

CHECKSPECIFICATION - standardize specification structure and auto-populate missing fields

 Usage:
   specification=dsCheckSpecification(specification)

 Input: DynaSim specification structure or equations

 Output:
   - DynaSim specification structure (standardized)
     .populations(i) (required): contains info for defining independent population models
       .name (default: 'pop1')      : name of population
       .size (default: 1)           : number of elements in population (i.e., # cells)
       .equations (required)        : string listing equations (see NOTE 1)
       .mechanism_list (default: []): cell array listing mechanisms (see NOTE 2)
       .parameters (default: [])    : parameters to assign across all equations in
                                      the population. provide as cell array list of
                                      key/value pairs, like
                                      {'param1',value1,'param2',value2,...}
       .conditionals (default: [])  : (see NOTE 3) if-then conditional actions
         .namespace (auto)    : string giving the namespace of the condition 
                                (pop_ or pop_mech_)
         .condition (required): string giving the condition to check
         .action (required)   : what to do if the condition is met
         .else (default: [])  : what to do if the condition is not met
       .monitors (default: [])      : (see NOTE 3) substructure with fields specifying
                                      what to record on each step of numerical
                                      integration in addition to state
                                      variables.
       .model (default: [])   : optional DynaSim model structure
     .connections(i) (default: []): contains info for linking population models
       .source (required if >1 pops): name of source population (see NOTE 7)
       .target (required if >1 pops): name of target population
       .mechanism_list (required)   : list of mechanisms that link two populations
       .parameters (default: [])    : parameters to assign across all equations in
                                      mechanisms in this connection's mechanism_list.

 Notes:
   - NOTE 1: .equations can be an equation string, cell array listing equation
       strings, or a file name pointing to a model / equations stored on disk
       (accepted file types: .eqns (equations of population), .m (function defining
       a model structure), ...)

   - NOTE 2: .mechanism_list is a cell array listing names of mechanisms to be
       included in the population or used to connect two populations. each mechanism
       name must have a mechanism file with the same name somewhere in the search
       path (the file should have extension .mech). The search path starts with the
       current directory, then the subdirectories of [dynasim]/models, and lastly
       the full matlab search path.

   - NOTE 3: conditionals and monitors are most easily specified by including
       them in .equations.
     - Example:
         spec.populations.equations='dv/dt=-v; if(v<eps)(v=10); monitor o=v^2'
         data=dsSimulate(spec); figure; plot(data.time,data.pop1_o);

   - NOTE 4: "pops" can be used instead of "populations". "cons" can be used
       instead of "connections".

   - NOTE 5: all population info can be embedded in the equation string.
       Specify name by starting the string with 'NAME: *' (e.g., 'E: dv/dt=-v').
       Specify size by including [SIZE] after the state variable (e.g., 'dv[5]/dt=-v').
       Specify mechanism_list by including cell array listing mechanism names
       without single quotes (e.g., 'dv/dt=@current; {iNa,iK}').

   - NOTE 6: the mechanism linker target IDENTIFIER used to link mechanism
       variables and functions to population equations can be overriden by including
       @NEWIDENTIFIER in the equations and after the mechanism name. (e.g.,
       'dv/dt=@M; {iNa,iK}@M'; or .mechanism_list={'iNa@M','iK@M'}).

   - NOTE 7: "direction" can be used instead of "source" and "target". The
       syntax is "SOURCE->TARGET" or "TARGET<-SOURCE". Either will be properly split
       into .source and .target fields. SOURCE and TARGET must be existing
       population names.

 Examples:
   - Example 1: obtain empty specification structure with all fields
       specification=dsCheckSpecification([]);

   - Example 2: standardize existing specification
       specification=dsCheckSpecification(specification)

   - Example 3: standardize equations in cell array
       eqns={
         's=10; r=27; b=2.666';
         'dx/dt=s*(y-x)';
         'dy/dt=r*x-y-x*z';
         'dz/dt=-b*z+x*y';
       specification=dsCheckSpecification(eqns);

   - Example 4: standardize equations in character array
       eqns='tau=10; R=10; E=-70; dV/dt=(E-V+R*1.55)/tau; if(V>-55)(V=-75)';
       specification=dsCheckSpecification(eqns);

   - Example 5: standardize specification with compact field names
       s.pops.size=10;
       s.pops.equations='dv/dt=-v';
       s.cons.mechanism_list='iGABAa';
       s=dsCheckSpecification(s)

   - Example 6: standardize specification with everything in equation string
       s.pops.equations='E:dv[10]/dt=@M+I; {iNa,iK}@M; I=10';
       s=dsCheckSpecification(s)
 
     Example 7: standardize equations from predefined population model
       specification=dsCheckSpecification('HH');

 See also: dsGenerateModel, dsCheckModel

 Author: Jason Sherfey, PhD <jssherfey@gmail.com>
 Copyright (C) 2016 Jason Sherfey, Boston University, USA

0001 function spec = dsCheckSpecification(specification, varargin)
0002 %CHECKSPECIFICATION - standardize specification structure and auto-populate missing fields
0003 %
0004 % Usage:
0005 %   specification=dsCheckSpecification(specification)
0006 %
0007 % Input: DynaSim specification structure or equations
0008 %
0009 % Output:
0010 %   - DynaSim specification structure (standardized)
0011 %     .populations(i) (required): contains info for defining independent population models
0012 %       .name (default: 'pop1')      : name of population
0013 %       .size (default: 1)           : number of elements in population (i.e., # cells)
0014 %       .equations (required)        : string listing equations (see NOTE 1)
0015 %       .mechanism_list (default: []): cell array listing mechanisms (see NOTE 2)
0016 %       .parameters (default: [])    : parameters to assign across all equations in
0017 %                                      the population. provide as cell array list of
0018 %                                      key/value pairs, like
0019 %                                      {'param1',value1,'param2',value2,...}
0020 %       .conditionals (default: [])  : (see NOTE 3) if-then conditional actions
0021 %         .namespace (auto)    : string giving the namespace of the condition
0022 %                                (pop_ or pop_mech_)
0023 %         .condition (required): string giving the condition to check
0024 %         .action (required)   : what to do if the condition is met
0025 %         .else (default: [])  : what to do if the condition is not met
0026 %       .monitors (default: [])      : (see NOTE 3) substructure with fields specifying
0027 %                                      what to record on each step of numerical
0028 %                                      integration in addition to state
0029 %                                      variables.
0030 %       .model (default: [])   : optional DynaSim model structure
0031 %     .connections(i) (default: []): contains info for linking population models
0032 %       .source (required if >1 pops): name of source population (see NOTE 7)
0033 %       .target (required if >1 pops): name of target population
0034 %       .mechanism_list (required)   : list of mechanisms that link two populations
0035 %       .parameters (default: [])    : parameters to assign across all equations in
0036 %                                      mechanisms in this connection's mechanism_list.
0037 %
0038 % Notes:
0039 %   - NOTE 1: .equations can be an equation string, cell array listing equation
0040 %       strings, or a file name pointing to a model / equations stored on disk
0041 %       (accepted file types: .eqns (equations of population), .m (function defining
0042 %       a model structure), ...)
0043 %
0044 %   - NOTE 2: .mechanism_list is a cell array listing names of mechanisms to be
0045 %       included in the population or used to connect two populations. each mechanism
0046 %       name must have a mechanism file with the same name somewhere in the search
0047 %       path (the file should have extension .mech). The search path starts with the
0048 %       current directory, then the subdirectories of [dynasim]/models, and lastly
0049 %       the full matlab search path.
0050 %
0051 %   - NOTE 3: conditionals and monitors are most easily specified by including
0052 %       them in .equations.
0053 %     - Example:
0054 %         spec.populations.equations='dv/dt=-v; if(v<eps)(v=10); monitor o=v^2'
0055 %         data=dsSimulate(spec); figure; plot(data.time,data.pop1_o);
0056 %
0057 %   - NOTE 4: "pops" can be used instead of "populations". "cons" can be used
0058 %       instead of "connections".
0059 %
0060 %   - NOTE 5: all population info can be embedded in the equation string.
0061 %       Specify name by starting the string with 'NAME: *' (e.g., 'E: dv/dt=-v').
0062 %       Specify size by including [SIZE] after the state variable (e.g., 'dv[5]/dt=-v').
0063 %       Specify mechanism_list by including cell array listing mechanism names
0064 %       without single quotes (e.g., 'dv/dt=@current; {iNa,iK}').
0065 %
0066 %   - NOTE 6: the mechanism linker target IDENTIFIER used to link mechanism
0067 %       variables and functions to population equations can be overriden by including
0068 %       @NEWIDENTIFIER in the equations and after the mechanism name. (e.g.,
0069 %       'dv/dt=@M; {iNa,iK}@M'; or .mechanism_list={'iNa@M','iK@M'}).
0070 %
0071 %   - NOTE 7: "direction" can be used instead of "source" and "target". The
0072 %       syntax is "SOURCE->TARGET" or "TARGET<-SOURCE". Either will be properly split
0073 %       into .source and .target fields. SOURCE and TARGET must be existing
0074 %       population names.
0075 %
0076 % Examples:
0077 %   - Example 1: obtain empty specification structure with all fields
0078 %       specification=dsCheckSpecification([]);
0079 %
0080 %   - Example 2: standardize existing specification
0081 %       specification=dsCheckSpecification(specification)
0082 %
0083 %   - Example 3: standardize equations in cell array
0084 %       eqns={
0085 %         's=10; r=27; b=2.666';
0086 %         'dx/dt=s*(y-x)';
0087 %         'dy/dt=r*x-y-x*z';
0088 %         'dz/dt=-b*z+x*y';
0089 %       specification=dsCheckSpecification(eqns);
0090 %
0091 %   - Example 4: standardize equations in character array
0092 %       eqns='tau=10; R=10; E=-70; dV/dt=(E-V+R*1.55)/tau; if(V>-55)(V=-75)';
0093 %       specification=dsCheckSpecification(eqns);
0094 %
0095 %   - Example 5: standardize specification with compact field names
0096 %       s.pops.size=10;
0097 %       s.pops.equations='dv/dt=-v';
0098 %       s.cons.mechanism_list='iGABAa';
0099 %       s=dsCheckSpecification(s)
0100 %
0101 %   - Example 6: standardize specification with everything in equation string
0102 %       s.pops.equations='E:dv[10]/dt=@M+I; {iNa,iK}@M; I=10';
0103 %       s=dsCheckSpecification(s)
0104 %
0105 %     Example 7: standardize equations from predefined population model
0106 %       specification=dsCheckSpecification('HH');
0107 %
0108 % See also: dsGenerateModel, dsCheckModel
0109 %
0110 % Author: Jason Sherfey, PhD <jssherfey@gmail.com>
0111 % Copyright (C) 2016 Jason Sherfey, Boston University, USA
0112 
0113 %% localfn output
0114 if ~nargin
0115   spec = localfunctions; % output var name specific to this fn
0116   return
0117 end
0118 
0119 %% auto_gen_test_data_flag argin
0120 options = dsCheckOptions(varargin,{'auto_gen_test_data_flag',0,{0,1}},false);
0121 if options.auto_gen_test_data_flag
0122   varargs = varargin;
0123   varargs{find(strcmp(varargs, 'auto_gen_test_data_flag'))+1} = 0;
0124   varargs(end+1:end+2) = {'unit_test_flag',1};
0125   argin = [{specification}, varargs]; % specific to this function
0126 end
0127 
0128 % check if input is a string or cell with equations and package in spec structure
0129 if ischar(specification) || iscell(specification)
0130   spec.populations.equations=specification;
0131 elseif isstruct(specification)
0132   spec=specification;
0133 elseif isempty(specification)
0134   spec=struct;
0135 else
0136   error('specification must be a DynaSim specification structure or a string with equations or sub-model filename.');
0137 end
0138 
0139 spec=backward_compatibility(spec, varargin{:});
0140 pop_field_order={'name','size','equations','mechanism_list','mechanisms','parameters',...
0141   'conditionals','monitors','model'};
0142 con_field_order={'source','target','mechanism_list','mechanisms','parameters'};
0143 
0144 if ~isfield(spec,'populations')
0145   spec.populations.name='pop1';
0146 end
0147 
0148 if ~isfield(spec,'connections')
0149   spec.connections=[];
0150 end
0151 
0152 if ~isfield(spec,'mechanisms')
0153   spec.mechanisms=[];
0154 end
0155 
0156 % 1.0 standardize populations
0157 if ~isfield(spec.populations,'name')
0158   spec.populations(1).name='pop1';
0159 end
0160 
0161 if ~isfield(spec.populations,'size')
0162   spec.populations(1).size=1;
0163 end
0164 
0165 if ~isfield(spec.populations,'equations')
0166   spec.populations(1).equations=[];
0167 end
0168 
0169 if ~isfield(spec.populations,'mechanism_list')
0170   spec.populations(1).mechanism_list=[];
0171 end
0172 
0173 if ~isfield(spec.populations,'mechanisms')
0174   spec.populations(1).mechanisms=[];
0175 end
0176 
0177 if ~isfield(spec.populations,'parameters')
0178   spec.populations(1).parameters={};
0179 end
0180 
0181 if ~isfield(spec.populations,'conditionals')
0182   spec.populations(1).conditionals=[];
0183 end
0184 
0185 if ~isfield(spec.populations,'monitors')
0186   spec.populations(1).monitors=[];
0187 end
0188 
0189 if ~isfield(spec.populations,'model')
0190   spec.populations(1).model=[];
0191 end
0192 
0193 % move compartments into populations if present
0194 if isfield(spec,'compartments')
0195   npops=length(spec.populations);
0196   fields=fieldnames(spec.compartments);
0197   for i=1:length(spec.compartments)
0198     for f=1:length(fields)
0199       spec.populations(npops+i).(fields{f})=spec.compartments(i).(fields{f});
0200     end
0201   end
0202   spec=rmfield(spec,'compartments');
0203 end
0204 
0205 % special cases of equation specification:
0206 for i=1:length(spec.populations)
0207   eqn=spec.populations(i).equations;
0208   if ~isempty(eqn) && ischar(eqn)
0209     % check for predefined population equations
0210     if exist([eqn '.eqns'],'file')
0211       eqn=[eqn '.eqns'];
0212     elseif exist([eqn '.pop'],'file')
0213       eqn=[eqn '.pop'];
0214     end
0215     if exist(eqn,'file')
0216       % load equations from file
0217       spec.populations(i).equations=dsReadText(eqn);
0218     elseif ~isempty(regexp(eqn,'\[[a-z_A-Z].*\]','match','once'))
0219       % split equations with '[...][...]...[...]' into multiple populations
0220       % create extra population
0221       tmp=regexp(eqn(2:end-1),'\],?\s*\[','split');
0222       spec.populations(i).equations=tmp{1};
0223       for j=2:length(tmp)
0224         spec.populations(end+1)=spec.populations(i);
0225         spec.populations(end).equations=tmp{j};
0226         spec.populations(end).name=sprintf('pop%g',length(spec.populations));
0227       end
0228     end
0229   end
0230 end
0231 
0232 % standardize each population separately
0233 for i=1:length(spec.populations)
0234   % population names
0235   if isempty(spec.populations(i).name)
0236     spec.populations(i).name=sprintf('pop%g',i);
0237   end
0238   
0239   % population sizes
0240   if isempty(spec.populations(i).size)
0241     spec.populations(i).size=1;
0242   end
0243   
0244   % make mechanism list a cell array of mechanism names
0245   if ischar(spec.populations(i).mechanism_list)
0246     spec.populations(i).mechanism_list={spec.populations(i).mechanism_list};
0247   end
0248   
0249   % parameter cell arrays
0250   if ~iscell(spec.populations(i).parameters)
0251     spec.populations(i).parameters={};
0252   end
0253   
0254   % standardize equations
0255   if ~isempty(spec.populations(i).equations)
0256     % convert cell array of equations into character array
0257     if iscell(spec.populations(i).equations)
0258       eqns=spec.populations(i).equations;
0259       for k=1:length(eqns)
0260         eqn=eqns{k};
0261         if ~isempty(eqn) && ~strcmp(eqn(end),';')
0262           eqns{k}(end+1)=';';
0263         end
0264       end
0265       spec.populations(i).equations=[eqns{:}];
0266     end
0267     
0268     % extract name from equations (e.g., TC:...)
0269     eqn=spec.populations(i).equations;
0270     name=regexp(eqn,'^\w+:','match','once');
0271     if ~isempty(name)
0272       % remove name indicator from equation
0273       eqn=strrep(eqn,name,'');
0274       
0275       % store name in specification
0276       name=regexp(name,'^(\w+):','tokens','once');
0277       spec.populations(i).equations=eqn;
0278       if strcmp(spec.populations(i).name,sprintf('pop%g',i))
0279         % replace default name with the name from equations
0280         spec.populations(i).name=name{1};
0281       end      
0282     end
0283     
0284     % extract size from equations if present (eg, v[4]'=.., dv[4]/dt=...)
0285     eqn=spec.populations(i).equations;
0286     %pattern='((\w+(\[\d+\])'')|(d\w+(\[\d+\])/dt))\s*='; % support size spec, dv[4]/dt
0287     pattern='((\w+(\[[\d,]+\])'')|(d\w+(\[[\d,]+\])/dt))\s*='; % support size spec, dv[4]/dt and dv[4,5]/dt
0288     
0289     % extract all differentials with size specification
0290     LHSs=regexp(eqn,pattern,'match');
0291     if ~isempty(LHSs)
0292       % extract sizes from all differentials (eg, 4 from v[4] or dv[4]/dt)
0293       for k=1:length(LHSs)
0294         tmp=regexp(LHSs{k},'\w+\[([\d,]+)\]''','tokens','once');
0295         if isempty(tmp)
0296           tmp=regexp(LHSs{k},'d\w+\[([\d,]+)\]/dt','tokens','once');
0297         end
0298         sz=cellfun(@str2double,regexp(tmp{1},',','split'));
0299         % check that all vars in same population have same size
0300         if k==1
0301           sz_first=sz;
0302         elseif sz~=sz_first
0303           error('all variables in same population must have same size. split ODEs with different sizes into different populations.');
0304         end        
0305         % remove size from ODE in population equations
0306         old=LHSs{k};
0307         new=strrep(LHSs{k},['[' tmp{1} ']'],'');
0308         eqn=strrep(eqn,old,new);
0309       end      
0310       spec.populations(i).equations=eqn;
0311       spec.populations(i).size=sz;
0312     end
0313     
0314     % add mechanisms embedded in equations to mechanism_list ({M1,M2,...})
0315     % ----------
0316     % todo: make the following a subfunction and apply it also to connection
0317     % mechanism lists (eg, for cons.mechanism_list='{AMPA,NMDA}@M')
0318     % ----------
0319     % extract mechanism list from equations
0320      mech_lists=regexp(spec.populations(i).equations,'\s*(\w+:)?{[\w\d@:,]*}\s*(@\w+)?;?\s*','match');
0321     % test: mech_list=regexp('v''=@M+sin(2*pi*t); {iNa, iK}','{.*}','match');
0322     if ~isempty(mech_lists)
0323       for k=1:length(mech_lists)
0324         mech_list=strtrim(mech_lists{k});
0325         
0326         % remove mechanism list from equations
0327         spec.populations(i).equations=strtrim(strrep(spec.populations(i).equations,mech_list,''));
0328         
0329         % append external link alias to each internal mechanism name (eg, {a,b}@M, alias @M)
0330         external_link=regexp(mech_list,'}(@[\w\d]+;?)','tokens');
0331         if ~isempty(external_link)
0332           % get external link alias
0333           external_link=[external_link{:}];
0334           
0335           % remove external link alias from mech_list
0336           mech_list=strrep(mech_list,external_link{1},'');
0337           
0338           % remove ';' from alias before appending to mech names
0339           external_link=strrep(external_link{1},';','');
0340           
0341           % get list of mechanism names in cell array
0342           words=regexp(mech_list(2:end-1),',','split');
0343           
0344           % append external link alias to each mechanism name
0345           for w=1:length(words)
0346             mech_list=strrep(mech_list,words{w},[words{w} external_link]);
0347           end
0348         end
0349         % prepend host name to each internal mechanism name (eg,
0350         % infbrain:{a,b} -> {infbrain:a,infbrain:b}
0351         host_name=regexp(mech_list,';?\s*([\w\d]+):{','tokens','once');
0352         if ~isempty(host_name)
0353           % get external link alias
0354           host_name=[host_name{:}];
0355           
0356           % remove external link alias from mech_list
0357           mech_list=strrep(mech_list,[host_name ':'],'');
0358           
0359           % get list of mechanism names in cell array
0360           words=regexp(mech_list(2:end-1),',','split');
0361           
0362           % append external link alias to each mechanism name
0363           for w=1:length(words)
0364             mech_list=strrep(mech_list,words{w},[host_name ':' words{w}]);
0365           end
0366         end
0367         % split into list of mechanism names
0368         mechanisms=regexp(mech_list,'[\w:@]+','match');
0369         
0370         % append mechanism from equations to mechanism_list
0371         if iscell(spec.populations(i).mechanism_list)
0372           spec.populations(i).mechanism_list=cat(2,mechanisms,spec.populations(i).mechanism_list);
0373         else
0374           spec.populations(i).mechanism_list=mechanisms;
0375         end
0376       end
0377     end
0378     % extract population-level parameters from equations
0379     param_name={};
0380     param_value={};
0381     eqn=spec.populations(i).equations;
0382     p=getfield(dsParseModelEquations(eqn, varargin{:}),'parameters');
0383     if ~isempty(p)
0384       param_name=cat(1,param_name,fieldnames(p));
0385       param_value=cat(1,param_value,struct2cell(p));
0386     end
0387     % extract mechanism-specific parameters defined in master equations
0388     % eg) eqn='dv/dt=@current+10; monitor iAMPA.functions; iNa.IC_noise=10; iK.g=g; g=3';
0389     o=regexp(eqn,';\s*[a-zA-Z]+\w*\.[a-zA-Z]+\w*\s*=[a-z_A-Z0-9\.]+','match'); 
0390       % eg) '; MECH.PARAM=VALUE', assumes param is not defined at start of equation string
0391     % add mechanism-specific keys (MECH.PARAM) and vals to p
0392     if ~isempty(o)
0393       % remove leading semicolons
0394       oo=regexprep(o,';',''); % {'MECH1.PARAM1=VAL1','MECH2.PARAM2=VAL2',..}
0395       for l=1:length(oo)
0396         tmp=strtrim(regexp(oo{l},'=','split')); % {'MECH1.PARAM1','VAL1'}
0397         param_name{end+1}=tmp{1}; % 'MECH1.PARAM1'
0398         param_value{end+1}=tmp{2}; % 'VAL1'
0399         % remove from equations
0400         eqn=strrep(eqn,o{l},'');
0401       end
0402       spec.populations(i).equations=eqn;
0403     end    
0404     % move user-defined parameters from equations to the parameters field
0405 %     if ~isempty(p)
0406 %       param_name=fieldnames(p);
0407 %       param_value=struct2cell(p);
0408     if ~isempty(param_name)
0409       for l=1:length(param_name)
0410         try
0411           value=eval(param_value{l});
0412         catch
0413           error('Values of this type are not supported for parameters set in equations.');
0414         end
0415         if isempty(spec.populations(i).parameters)
0416           spec.populations(i).parameters={param_name{l},value};
0417         elseif ~ismember(param_name{l},spec.populations(i).parameters(1:2:end))
0418           spec.populations(i).parameters{end+1}=param_name{l};
0419           spec.populations(i).parameters{end+1}=value;
0420         end
0421       end
0422     end    
0423     % TODO: remove support for MECH.PARAM from dsParseModelEquations,
0424     % because that returns MECH_PARAM, whereas we now support MECH.PARAM.
0425     
0426     % incorporate user-supplied parameters in pop equations if used in them
0427     % note: this is necessary for proper substitution when the master
0428     % equations are parsed in dsGenerateModel.
0429     if ~isempty(spec.populations(i).parameters)
0430       keys=spec.populations(i).parameters(1:2:end);
0431       vals=spec.populations(i).parameters(2:2:end);
0432       
0433       % add user-supplied params to pop equations if present in them
0434       % approach: look for populations.parameters in population.equations that are not explicitly
0435       % defined in population.equations and append their definition explicitly to pop.eqns
0436       eqn=spec.populations(i).equations;
0437      
0438       % get list of parameters/variables/functions in population equations
0439       words=unique(regexp(eqn,'[a-zA-Z]+\w*','match'));
0440       
0441       % find those in user-supplied parameters
0442       found_words=words(ismember(words,keys));
0443       if ~isempty(found_words)
0444         % set in population equations if not already defined there
0445         for ff=1:length(found_words)
0446           found_word=found_words{ff};
0447           precision=8; % number of digits allowed for user-supplied values
0448           found_value = toString(vals{strcmp(found_word,keys)},precision);
0449           
0450           % check if not explicitly set in population equations
0451           if isempty(regexp(eqn,[';\s*' found_word '\s*='],'once')) && ... % not in middle or at end
0452              isempty(regexp(eqn,['^' found_word '\s*='],'once')) % not at beginning
0453             % append and explicitly set in population equations
0454             if eqn(end)~=';', eqn(end+1)=';'; end % add semicolon if necessary
0455             eqn=[eqn sprintf(' %s=%s;',found_word,found_value)];
0456           else
0457             % update values in population equations
0458             old=regexp(eqn,[';\s*' found_word '\s*=\s*[\w\.'']+'],'match','once'); % in middle or at end
0459             if ~isempty(old)
0460 %               % remove semicolon for proper substitution using dsStrrep
0461 %               old=old(2:end);
0462               % replace value in middle or at end
0463               new=['; ' found_word '=' found_value];
0464 %               new=[' ' found_word '=' found_value];
0465             else
0466               % replace value at the beginning
0467               old=regexp(eqn,['^' found_word '\s*=\s*[\w\.'']+'],'match','once');
0468               new=[found_word '=' found_value];
0469             end
0470             eqn=strrep(eqn,old,new); % update value in equations
0471 %             eqn=dsStrrep(eqn,old,new); % update value in equations
0472           end
0473         end
0474         spec.populations(i).equations=eqn;
0475       end
0476     end
0477   end
0478   % expand mechanism list if any element is itself a list of mechanisms (eg, {'iCa','{CaBuffer,iCan}'} or '{CaBuffer,iCan}')
0479   spec.populations(i).mechanism_list=expand_list(spec.populations(i).mechanism_list, varargin{:});
0480 end
0481 
0482 % 2.0 standardize connections
0483 if ~isempty(spec.connections)
0484   % check for proper fields in connections substructure
0485   if ~isfield(spec.connections,'source')
0486     spec.connections(1).source=[];
0487   end
0488   
0489   if ~isfield(spec.connections,'target')
0490     spec.connections(1).target=[];
0491   end
0492   
0493   if ~isfield(spec.connections,'mechanism_list')
0494     spec.connections(1).mechanism_list=[];
0495   end
0496 
0497   if ~isfield(spec.connections,'mechanisms')
0498     spec.connections(1).mechanisms=[];
0499   end
0500 
0501   if ~isfield(spec.connections,'parameters')
0502     spec.connections(1).parameters={};
0503   end
0504 end
0505 for i=1:length(spec.connections)
0506   if isempty(spec.connections(i).source) && length(spec.populations)==1
0507     spec.connections(i).source=spec.populations(1).name;
0508     spec.connections(i).target=spec.populations(1).name;
0509   elseif isempty(spec.connections(i).source) && length(spec.populations)>1
0510     error('connection source and target populations must be specified in specification.connections when the model contains more than one population.');
0511   end
0512   
0513   % make mechanism list a cell array of mechanism names
0514   if ischar(spec.connections(i).mechanism_list)
0515     spec.connections(i).mechanism_list={spec.connections(i).mechanism_list};
0516   end
0517   
0518   % expand mechanism list if any element is itself a list of mechanisms (eg, {'AMPA','{GABAa,GABAb}'} or '{GABAa,GABAb}')
0519   spec.connections(i).mechanism_list=expand_list(spec.connections(i).mechanism_list, varargin{:});
0520   
0521   % parameter cell arrays
0522   if ~iscell(spec.connections(i).parameters)
0523     spec.connections(i).parameters={};
0524   end
0525 end
0526 
0527 % remove populations with size==0
0528 sizes=[spec.populations.size];
0529 if any(sizes==0)
0530   % find null populations
0531   null_pops=find(sizes==0);
0532   null_names={spec.populations(null_pops).name};
0533   
0534   % remove from .populations
0535   spec.populations(null_pops)=[];
0536   
0537   % remove from connections
0538   if ~isempty(spec.connections)
0539     sources={spec.connections.source};
0540     targets={spec.connections.target};
0541     null_conns=ismember(sources,null_names) | ismember(targets,null_names);
0542     spec.connections(null_conns)=[];
0543   end
0544 end
0545 
0546 % 3.0 sort fields
0547 % remove extra fields
0548 otherfields=setdiff(fieldnames(spec.populations),pop_field_order);
0549 spec.populations=rmfield(spec.populations,otherfields);
0550 
0551 % sort standardized fields
0552 spec.populations=orderfields(spec.populations,pop_field_order);
0553 if isstruct(spec.connections)
0554   otherfields=setdiff(fieldnames(spec.connections),con_field_order);
0555   spec.connections=rmfield(spec.connections,otherfields);
0556   spec.connections=orderfields(spec.connections,con_field_order);
0557 end
0558 spec=orderfields(spec,{'populations','connections','mechanisms'});
0559 
0560 % 4.0 standardize mechanism equations
0561 [~,files]=dsLocateModelFiles(spec);
0562 % read mechanism files
0563 for f=1:length(files)
0564   [~,name]=fileparts(files{f}); % name of mechanism
0565   if isempty(spec.mechanisms) || ~ismember(name,{spec.mechanisms.name})
0566     spec.mechanisms(end+1).name=name;
0567     spec.mechanisms(end).equations=read_mechanism_file(files{f});
0568   end
0569 end
0570 % make sure equations are all in a single string
0571 for m=1:length(spec.mechanisms)
0572   if iscellstr(spec.mechanisms(m).equations)
0573     % make sure each line ends with a semicolon and space
0574     eqn=spec.mechanisms(m).equations;
0575     idx=cellfun(@isempty,regexp(eqn,';$')); % lines that need semicolons
0576     eqn(idx)=cellfun(@(x)[x ';'],eqn(idx),'uni',0);
0577     idx=cellfun(@isempty,regexp(eqn,'\s$')); % lines that need a space
0578     eqn(idx)=cellfun(@(x)[x ' '],eqn(idx),'uni',0);
0579     % concatenate lines into a single string
0580     spec.mechanisms(m).equations=[eqn{:}];
0581   end
0582 end
0583 
0584 % 4.1 replace mechanism names by full file names
0585 % this is necessary so that regenerated models will use the same mechanism
0586 % files to recreate the model (e.g., when a cluster job simulates a
0587 % modified version of an original base model).
0588 % also: add mech equations to pop specification, and apply global population
0589 % parameters to pop-specific mechanism equations.
0590 
0591 fnames={};
0592 if ~isempty(files)
0593   for f=1:length(files)
0594     [~,name]=fileparts2(files{f});
0595     fnames{f}=name;
0596   end
0597 end
0598 if ~isempty(spec.mechanisms)
0599   mnames={spec.mechanisms.name};
0600 else
0601   mnames={};
0602 end
0603 fields={'populations','connections'};
0604 % update population and connection mechanism lists
0605 for f=1:length(fields)
0606   object=fields{f};
0607   for i=1:length(spec.(object))
0608     for j=1:length(spec.(object)(i).mechanism_list)
0609       mech=spec.(object)(i).mechanism_list{j};
0610       if ismember(mech,fnames)
0611         % replace mechanism names by full file names
0612         % this is necessary so that regenerated models will use the same mechanism
0613         % files to recreate the model (e.g., when a cluster job simulates a
0614         % modified version of an original base model).
0615         index=find(ismember(fnames,mech),1,'first');
0616         spec.(object)(i).mechanism_list{j}=files{index};
0617       end
0618       if ismember('@',mech) % remove linker alias
0619         mech=regexp(mech,'@','split');
0620         mech=mech{1};
0621       end
0622       if isfield(spec.(object),'mechanisms') && ~isempty(spec.(object)(i).mechanisms) && ismember(mech,{spec.(object)(i).mechanisms.name})
0623         % mechanism already defined for this population, nothing else to do
0624         continue;
0625       end
0626       if ismember(mech,mnames)
0627         % store mechanism equations in object-level specification
0628         index=find(ismember(mnames,mech),1,'first');
0629         if isempty(spec.(object)(i).mechanisms)
0630           spec.(object)(i).mechanisms=spec.mechanisms(index);
0631         elseif ~ismember(mech,{spec.(object)(i).mechanisms.name})
0632           % store if mechanism is not already defined at population-level
0633           spec.(object)(i).mechanisms(j)=spec.mechanisms(index);
0634         end
0635         % apply global population parameters to pop-specific mechanism equations
0636         if ~isempty(spec.(object)(i).parameters)
0637           % approach: set key=val for all keys in eqns
0638           eqns=spec.(object)(i).mechanisms(j).equations;
0639           keys=spec.(object)(i).parameters(1:2:end);
0640           vals=spec.(object)(i).parameters(2:2:end);
0641           % get list of parameters/variables/functions in population equations
0642           words=unique(regexp(eqns,'[a-zA-Z]+\w*','match'));
0643           % find words in user-supplied parameters (keys)
0644           found_words=words(ismember(words,keys));
0645           if ~isempty(found_words)
0646             for ff=1:length(found_words)
0647               found_word=found_words{ff};
0648               % new parameter assignment
0649               precision=8; % number of digits allowed for user-supplied values
0650               found_value=toString(vals{strcmp(found_word,keys)},precision);
0651               rep=sprintf(' %s=%s;',found_word,found_value);
0652               % replace old parameter assignment in the middle of equations
0653               pat=['([^\w]{1})' found_word '\s*=\s*\w+;']; % find in the middle
0654               eqns=regexprep(eqns,pat,['$1' rep]);
0655               % replace old parameter assignment at the beginning of equations
0656               pat=['^' found_word '\s*=\s*\w+;']; % find at the beginning
0657               eqns=regexprep(eqns,pat,rep);
0658             end
0659           end
0660           spec.(object)(i).mechanisms(j).equations=eqns;
0661         end        
0662       end
0663     end
0664   end
0665 end
0666 
0667 %% auto_gen_test_data_flag argout
0668 if options.auto_gen_test_data_flag
0669   argout = {spec}; % specific to this function
0670   
0671   dsUnitSaveAutoGenTestData(argin, argout);
0672 end
0673 
0674 end % main fn
0675 
0676 
0677 %% local fns
0678 function txt=read_mechanism_file(file)
0679   fid=fopen(file,'rt');
0680   % read all text
0681   txt=textscan(fid,'%s','Delimiter','\n');
0682   if ~isempty(txt)
0683     % remove empty lines
0684     txt=txt{1}(~cellfun(@isempty,txt{1}));
0685     % remove comments
0686     txt=txt(~cellfun(@isempty,regexp(txt,'^[^#%]')));
0687     txt=regexp(txt,'^[^%#]*','match');
0688     txt=[txt{:}];
0689     % remove leading/trailing white space
0690     txt=strtrim(txt);
0691     % make sure each line ends with a semicolon and space
0692     idx=cellfun(@isempty,regexp(txt,';$')); % lines that need semicolons
0693     txt(idx)=cellfun(@(x)[x ';'],txt(idx),'uni',0);
0694     idx=cellfun(@isempty,regexp(txt,'\s$')); % lines that need a space
0695     txt(idx)=cellfun(@(x)[x ' '],txt(idx),'uni',0);
0696     % concatenate lines into a single string
0697     txt=[txt{:}];
0698   end
0699   % close text file
0700   fclose(fid);
0701 end
0702 
0703 function list = expand_list(list, varargin)
0704   %% auto_gen_test_data_flag argin
0705   options = dsCheckOptions(varargin,{'auto_gen_test_data_flag',0,{0,1}},false);
0706   if options.auto_gen_test_data_flag
0707     varargs = varargin;
0708     varargs{find(strcmp(varargs, 'auto_gen_test_data_flag'))+1} = 0;
0709     varargs(end+1:end+2) = {'unit_test_flag',1};
0710     argin = [{list}, varargs];
0711   end
0712 
0713   % expand mechanism list if any element is itself a list of mechanisms (eg, {'AMPA','{GABAa,GABAb}'} or '{GABAa,GABAb}')
0714   if isempty(list)
0715     return;
0716   end
0717 
0718   if any(~cellfun(@isempty,regexp(list,'[{,}]+')))
0719     mechs={};
0720     for k=1:length(list)
0721       tmp=regexp(list{k},'\w+','match');
0722       mechs=cat(2,mechs,tmp{:});
0723     end
0724     list=mechs;
0725   end
0726 
0727   %% auto_gen_test_data_flag argout
0728   if options.auto_gen_test_data_flag
0729     argout = {list};
0730 
0731     dsUnitSaveAutoGenTestDataLocalFn(argin, argout); % localfn
0732   end
0733 
0734 end
0735 
0736 
0737 function spec = backward_compatibility(spec, varargin)
0738   % purpose: change name of fields from old to new convention
0739   % rename "nodes" or "entities" to "populations"
0740 
0741   %% auto_gen_test_data_flag argin
0742   options = dsCheckOptions(varargin,{'auto_gen_test_data_flag',0,{0,1}},false);
0743   if options.auto_gen_test_data_flag
0744     varargs = varargin;
0745     varargs{find(strcmp(varargs, 'auto_gen_test_data_flag'))+1} = 0;
0746     varargs(end+1:end+2) = {'unit_test_flag',1};
0747     argin = [{spec}, varargs];
0748   end
0749 
0750   if isfield(spec,'nodes')
0751     spec.populations=spec.nodes;
0752     spec=rmfield(spec,'nodes');
0753   end
0754 
0755   if isfield(spec,'cells')
0756     spec.populations=spec.cells;
0757     spec=rmfield(spec,'cells');
0758   end
0759 
0760   if isfield(spec,'entities')
0761     spec.populations=spec.entities;
0762     spec=rmfield(spec,'entities');
0763   end
0764 
0765   if isfield(spec,'pops')
0766     spec.populations=spec.pops;
0767     spec=rmfield(spec,'pops');
0768   end
0769 
0770   if isfield(spec,'cons')
0771     spec.connections=spec.cons;
0772     spec=rmfield(spec,'cons');
0773   end
0774   if isfield(spec,'links')
0775     spec.connections=spec.links;
0776     spec=rmfield(spec,'links');
0777   end
0778   if isfield(spec,'edges')
0779     spec.connections=spec.edges;
0780     spec=rmfield(spec,'edges');
0781   end
0782 
0783   if isfield(spec,'edges')
0784     spec.connections=spec.edges;
0785     spec=rmfield(spec,'edges');
0786   end
0787   
0788   if isfield(spec,'links')
0789     spec.connections=spec.links;
0790     spec=rmfield(spec,'links');
0791   end
0792   
0793   if isfield(spec,'comps')
0794     spec.compartments=spec.comps;
0795     spec=rmfield(spec,'comps');
0796   end
0797   
0798   if isfield(spec,'populations')
0799     % rename population "label" to "name"
0800     if isfield(spec.populations,'label')
0801       for i=1:length(spec.populations)
0802         spec.populations(i).name=spec.populations(i).label;
0803       end
0804       spec.populations=rmfield(spec.populations,'label');
0805     end
0806 
0807     % rename population "multiplicity" to "size"
0808     if isfield(spec.populations,'multiplicity')
0809       for i=1:length(spec.populations)
0810         spec.populations(i).size=spec.populations(i).multiplicity;
0811       end
0812       spec.populations=rmfield(spec.populations,'multiplicity');
0813     end
0814 
0815     % rename population "dynamics" to "equations"
0816     if isfield(spec.populations,'dynamics')
0817       for i=1:length(spec.populations)
0818         spec.populations(i).equations=spec.populations(i).dynamics;
0819       end
0820       spec.populations=rmfield(spec.populations,'dynamics');
0821     end
0822 
0823   end
0824 
0825   if isfield(spec,'connections') && isfield(spec.connections,'direction')
0826     for i=1:length(spec.connections)
0827       if ischar(spec.connections(i).direction)
0828         str=spec.connections(i).direction;
0829         if any(regexp(str,'->','once'))
0830           pops=regexp(str,'->','split');
0831           spec.connections(i).source=pops{1};
0832           spec.connections(i).target=pops{2};
0833         elseif any(regexp(str,'<-','once'))
0834           pops=regexp(str,'<-','split');
0835           spec.connections(i).source=pops{2};
0836           spec.connections(i).target=pops{1};
0837         end
0838       end
0839     end
0840     spec.connections=rmfield(spec.connections,'direction');
0841   end
0842 
0843   %% auto_gen_test_data_flag argout
0844   if options.auto_gen_test_data_flag
0845     argout = {spec};
0846 
0847     dsUnitSaveAutoGenTestDataLocalFn(argin, argout); % localfn
0848   end
0849 
0850 end