// !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
// ! This script is provided "as is" by Applied Maths.              !
// ! You are free to use and modify this script for your own needs. !
// ! Redistribution or reproduction of the script is prohibited.    !
// ! DISCLAIMER:                                                    !
// ! Improper use of scripts may corrupt your database.             !
// ! Running this script is entirely at your own responsibility.    !
// ! Applied Maths accepts no lialibility for any consequences      !
// ! resulting from its use.                                        !
// !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

integer i,j,ok,ok2,ps,ct,filenr,resol;
string homedir,dbasepath,filename,curfilename,st,st1,str_resol;
integer enzlist_count;
string enzlist_name[],enzlist_seq[];
integer enzlist_cutpos[];
string enzlist_all,enzlist_sel;
string seqtypelist,seqtypename,fprintlist,fprintname;
FILE fl;
DIALOG dlg;
string key;
integer lanecount;
CRFPRINT fp;
SEQUENCE sq;
TIME tm;
integer len,startpos,endpos,cutarray[];

integer enzsel_count,enzsel_idx[];
string primer1,primer2;
integer primer1_mismatch,primer2_mismatch;
string primer1_strmismatches,primer2_strmismatches;


//Check if at least some entries are selected
j=0;
for i=1 To DbGetEntryCount do
   if DbGetSel(DbGetEntryKey(i)) then j=j+1;
if j<=0 then {
   message("You should select some entries in the database");
   stop;
}

dbasepath=DbGetPath;
for i=1 to length(dbasepath) do
   if substring(dbasepath,i,i)="\" then j=i;
homedir=substring(dbasepath,1,j-1);


//Sequence type list
DlgAddText(dlg,"Sequence type:",20,20,120,15);
seqtypelist="";seqtypename="";
for i=1 to DbGetExperCount do
   if DbGetExperclass(i)="SEQ" then {
      seqtypelist=fprintlist+DbGetExperName(i)+"	";
      if length(seqtypename)=0 then seqtypename=DbGetExperName(i);
   }
DlgAddList(dlg,seqtypelist,seqtypename,20,40,120,110,"LIST");





//Fingerprint type list
DlgAddText(dlg,"Fingerprint type:",20,150,120,15);
fprintlist="";fprintname="";
for i=1 to DbGetExperCount do
   if DbGetExperclass(i)="FPR" then {
      fprintlist=fprintlist+DbGetExperName(i)+"	";
      if length(fprintname)=0 then fprintname=DbGetExperName(i);
   }
DlgAddList(dlg,fprintlist,fprintname,20,170,120,110,"LIST");






//Create list with restriction enzymes
filename=homedir+"\restr.txt";
if not(FileOpenRead(fl,filename)) then {
   message("Unable to find the restriction enzymes file~n'"+filename+"'");
   stop;
}
enzlist_count=0;
while not(FileIsEnd(fl)) do {
   st=FileRead(fl,1000);
   if length(st)>0 then {
      enzlist_count=enzlist_count+1;
      enzlist_name[enzlist_count]=splitstring(st,",");
      enzlist_seq[enzlist_count]=splitstring(st,",");
      enzlist_cutpos[enzlist_count]=val(splitstring(st,","));
   }
}
FileClose(fl);

DlgAddText(dlg,"Restriction enzyme(s)",170,20,110,15);
enzlist_all="";
for i=1 to enzlist_count do enzlist_all=enzlist_all+enzlist_name[i]+"	";
DlgAddList(dlg,enzlist_all,enzlist_sel,170,40,100,180,"MULTI");





//Create primer edit boxes

primer1="";
primer2="";

DlgAddText(dlg,"Primer 1 sequence",300,20,200,15);
DlgAddEdit(dlg,primer1,300,40,200,18);
DlgAddText(dlg,"Allow mismatches:",300,70,100,15);
primer1_strmismatches="1";
DlgAddEdit(dlg,primer1_strmismatches,420,70,30,18);

DlgAddText(dlg,"Primer 2 sequence",300,120,200,15);
DlgAddEdit(dlg,primer2,300,140,200,18);
DlgAddText(dlg,"Allow mismatches:",300,170,100,15);
primer2_strmismatches="1";
DlgAddEdit(dlg,primer2_strmismatches,420,170,30,18);



//Create file name box
filename="ARDRA";
DlgAddText(dlg,"Ouput file name:",170,230,80,15);
DlgAddEdit(dlg,filename,260,230,130,20);

//Create resolution box
str_resol="1000";
DlgAddText(dlg,"Curve resolution:",170,260,85,15);
DlgAddEdit(dlg,str_resol,260,260,50,20);


//prompt for dialog box
ok=0;
while not(ok) do {
   if not(DlgShow(dlg,"ARDRA profile generation",530,320)) then stop;
   ok=1;
   enzsel_count=0;
   while length(enzlist_sel)>0 do {
      st1=splitstring(enzlist_sel,"	");
      j=0;for i=1 to enzlist_count do if st1=enzlist_name[i] then j=i;
      if j>0 then { enzsel_count=enzsel_count+1;enzsel_idx[enzsel_count]=j; }
   }
   if enzsel_count<=0 then { message("You should select at least one restriction enzyme"); ok=0; }
   if enzsel_count> 5 then { message("You should select no more than 5 restriction enzymes"); ok=0; }
   
   if (length(primer1)<5) or (length(primer1)>50) or
      (length(primer2)<5) or (length(primer2)>50) then {
         Message("The primer sequences should be between 5 and 50 base pairs");
         ok=0;
      }

   primer1_mismatch=val(primer1_strmismatches);
   primer2_mismatch=val(primer2_strmismatches);
   resol=val(str_resol);
            
   if length(seqtypename)<=0 then { message("You should select a sequence type"); ok=0; }
   if length(fprintname)<=0 then { message("You should select a fingerprint name"); ok=0; }
   if length(filename)<=0 then { message("You should enter a file name"); ok=0; }
   if (resol<100) or (resol>15000) then { message("The resolution should be between 100 and 15000"); ok=0; }
   
}

//Checks that all entries have that sequence type
j=0;
for ct=1 to DbGetEntryCount do  if DbGetSel(DbGetEntryKey(ct)) then {
   key=DbGetEntryKey(ct);
   if length(DbGetEntryLink(key,seqtypename))<=0 then j=j+1;
}
if j>0 then {
   message("Some of the selected entries do not have the sequence type '"+seqtypename+"'");
   stop;
}

//Checks that no entries have already that fingerprint type
j=0;
for ct=1 to DbGetEntryCount do  if DbGetSel(DbGetEntryKey(ct)) then {
   key=DbGetEntryKey(ct);
   if length(DbGetEntryLink(key,fprintname))>0 then j=j+1;
}
if j>0 then {
   message("Some of the selected entries have already an experiment that belongs to '"+fprintname+"'");
   stop;
}

//Loops through all database entries & creates fingerprints
curfilename="";lanecount=0;filenr=0;
setbusy("Processing - 0");
for ct=1 to DbGetEntryCount do  if DbGetSel(DbGetEntryKey(ct)) then {
   key=DbGetEntryKey(ct);
   if floor(ct/10)=ct/10 then setbusy("Processing - "+str(ct,0,0));
   
//If necessary, create a new fingerprint file
   if (lanecount>=250) or length(curfilename)<=0 then {
      if lanecount>0 then {
         CrFprBandsToCurve(fp,resol);
         CrFprSave(fp);
         lanecount=0;
      }
      ok=0;
      while (not(ok)) and (filenr<2000) do {
         filenr=filenr+1;curfilename=filename+str(filenr,0,0);
         ok=1;
         st=DbGetPath+"\fprint\"+curfilename+".DEF";
         if TmGetFileTime(st,"ACCESS",tm) then ok=0;
         st=DbGetPath+"\images\"+curfilename+".CRV";
         if TmGetFileTime(st,"ACCESS",tm) then ok=0;
      }
      CrFprCreate(fp,curfilename,fprintname);lanecount=0;
   }

//Load & process a sequence
   SeqLoad(sq,key,seqtypename);
   SeqGet(sq,st);len=length(st);
   ok=0;
   if len>10 then ok=1;
   
//Determination of forward primer position
   j=val(primer1_strmismatches);
   startpos=SeqFind(sq,primer1,0,j);

//Determination of reverse primer position
   j=val(primer2_strmismatches);
   endpos=SeqFind(sq,primer2,0,j);
   if endpos>0 then endpos=endpos+length(primer2);
   
   if startpos<1 then ok=0;
   if endpos<1 then ok=0;
   if startpos>=endpos then ok=0;
   
   if ok then {
      SeqGet(sq,st);
      st=substring(st,startpos,endpos);
      SeqSet(sq,st);
      len=length(st);
   }

//Determination of enzyme cutting positions
   for i=1 to len do cutarray[i]=0;
   for i=1 to enzsel_count do {
      j=enzsel_idx[i];st1=enzlist_seq[j];
      ok2=1;ps=0;
      while ok2 do {
         ps=SeqFind(sq,st1,ps+1,0);ok2=0;
         if ps>0 then {
            ok2=1;
            ps=ps+enzlist_cutpos[j];
            cutarray[ps]=1;
         }
      }
   }
   
//Determination of fragments
   if ok then {
      lanecount=CrFprAddLane(fp,fprintname,key);
      i=1;
      while i<len do {
         j=i+1;
         while (j<len) and (cutarray[j]=0) do j=j+1;
         CrFprAddBand(fp,lanecount,j-i,255,2);
         i=j+1;
      }
   }
   
      
}


//If necessary, saves the last fingerprint file
if lanecount>0 then {
   CrFprBandsToCurve(fp,1000);
   CrFprSave(fp);
   lanecount=0;
}


setbusy("");