/home/peifer/diffit/modules/readData.cc

Go to the documentation of this file.
00001 
00002 
00004 #include <stdlib.h>
00005 #include <stdio.h>
00006 #include <iostream>
00007 #include <fstream>
00008 #include <ctype.h>
00009 #include <math.h>
00010 #include <string.h>
00011 
00012 #include "../nr.h"
00013 #include "../def.h"
00014 
00015 using namespace std;
00016 
00017 void readData(GlobExp *ex,Glob *globs,long expNr)
00018 {
00019   int lineNr=0;
00020   long k,l,index;
00021   long nLines=0,maxLeng=0;
00022   long nobsmax=ex[expNr].nobs;
00023   char c,*line,*arg;
00024   double *x;
00025   double **y;
00026   double **sig;
00027 
00028   ifstream in;
00029 
00030   FILE *fp;
00031 
00032   fp=fopen(ex[expNr].fileName,"r");
00033   if (fp==NULL) {
00034     cerr << "Error opening data file " << ex[expNr].fileName << "\n";
00035     exit(1);
00036   }
00037   
00038   //read total number of lines
00039   while((c=fgetc(fp))!=EOF)
00040         {
00041           for(k=1;c!='\n';k++)
00042             {
00043               if(k>maxLeng)
00044                 maxLeng=k;
00045               c=fgetc(fp);
00046             }
00047           nLines++;
00048         }
00049 #ifdef DEBUGREADDATA
00050   *dbg << "Total number of lines in " <<  ex[expNr].fileName << " : " << nLines << endl;
00051   *dbg << "Maximal length of line in " << ex[expNr].fileName << " : " << maxLeng << endl;
00052 #endif
00053   fclose(fp);
00054 
00055   in.open(ex[expNr].fileName);
00056 
00057   //allocate memory
00058   line=new char[maxLeng+2];
00059   arg=new char[maxLeng+2];
00060   x=dvector(1,nLines);
00061   y=dmatrix(1,nLines,1,ex[expNr].nobs);
00062   sig=dmatrix(1,nLines,1,ex[expNr].nobs);
00063 
00064   //extracting the data
00065   while(!in.eof())
00066         {
00067           in.getline(line,maxLeng+2,'\n');
00068           if(isspace(line[0])==FALSE && line[0]!= '#' && line[0] != '\0' )
00069             {
00070               lineNr++;
00071 #ifdef DEBUGREADDATA
00072               *dbg << "Extracting data line no. : " << lineNr << endl;
00073 #endif
00074               l=0;
00075               index=0;
00076               for(k=0;k<=strlen(line);k++)
00077                 {
00078                   if(isspace(line[k])!=FALSE)
00079                     {
00080                       arg[l]='\0';
00081                       l=0;
00082                       if(index==0)
00083                         x[lineNr]=atof(arg);
00084                       else if(index <=2*ex[expNr].nobs && (index % 2)==1)
00085                         y[lineNr][(index+1)/2]=atof(arg);
00086                       else if(index <=2*ex[expNr].nobs && (index % 2)==0)
00087                         sig[lineNr][index/2]=atof(arg);
00088                       index++;
00089                     }
00090                   else
00091                     {
00092                       arg[l]=line[k];
00093                       l++;
00094                     }
00095                 }
00096               arg[l]='\0';      
00097               if(index==0)
00098                 x[lineNr]=atof(arg);
00099               else if(index <=2*ex[expNr].nobs && (index % 2)==1)
00100                 y[lineNr][(index+1)/2]=atof(arg);
00101               else if(index <=2*ex[expNr].nobs && (index % 2)==0)
00102                 sig[lineNr][index/2]=atof(arg);
00103               
00104               //determ. max. complete obs.
00105               if((index % 2)==1)
00106                 index--;
00107               if(index==0)
00108                 {
00109                   cerr << "Data line " << lineNr << " contains no data." << endl;
00110                   exit(1);
00111                 }
00112               if(index/2 <=  nobsmax)
00113                 nobsmax=index/2;
00114             }
00115                     
00116         }
00117   in.close();
00118 
00119   ex[expNr].nobs=nobsmax;
00120   ex[expNr].nMeasure=lineNr;
00121   if(ex[expNr].nMeasure==0)
00122     {
00123       cerr << "File " << ex[expNr].fileName << " contains no data." << endl;
00124       exit(1);
00125     }
00126 
00127   //check if the time is in ascending order
00128     for(k=2;k<=ex[expNr].nMeasure;k++)
00129       {
00130         if(x[k-1] >= x[k])
00131           {
00132             cerr << "Data not in temporal ascending order at dataline " << k << "\n";
00133             exit(1);
00134           }
00135       }
00136     //fit range stuff
00137     if(ex[expNr].fitstart == -1e99)
00138       ex[expNr].fitstart=x[1];
00139     if(ex[expNr].fitend ==  1e99)
00140       ex[expNr].fitend=x[ex[expNr].nMeasure];
00141     
00142     //for memory allocation
00143     index=1;
00144     for(k=1;k<=ex[expNr].nMeasure;k++)
00145       {
00146         if((x[k] >= ex[expNr].fitstart) && (x[k] <= ex[expNr].fitend))
00147           { 
00148             index++;
00149           }
00150       }
00151     index--; 
00152   
00153     //the actual amount of data
00154     ex[expNr].nMeasure=index;
00155 
00156     //memory allocation
00157     ex[expNr].xMeasure=dvector(1,ex[expNr].nMeasure);
00158     ex[expNr].yMeasure=dmatrix(1,ex[expNr].nMeasure,1,ex[expNr].nobs);
00159     ex[expNr].sigma=dmatrix(1,ex[expNr].nMeasure,1,ex[expNr].nobs);
00160 
00161     // copy data to data structure
00162     index=1;
00163     for(k=1;k<=lineNr;k++)
00164       {
00165         if((x[k] >= ex[expNr].fitstart) && (x[k] <= ex[expNr].fitend))
00166           {
00167             ex[expNr].xMeasure[index]=x[k];
00168             for(l=1;l<=ex[expNr].nobs;l++)
00169               {
00170                 ex[expNr].yMeasure[index][l]=y[k][l];
00171                 ex[expNr].sigma[index][l]=sig[k][l];
00172               }
00173             index++;
00174           }
00175       }
00176     index--;
00177 
00178     //setting some variables
00179     ex[expNr].firstMeasure=1;
00180     ex[expNr].lastMeasure=ex[expNr].nMeasure;
00181     
00182     //free memory
00183     free_dvector(x,1,nLines);
00184     free_dmatrix(y,1,nLines,1,ex[expNr].nobs);
00185     free_dmatrix(sig,1,nLines,1,ex[expNr].nobs);
00186 }

Generated on Mon Jan 29 17:09:13 2007 for Diffit by  doxygen 1.4.6