%Modelbrick Assembler, version 1.0, created: August 22, 2008
%Developed by the 2008 iGEM team from the University of Groningen
%Author: Sander Onur (contact: a 'dot' r 'dot' onur 'at' gmail 'dot' com) 
%
%This script is meant to be self-explanatory, for any questions regarding its
%workings, please contact: igemgroningen@googlegroups.com or contact the
%author directly.
%
%Before running, make sure that your library file is entitled 'partslibrary.mat' and that
%it is located in the same directory as this file. Note that partnames are
%CaSe-CeNsItIvE.
%
%This script is meant to be a proof of concept. It has only been tested on the University of Groningen 2008 iGEM
%project in combination with the 'partslibrary.mat' file as it is available on 2008.igem.org/Teams:Groningen (October 28, 2008).
%
home

sbiodesktop

load partslibrary

home

fprintf('This is the iGEM Groningen Modelbrick Assembler, version 1.0, created: August 22, 2008 \n')
fprintf('For help on using this product visit 2008.igem.org/Teams:Groningen \n')
fprintf('If you have any questions, remarks or suggestions for improving the product, please contact igemgroningen@googlegroups.com \n')

fprintf('\r\r');
%Choose a name
name = input('Choose a name for your model: ','s');
fprintf('\r');
home
%Choose a chassis to embed the parts in
chassis = input('Please type the name of the chassis here: ');
fprintf('\r');
home

i=0;
ans1 = 'y';

while ans1 == 'y'
    %Add parts to the model, partnames are case-sensitive
    ans2 = 'y';
    j=0;
    fprintf('\r');
    ans1 = input('Would you like to add a(nother) transcription region? (y=yes, n=no) ','s');
    fprintf('\r');
    home
    if ans1 == 'y'
        i=i+1;
        while ans2 =='y'
            j=j+1;
            plist{i,j} = input('Please type the name of the part you want to add: ');
            fprintf('\r');
            home
            ans2 = input('Would you like to add another part to this transcription region? (y=yes, n=no) ','s');
            fprintf('\r');
            home
        end
        cntparts(i)=j;
    elseif ans1 == 'n'
        fprintf('Assembling device, please wait... \n \n');
        %Create new model
        newmodel = sbiomodel(name);
        %Use chosen chassis to define a compartment
        copyobj(chassis.compartments(1), newmodel);

        %Rename compartments to match chassis
        for i=1:size(plist,1)
            for j=1:size(plist,2)
                plist{i,j}.compartments(1).name = newmodel.compartments(1).name;
            end
        end

        sname{1} = regexp(newmodel.species(1).name, '\_', 'split');
        sname{2} = regexp(newmodel.species(2).name, '\_', 'split');
        
        %Copy components to chassis
        scnt=1;
        rcnt=0;
        rcnt2=0;
        for i=1:size(plist,1)
            
            
            for j=1:cntparts(i)
                rcnt3=rcnt2+1;
                rcnt2=rcnt;
                for k=1:length(plist{i,j}.reactions)
                    rcnt=rcnt+1;
                    copyobj(plist{i,j}.reactions(k), newmodel);
                    newmodel.reactions(rcnt).name = strcat(num2str(i),'_',num2str(j),'_',newmodel.reactions(rcnt).name);
                    rname{rcnt} = regexp(newmodel.reaction(rcnt).name, '\_', 'split');
                    
                end
                for k=scnt+2:length(newmodel.species)
                    newmodel.species(k).name = strcat(num2str(i),'_',num2str(j),'_',newmodel.species(k).name);
                    sname{k} = regexp(newmodel.species(k).name, '\_', 'split');
                    if strcmp(sname{k}{3},'pops-out')==1
                        for l=rcnt3:rcnt2
                            if strcmp(rname{l}{3},'pops-out')==1
                                addproduct(newmodel.reactions(l),newmodel.species(k));
                            end
                        end
                    end
                    if strcmp(sname{k}{3},'rips-out')==1
                        for l=rcnt3:rcnt2
                            if strcmp(rname{l}{3},'rips-out')==1
                                addproduct(newmodel.reactions(l),newmodel.species(k));
                            end
                        end
                    end
                    scnt=scnt+1;
                end
            end
        end

        %Add pools
        p=1;
        pools={};
        for i=1:size(plist,1)
            for j=1:cntparts(i)
                tags = regexp(plist{i,j}.name, '\_', 'split');
                if strcmp(tags(2), 'pr')==1
                    nops= str2num(tags{3});
                    for k=4:3+nops
                        pools{p}=eval(tags{k});
                        p=p+1;
                    end
                end
                if strcmp(tags(2), 'cr')==1
                    pools{p}=eval(tags{3});
                    p=p+1;
                end
            end
        end
        p=1;

        poolnm{1}='a';
        for i=1:length(pools)
                if strcmp(pools{i}.name, poolnm) == zeros(1,length(poolnm))
                    pools{i}.compartments(1).name = newmodel.compartments(1).name;
                    for k=1:length(pools{i}.reactions)
                        copyobj(pools{i}.reactions(k), newmodel);
                        rcnt=rcnt+1;
                        rname{rcnt} = regexp(pools{i}.reactions(k).name, '\_', 'split');
                        
                    end
                    poolnm{p}=pools{i}.name;
                    for l=1:length(pools{i}.species)
                        sname{scnt+2}=regexp(pools{i}.species(l).name, '\_', 'split');
                        scnt=scnt+1;
                    end
                    p=p+1;
                end
        end
     
        %Make connections
        poolnm{p}='rnap';
        poolnm{p+1}='rib';
         for i=1:length(newmodel.reactions)
            if length(rname{i})>=4
            if strcmp(rname{i}{4},'in')==1
               for j=1:length(sname)
                   if length(sname{j})==2
                      if strcmp(strcat(sname{j}{1},sname{j}{2}),strcat(rname{i}{3},rname{i}{4}))==1
                         addreactant(newmodel.reactions(i),strcat(sname{j}{1},'_in'));
                      elseif (strcmp(sname{j}{1},rname{i}{3})==1) && (strcmp(sname{j}{2},'inout')==1)
                         addreactant(newmodel.reactions(i),strcat(sname{j}{1},'_inout'));
                      end
                   end
                end
            end
            if strcmp(rname{i}{4},'out')==1
                for j=1:length(newmodel.species)
                    if length(sname{j})==2
                    if strcmp(strcat(sname{j}{1},sname{j}{2}),strcat(rname{i}{3},rname{i}{4}))==1
                        addproduct(newmodel.reactions(i),strcat(sname{j}{1},'_out'));
                    elseif (strcmp(sname{j}{1},rname{i}{3})==1) && (strcmp(sname{j}{2},'inout')==1)
                        addproduct(newmodel.reactions(i),strcat(sname{j}{1},'_inout'));
                    end
                    end
                end
            end
            end
         end

        %Repair reactions
        for i=1:length(newmodel.reactions)
            newmodel.reactions(i).kineticlaw(1).kineticlawname='MassAction';
            for j=1:length(newmodel.reactions(i).kineticlaw(1).Parameters)
                pars{j}=newmodel.reactions(i).kineticlaw(1).Parameters(j).Name;
            end
            set (newmodel.reactions(i).kineticlaw(1), 'ParameterVariablenames', pars);
            clear pars;
        end
            
        %Save and finish
        newmodel
        sbiodesktop(newmodel)
        fprintf('Your device has been exported to the workspace as ''newmodel'' \r \r')
    else
        fprintf('unrecognized input \n');
        break
    end
end

%Clear redundant variables
clear plist i j k name ans1 ans2 chassis cntparts nops p poolnm pools tags rname sname rcnt rcnt2 rcnt3 scnt l