%**************************************************************************
% CLUSTER project: Curlometer method               Autor: Christian Maszl *
% Formulas: y0*Jav*(r1a x r1b) = dB1a*r1b-dB1b*r1a                        *
%           <div B>*|r1a*(r1b x r1c)| = |Sum(cyclic) dB1a*(r1b x r1c)|    *
% Unit System = [SI]                  Mail: e0025754@student.tuwien.ac.at *
%**************************************************************************

clear variables;
tic
%**************************************************************************
% Load configuration data like paths, constants, date, ...
%**************************************************************************

% Load configuration data and check the entries.
syscheck;
% Initial length guess for array preallocating
dimension=270000; 
% Preallocate space for the input data..
for i=1:1:4
    eval(['D',int2str(i),'=zeros(dimension,4);']);
end 
%**************************************************************************
%% Read data from Datafiles / [B]=nT [r]=km / Coord. system=GSE 
% [OK] Write a subroutine that detects the files. Style of the filenames: 
% [yymmdd]0000[type].[satellite].[t1]-[t2].[coordsys].[timeformat].[igmasci]
% Example: 0203040000fn.01.0700-1019.gse.hr.igmasc
filelist;
% Output processing status, remove in the real program
    text1=sprintf('%s \n','File list ready. Proceed in execution!');
          disp(text1);
% Ask User a few questions.
userinput;
% The "BIG" loop for all packages or a manual range entry
for lauf=startloop:1:endloop
% Store data in memory (binary datafiles=5x faster)
for i=1:1:4
    filename=[datapath,myfilelist{lauf,i}];
    eval(['D',int2str(i),'=load(filename,datatype);']);
end
% Output processing status, remove in the real program
    text1=sprintf('\n%s\n','Data loaded to RAM. Proceed in execution!');
          disp(text1);
% [Ok] Sort the data of the four satellites according to the date, time 
% and the satellite numbers
% [Ok] Write a function for the sort routine
[D1, D2, D3, D4] = sorting(D1,D2,D3,D4);
% Output processing status, remove in the real program
    text1=sprintf('%u%s\n',lauf,'. data-package sorted. Proceed in'...
                                ' execution!');
          disp(text1);
% [Ok] Split the input data into B1,B2,B3,B4 and r1,r2,r3,r4;
for i=1:1:4    
    eval(['B',int2str(i),'=D',int2str(i),'(:,2:4);']);
    eval(['r',int2str(i),'=D',int2str(i),'(:,5:7);']);
end
% Store the time information for the output-datafile (time accuracy 8E-6h)
    T = D1(:,1); 
% Clear dispensable variables
    clear dimension filename readin D1 D2 D3 D4 i;
%**************************************************************************
% Perform calculations (see tetrahedron.jpg for details about
% geometrical coherences).
%
% [Ok] Write a function faster execution 
[jav,divB] = curlB(B1,B2,B3,B4,r1,r2,r3,r4,y0);
% Output processing status, remove in the real program
    text1=sprintf('%u%s\n',lauf,'. data-package ready. Proceed in'...
                                ' execution!');
          disp(text1)
% Store data to harddisk
try
     datatodisk;
catch
     text1=['Could not write data to HD.The disk quota is either'...
           ' exceeded or you do not have write permission for the'...
           ' directory!'];
     disp(text1);
     disp('Please restart the program!');
     break;
end
     if (lauf~=endloop)
        disp('Data stored to harddisk. Proceed in execution!')
     else
        disp('Data stored to harddisk.')
	counter=endloop-startloop+1;
        text1=sprintf('\n%s %u%s\n','Computation of',counter, ...
                      ' data-packages ready. End of execution!');
              disp(text1)                
     end
end
t=toc;
% Show the computation time
t=num2str(sec2hm(t));
l=size(t,2);
if (l<=2)
    h='0';
    m=t;
else
    h=t(1,1:l-2);
    m=t(1,l-1:l);
end
    text1=sprintf('\n%s%s%s%s\n','Computation time: ',h,':',m);
          disp(text1);
% Free memory
clear variables;
