src/MEDMEM/MEDMEM_AsciiFieldDriver.hxx

   1 #ifndef ASCII_FIELD_DRIVER_HXX
   2 #define ASCII_FIELD_DRIVER_HXX
   3
   4 #include "MEDMEM_Mesh.hxx"
   5 #include "MEDMEM_Support.hxx"
   6 #include "MEDMEM_GenDriver.hxx"
   7 #include "MEDMEM_Exception.hxx"
   8 #include "MEDMEM_Unit.hxx"
   9
  10 #include <list>
  11 #include <string>
  12 #include <ctype.h>
  13 #include <iomanip>
  14 #include <stdlib.h>
  15 #include <string.h>
  16 #include <fstream.h>
  17
  18 #define PRECISION_IN_ASCII_FILE 10
  19 #define PRECISION_IN_COMPARE 1e-10
  20 #define SPACE_BETWEEN_NBS 19
  21
  22 namespace MEDMEM {
  23
  24   template<class T>
  25   class FIELD;
  26
  27   template<int N,unsigned int CODE>
  28   void fill(double *a, const double *b)
  29   {
  30     a[N]=b[CODE & 0x3 ];
  31     fill<N-1, (CODE>>2) >(a,b);
  32   }
  33
  34   template<int N>
  35   bool compare(const double* a, const double* b)
  36   {
  37     double sign=b[N]<0 ? -1 : 1;
  38     if(a[N]<b[N]*(1-sign*PRECISION_IN_COMPARE))
  39       return true;
  40     if(a[N]>b[N]*(1+sign*PRECISION_IN_COMPARE))
  41       return false;
  42     return compare<N-1>(a,b);
  43   }
  44
  45   template<>
  46   void fill<-1,0x3>(double *a, const double *b);
  47
  48   template<>
  49   bool compare<-1>(const double *a, const double *b);
  50
  51   template <class T, int SPACEDIMENSION, unsigned int SORTSTRATEGY>
  52   class SDForSorting
  53   {
  54   private:
  55     double _coords[SPACEDIMENSION];
  56     T *_components;
  57     int _nbComponents;
  58   public:
  59     SDForSorting(const double *coords, const T* comp, int nbComponents);
  60     SDForSorting(const SDForSorting& other);
  61     ~SDForSorting();
  62     bool operator< (const SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY>& other) const;
  63     void writeLine(ofstream& file) const;
  64   };
  65
  66   template <class T>
  67   class ASCII_FIELD_DRIVER : public GENDRIVER
  68   {
  69   private:
  70     MESH          *_mesh;
  71     SUPPORT       *_support;
  72     mutable FIELD<T> *     _ptrField;
  73     std::string         _fileName;
  74     mutable ofstream             _file;
  75     unsigned int         _code;
  76     MED_EN::med_sort_direc _direc;
  77     int                  _nbComponents;
  78     int                  _spaceDimension;
  79     //static int           _nbComponentsForCpyInfo;
  80   public:
  81     ASCII_FIELD_DRIVER(const string & fileName, FIELD<T> * ptrField,
  82                        MED_EN::med_sort_direc direction=MED_EN::ASCENDING,
  83                        const char *priority="");
  84     ASCII_FIELD_DRIVER(const ASCII_FIELD_DRIVER<T>& other);
  85     void open() throw (MEDEXCEPTION);
  86     void close();
  87     void read ( void ) throw (MEDEXCEPTION);
  88     void write( void ) const throw (MEDEXCEPTION);
  89     GENDRIVER* copy() const;
  90   private:
  91     void buildIntroduction() const;
  92     template<int SPACEDIMENSION, unsigned int SORTSTRATEGY>
  93     void sortAndWrite() const;
  94     //template<int SPACEDIMENSION, unsigned int SORTSTRATEGY>//, std::list< SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY > > >
  95     //static void copyInfo(const double *a,T *b);
  96     //static void copyInfo2(const double *,T *);
  97   };
  98 }
  99
 100 #include "MEDMEM_Field.hxx"
 101
 102 namespace MEDMEM {
 103
 104   template <class T, int SPACEDIMENSION, unsigned int SORTSTRATEGY>
 105   SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY>::SDForSorting(const double *coords, const T* comp, int nbComponents):_nbComponents(nbComponents)
 106   {
 107     fill<SPACEDIMENSION-1,SORTSTRATEGY>(_coords,coords);
 108     _components=new T[_nbComponents];
 109     memcpy(_components,comp,sizeof(T)*_nbComponents);
 110   }
 111
 112   template <class T, int SPACEDIMENSION, unsigned int SORTSTRATEGY>
 113   SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY>::SDForSorting(const SDForSorting& other):_nbComponents(other._nbComponents)
 114   {
 115     memcpy(_coords,other._coords,sizeof(double)*SPACEDIMENSION);
 116     _components=new T[_nbComponents];
 117     memcpy(_components,other._components,sizeof(T)*_nbComponents);
 118   }
 119
 120   template <class T, int SPACEDIMENSION, unsigned int SORTSTRATEGY>
 121   SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY>::~SDForSorting()
 122   {
 123     delete [] _components;
 124   }
 125
 126   template <class T, int SPACEDIMENSION, unsigned int SORTSTRATEGY>
 127   bool SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY>::operator< (const SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY>& other) const
 128   {
 129     return compare<SPACEDIMENSION-1>(_coords,other._coords);
 130   }
 131
 132   template <class T, int SPACEDIMENSION, unsigned int SORTSTRATEGY>
 133   void SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY>::writeLine(ofstream& file) const
 134   {
 135     int i;
 136     double temp[SPACEDIMENSION];
 137     fill<SPACEDIMENSION-1,SORTSTRATEGY>(temp,_coords);
 138     for(i=0;i<SPACEDIMENSION;i++)
 139       file << setw(SPACE_BETWEEN_NBS) << temp[i];
 140     for(i=0;i<_nbComponents;i++)
 141       file << setw(SPACE_BETWEEN_NBS) << _components[i];
 142     file << endl;
 143   }
 144
 145   template <class T>
 146   ASCII_FIELD_DRIVER<T>::ASCII_FIELD_DRIVER(const string & fileName, FIELD<T> * ptrField,
 147                                             MED_EN::med_sort_direc direction,const char *priority)
 148     :GENDRIVER(fileName,MED_EN::MED_ECRI),_ptrField(ptrField),_fileName(fileName),_direc(direction)
 149   {
 150     _nbComponents=_ptrField->getNumberOfComponents();
 151     if(_nbComponents<=0)
 152       throw MEDEXCEPTION("ASCII_FIELD_DRIVER : No components in FIELD<T>");
 153     _support=(SUPPORT *)_ptrField->getSupport();
 154     _mesh=(MESH *)_support->getMesh();
 155     _spaceDimension=_mesh->getSpaceDimension();
 156     _code=3;
 157     int i;
 158     if(priority[0]=='\0')
 159       for(i=_spaceDimension-1;i>=0;i--)
 160         {
 161           _code<<=2;
 162           _code+=i;
 163         }
 164     else
 165       {
 166         if(_spaceDimension!=strlen(priority))
 167           throw MEDEXCEPTION("ASCII_FIELD_DRIVER : Coordinate priority invalid with spaceDim");
 168         for(i=_spaceDimension-1;i>=0;i--)
 169           {
 170             char c=toupper(priority[i]);
 171             if(int(c-'X')>(_spaceDimension-1))
 172               throw MEDEXCEPTION("ASCII_FIELD_DRIVER : Invalid priority definition");
 173             _code<<=2;
 174             _code+=c-'X';
 175           }
 176       }
 177   }
 178
 179   template <class T>
 180   ASCII_FIELD_DRIVER<T>::ASCII_FIELD_DRIVER(const ASCII_FIELD_DRIVER<T>& other)
 181     :_ptrField(other._ptrField),_fileName(other._fileName),_direc(other._direc),_mesh(other._mesh),_nbComponents(other._nbComponents),
 182      _code(other._code),_spaceDimension(other._spaceDimension),_support(other._support)
 183   {
 184   }
 185
 186   template <class T>
 187   void ASCII_FIELD_DRIVER<T>::open() throw (MEDEXCEPTION)
 188   {
 189     _file.open(_fileName.c_str(),ofstream::out | ofstream::app);
 190   }
 191
 192   template <class T>
 193   void ASCII_FIELD_DRIVER<T>::close()
 194   {
 195     _file.close();
 196   }
 197
 198   template <class T>
 199   void ASCII_FIELD_DRIVER<T>::read ( void ) throw (MEDEXCEPTION)
 200   {
 201     throw MEDEXCEPTION("ASCII_FIELD_DRIVER::read : Can't read with a WRONLY driver !");
 202   }
 203
 204   template <class T>
 205   GENDRIVER* ASCII_FIELD_DRIVER<T>::copy() const
 206   {
 207     return new ASCII_FIELD_DRIVER<T>(*this);
 208   }
 209
 210   template <class T>
 211   void ASCII_FIELD_DRIVER<T>::write( void ) const throw (MEDEXCEPTION)
 212   {
 213     buildIntroduction();
 214     switch(_spaceDimension)
 215       {
 216       case 2:
 217         {
 218           switch(_code)
 219             {
 220             case 52: //XY
 221               {
 222                 sortAndWrite<2,52>();
 223                 break;
 224               }
 225             case 49: //YX
 226               {
 227                 sortAndWrite<2,49>();
 228                 break;
 229               }
 230             default:
 231               MEDEXCEPTION("ASCII_FIELD_DRIVER : Invalid priority definition");
 232             }
 233           break;
 234         }
 235         case 3:
 236         {
 237           switch(_code)
 238             {
 239             case 228: //XYZ
 240               {
 241                 sortAndWrite<3,228>();
 242                 break;
 243               }
 244             case 216: //XZY
 245               {
 246                 sortAndWrite<3,216>();
 247                 break;
 248               }
 249             case 225://YXZ
 250               {
 251                 sortAndWrite<3,225>();
 252                 break;
 253               }
 254             case 201://YZX
 255               {
 256                 sortAndWrite<3,201>();
 257                 break;
 258               }
 259             case 210://ZXY
 260               {
 261                 sortAndWrite<3,210>();
 262                 break;
 263               }
 264             case 198://ZYX
 265               {
 266                 sortAndWrite<3,198>();
 267                 break;
 268               }
 269              default:
 270               MEDEXCEPTION("ASCII_FIELD_DRIVER : Invalid priority definition");
 271             }
 272           break;
 273           }
 274       default:
 275         MEDEXCEPTION("ASCII_FIELD_DRIVER : Invalid space dimension must be 2 or 3");
 276         }
 277   }
 278
 279   template <class T>
 280   void ASCII_FIELD_DRIVER<T>::buildIntroduction() const
 281   {
 282     int i;
 283     _file  << setiosflags(ios::scientific);
 284     _file << "#TITLE: table " << _ptrField->getName() << " TIME: " << _ptrField->getTime() << " IT: " << _ptrField->getIterationNumber() << endl;
 285     _file << "#COLUMN_TITLES: ";
 286     for(i=0;i<_spaceDimension;i++)
 287       _file << char('X'+i) << " | ";
 288     const std::string *compoNames=_ptrField->getComponentsNames();
 289     for(i=0;i<_nbComponents;i++)
 290       {
 291         if(!compoNames)
 292           _file << compoNames[i];
 293         else
 294           _file << "None";
 295         if(i<_nbComponents-1)
 296           _file << " | ";
 297       }
 298     _file << endl;
 299     _file << "#COLUMN_UNITS: ";
 300     compoNames=_mesh->getCoordinateptr()->getCoordinatesUnits();
 301     for(i=0;i<_spaceDimension;i++)
 302       {
 303         if(!compoNames)
 304           _file << compoNames[i];
 305         else
 306           _file << "None";
 307         _file << " | ";
 308       }
 309     const UNIT *compoUnits=_ptrField->getComponentsUnits();
 310     for(i=0;i<_nbComponents;i++)
 311       {
 312         if(!compoUnits)
 313           _file << compoUnits[i].getName();
 314         else
 315           _file << "None";
 316         if(i<_nbComponents-1)
 317           _file << " | ";
 318       }
 319     _file << endl;
 320   }
 321
 322   template <class T>
 323   template<int SPACEDIMENSION, unsigned int SORTSTRATEGY>
 324   void ASCII_FIELD_DRIVER<T>::sortAndWrite() const
 325   {
 326     int i,j;
 327     int numberOfValues=_ptrField->getNumberOfValues();
 328     std::list< SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY > > li;
 329     const double * coord;
 330     FIELD<double> * barycenterField=0;
 331     double * xyz[SPACEDIMENSION];
 332     bool deallocateXyz=false;
 333     if(_support->getEntity()==MED_EN::MED_NODE)
 334       {
 335         if (_support->isOnAllElements())
 336           {
 337             coord=_mesh->getCoordinates(MED_EN::MED_NO_INTERLACE);
 338             for(i=0; i<SPACEDIMENSION; i++)
 339               xyz[i]=(double *)coord+i*numberOfValues;
 340           }
 341         else
 342           {
 343             coord = _mesh->getCoordinates(MED_EN::MED_FULL_INTERLACE);
 344             const int * nodesNumber=_support->getNumber(MED_EN::MED_ALL_ELEMENTS);
 345             for(i=0; i<SPACEDIMENSION; i++)
 346               xyz[i]=new double[numberOfValues];
 347             deallocateXyz=true;
 348             for(i=0;i<numberOfValues;i++)
 349               {
 350                 for(j=0;j<SPACEDIMENSION;j++)
 351                   xyz[j][i]=coord[(nodesNumber[i]-1)*SPACEDIMENSION+j];
 352               }
 353           }
 354       }
 355     else
 356       {
 357         barycenterField = _mesh->getBarycenter(_support);
 358         coord=barycenterField->getValue(MED_EN::MED_NO_INTERLACE);
 359         for(i=0; i<SPACEDIMENSION; i++)
 360           xyz[i]=(double *)(coord+i*numberOfValues);
 361       }
 362     T* valsToSet=(T*)_ptrField->getValue(MED_EN::MED_FULL_INTERLACE);
 363     double temp[SPACEDIMENSION];
 364     for(i=0;i<numberOfValues;i++)
 365       {
 366         for(j=0;j<SPACEDIMENSION;j++)
 367           temp[j]=xyz[j][i];
 368         li.push_back(SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY>(temp,valsToSet+i*_nbComponents,_nbComponents));
 369       }
 370     if(barycenterField)
 371       delete barycenterField;
 372     if(deallocateXyz)
 373       for(j=0;j<SPACEDIMENSION;j++)
 374         delete [] xyz[j];
 375     li.sort();
 376     _file << setprecision(PRECISION_IN_ASCII_FILE);
 377     if(_direc==MED_EN::ASCENDING)
 378       {
 379         typename std::list< SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY > >::iterator iter;
 380         for(iter=li.begin();iter!=li.end();iter++)
 381           (*iter).writeLine(_file);
 382         _file << endl;
 383       }
 384     else if(_direc==MED_EN::DESCENDING)
 385       {
 386         typename std::list< SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY > >::reverse_iterator iter;
 387         for(iter=li.rbegin();iter!=li.rend();iter++)
 388           (*iter).writeLine(_file);
 389         _file << endl;
 390       }
 391     else
 392       MEDEXCEPTION("ASCII_FIELD_DRIVER : Invalid sort direction");
 393   }
 394   //{
 395     //_nbComponentsForCpyInfo=_nbComponents;
 396     //_ptrField->fillFromAnalytic <TEST<T>::copyInfo3> ();
 397     //_ptrField->fillFromAnalytic<  ASCII_FIELD_DRIVER<T>::copyInfo<SPACEDIMENSION,SORTSTRATEGY> > ();
 398     //li.sort();
 399     //typename std::list< SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY > >::iterator iter;
 400     //_file << setiosflags(ios::scientific) << setprecision(PRECISION_IN_ASCII_FILE);
 401     //for(iter=li.begin();iter!=li.end();iter++)
 402     //  {
 403     //(*iter).writeLine(_file);
 404     // }
 405
 406 //   template <class T>
 407 //   template<int SPACEDIMENSION, unsigned int SORTSTRATEGY>//, std::list< SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY> > lis>
 408 //   void ASCII_FIELD_DRIVER<T>::copyInfo(const double *a,T *b)
 409 //   {
 410 //     //lis.push_back(SDForSorting<T,SPACEDIMENSION,SORTSTRATEGY>(a,b,_nbComponentsForCpyInfo));
 411 //   }
 412
 413 //   template <class T>
 414 //   int ASCII_FIELD_DRIVER<T>::_nbComponentsForCpyInfo=0;
 415 }
 416
 417 #endif