FFmaGenericLoader.hpp 24.8 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
// ===================================================================================
// Copyright ScalFmm 2011 INRIA, Olivier Coulaud, Berenger Bramas, Matthias Messner
// olivier.coulaud@inria.fr, berenger.bramas@inria.fr
// This software is a computer program whose purpose is to compute the FMM.
//
// This software is governed by the CeCILL-C and LGPL licenses and
// abiding by the rules of distribution of free software.  
// 
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU General Public and CeCILL-C Licenses for more details.
// "http://www.cecill.info". 
// "http://www.gnu.org/licenses".
// ===================================================================================
COULAUD Olivier's avatar
COULAUD Olivier committed
16 17
// author Berenger Bramas and Olivier Coulaud
//
18 19 20 21 22 23 24
#ifndef FFmaGenericLoader_HPP
#define FFmaGenericLoader_HPP

#include <iostream>
#include <iomanip>
#include <fstream>
#include <string>
COULAUD Olivier's avatar
COULAUD Olivier committed
25 26
//
#include "Utils/FGlobal.hpp"
27
#include "FAbstractLoader.hpp"
COULAUD Olivier's avatar
COULAUD Olivier committed
28
#include "Utils/FPoint.hpp"
29

COULAUD Olivier's avatar
COULAUD Olivier committed
30 31 32 33 34 35 36 37 38 39 40 41 42 43 44

//


//! \class  FmaBasicParticle
//!
//! \brief Basic  Particle class used in FMA loader and writer
//!
//!  Here we consider the position, the physical value in the structure
//! but we read only the four first values the position, the physical value
//! this class is used in the generateDistributions example.
//!
class FmaBasicParticle {
public:
	FPoint position;            ///< position of the particle
45
	//	FReal  x,y,z;                 ///< position of the particle
COULAUD Olivier's avatar
COULAUD Olivier committed
46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130
	FReal  physicalValue;    ///<  its physical value
	/**
	 *  return a pointer on the first value of the structure
	 */
	FReal  * getPtrFirstData()
	{return position.getDataValue() ;}
	const  FReal * getPtrFirstData() const
	{return position.getDataValue() ;}
	/**
	 *  return The number of data to read. it is used in FMAGenericLoader  (here 4)
	 */
	int getReadDataNumber()
	{ return 4;}
	/**
	 *  return The number of data to write. it is used in FMAGenericWriter (here 4)
	 */
	int getWriteDataNumber() const
	{ return 4;}
	/**
	 *  return size in byte of the structure. it is used in FMAGenericWriter
	 */
	unsigned int getWriteDataSize() const
	{ return sizeof(FmaBasicParticle);}
};

//! \class  FmaRParticle
//!
//! \brief The Particle class used in FMA loader and writer
//!
//! Here we consider the position, the physical value, the potential and the force in the structure
//! but we read only the four first values the position, the physical value and we write all the data in a file.
//! This class is used if you create an array of particles and you read them from a file generate by gerateDistributions
//! and you want to store the result.
//!
class FmaRParticle {
public:
	FPoint position;            ///< position of the particle
	FReal  physicalValue;    ///< its physical value (mass or charge
	FReal  potential;           ///< the potential
	FReal  forces[3];            ///<the force

	/**
	 *  return a pointer on the first value of the structure
	 */
	FReal * getPtrFirstData()
	{return position.getDataValue() ;}
	const  FReal * getPtrFirstData() const
	{return position.getDataValue() ;}
	/**
	 *  return The number of data to read. it is used in FMAGenericLoader  (here 4)
	 */
	int getReadDataNumber()
	{ return 4;}
	/**
	 *  return The number of data to write. it is used in FMAGenericWriter (here 8)
	 */
	int getWriteDataNumber() const
	{ return 8;}
	/**
	 *  return size in byte of the structure. it is used in FMAGenericWriter
	 */
	unsigned int getWriteDataSize() const
	{ return sizeof(FmaRParticle);}
};
//
//! \class  FmaParticle
//!
//! \brief The Particle class used in FMA loader and writer
//!
//! Here we consider the position, the physical value, the potential and the force in the structure
//! and we read and write all the data in the structure
//! This class is used if you create an array of particles  from a file generate by a previous computation (directComputation or a FMM computation)
//! and you want to compare the result.

//!
class FmaParticle : public FmaRParticle {
public:
	/**
	 *  return The number of data to read. it is used in FMAGenericLoader  (here 8)
	 *  Useful to read the result of the output of DirectComputation and to compare the potential and the force
	 *  with the FMM computation
	 */
	int getReadDataNumber()
	{ return 8;}
};
131 132 133 134

typedef FmaBasicParticle FmaR4W4Particle ;
typedef FmaRParticle       FmaR4W8Particle ;
typedef FmaParticle         FmaR8W8Particle ;
COULAUD Olivier's avatar
COULAUD Olivier committed
135 136 137 138 139 140 141 142
//
//! \class  FFmaGenericLoader
//!
//! \brief Read a set of particles in FMA format
//!
//! The FMA format is a simplest format to store the particles in a file.
//!
//!  It is organized as follow<br>
COULAUD Olivier's avatar
COULAUD Olivier committed
143
//!   DatatypeSise Number_of_record_per_line  <br>
COULAUD Olivier's avatar
COULAUD Olivier committed
144
//!    NB_particles half_Box_width Center_X Center_Y Center_Z <br>
COULAUD Olivier's avatar
COULAUD Olivier committed
145 146 147 148 149 150 151 152 153 154
//!    Particle Values
//!
//!  if  Number_of_record_per_line is  <br>
//!      4    the Particle Values  are X Y Z Q <br>
//!      8   the Particle Values  are X X Y Z Q  P FX FY FZ<br>
//!
//! There is 3 methods to read the data from the file <br>
//!  1) 	fillParticle(FPoint*, FReal*);<br>
//!  2)   fillParticle(FReal*, int);<br>
//!  3   fillParticle(PartClass);<br>
COULAUD Olivier's avatar
COULAUD Olivier committed
155 156
//!
//!  \code
COULAUD Olivier's avatar
COULAUD Olivier committed
157
//!     FFmaGenericLoader  loader("../Data/unitCubeXYZQ20k.fma");    // extension fma --> ascii format
COULAUD Olivier's avatar
COULAUD Olivier committed
158
//!
COULAUD Olivier's avatar
COULAUD Olivier committed
159 160
//!	 FmaRParticle * const particles = new FmaRParticle[nbParticles];
//!	  memset(particles, 0, sizeof(FmaRParticle) * nbParticles) ;
COULAUD Olivier's avatar
COULAUD Olivier committed
161 162 163 164
//!
//!    FPoint position ;
//!    Freal   physicalValue ;
//!    for(int idx = 0 ; idx < loader.getNumberOfParticles() ; ++idx){
COULAUD Olivier's avatar
COULAUD Olivier committed
165
//!              loader.fillParticle(particles[idx]);
COULAUD Olivier's avatar
COULAUD Olivier committed
166 167 168
//!    }
//! \endcode
//!
COULAUD Olivier's avatar
COULAUD Olivier committed
169
//! \warning It works only in shared memory (doesn't work with MPI)
COULAUD Olivier's avatar
COULAUD Olivier committed
170
//!
COULAUD Olivier's avatar
COULAUD Olivier committed
171 172
//
class FFmaGenericLoader : public FAbstractLoader  {
173
protected:
COULAUD Olivier's avatar
COULAUD Olivier committed
174 175 176 177
	std::fstream *file;                   ///< the stream used to read the file
	bool binaryFile  ;                    ///< if true the file to read is in binary mode
	FPoint     centerOfBox;            ///< The center of box (read from file)
	FReal boxWidth;                     ///< the box width (read from file)
COULAUD Olivier's avatar
COULAUD Olivier committed
178 179 180 181 182
	FSize    nbParticles;                ///< the number of particles (read from file)
	unsigned int typeData[2];      ///< {Size of the data to read, number of data on 1 line}
private:
	FReal   *tmpVal ;                         ///  Temporary array to read data
	int       otherDataToRead  ;       ///< <<number of other data (>4)to read in a particle record
183 184 185 186
public:
	/**
	 * The constructor need the file name
	 * @param filename the name of the file to open
COULAUD Olivier's avatar
COULAUD Olivier committed
187 188
	 * @param binary   true if the file to open is in binary mode
	 *
COULAUD Olivier's avatar
COULAUD Olivier committed
189
	 *  This function also read the header of the fma file (Two first lines)
COULAUD Olivier's avatar
COULAUD Olivier committed
190
	 *  This means that we can obtain after the call the number of particles, ...
191 192
	 * you can test if file is successfully open by calling hasNotFinished()
	 */
COULAUD Olivier's avatar
COULAUD Olivier committed
193 194
	FFmaGenericLoader(const std::string & filename,const bool binary ) :file(nullptr),binaryFile(binary),
	centerOfBox(0.0,0.0,0.0),boxWidth(0.0),nbParticles(0),tmpVal(nullptr),otherDataToRead(0){
195
		if(binary) {
196
			this->file = new std::fstream (filename.c_str(),std::ifstream::in| std::ios::binary);
197 198
		}
		else {
199
			this->file = new std::fstream(filename.c_str(),std::ifstream::in) ;
200 201
		}
		// test if open
COULAUD Olivier's avatar
COULAUD Olivier committed
202 203 204
		if(! this->file->is_open()){
			std::cerr << "File "<< filename<<" not opened! " <<std::endl;
			exit(-1);
205
		}
COULAUD Olivier's avatar
COULAUD Olivier committed
206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233
		this->readHeader();
	}
	/**
	 * The constructor need the file name
	 * @param filename the name with the extension .fma, .bfma of the file to open. We check the extension to know if the file is in ASCII mode or binary mode
	 *
	 *  This function also read the header of the  file (Two first lines)
	 *  This means that we can obtain after the call the number of particles, ...
	 * you can test if file is successfully open by calling hasNotFinished()
	 */
	FFmaGenericLoader(const std::string & filename) : binaryFile(false) ,tmpVal(nullptr),otherDataToRead(0){
		std::string ext(".bfma");
		// open particle file
		if(filename.find(ext) != std::string::npos) {
			binaryFile = true;
			this->file = new std::fstream (filename.c_str(),std::ifstream::in| std::ios::binary);
		}
		else if(filename.find(".fma")!=std::string::npos ) {
			this->file = new std::fstream(filename.c_str(),std::ifstream::in) ;
		}
		else  {
			std::cout << "Input file not allowed only .fma or .bfma extensions" <<std::endl;
			exit (-1) ;
		}
		// test if open
		if(! this->file->is_open()){
			std::cerr << "File "<< filename<<" not opened! " <<std::endl;
			exit(-1);
234
		}
COULAUD Olivier's avatar
COULAUD Olivier committed
235
		this->readHeader();
236 237 238 239 240 241
	}
	/**
	 * Default destructor, simply close the file
	 */
	virtual ~FFmaGenericLoader(){
		file->close();
COULAUD Olivier's avatar
COULAUD Olivier committed
242
		delete file ;
COULAUD Olivier's avatar
COULAUD Olivier committed
243
		delete tmpVal;
244 245 246 247 248 249 250 251 252 253 254 255 256 257 258
	}

	/**
	 * To know if file is open and ready to read
	 * @return true if loader can work
	 */
	bool isOpen() const{
		return this->file->is_open() && !this->file->eof();
	}

	/**
	 * To get the number of particles from this loader
	 * @param the number of particles the loader can fill
	 */
	FSize getNumberOfParticles() const{
COULAUD Olivier's avatar
COULAUD Olivier committed
259
		return this->nbParticles;
260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275
	}

	/**
	 * The center of the box from the simulation file opened by the loader
	 * @return box center
	 */
	FPoint getCenterOfBox() const{
		return this->centerOfBox;
	}
	/**
	 * The box width from the simulation file opened by the loader
	 * @return box width
	 */
	FReal getBoxWidth() const{
		return this->boxWidth;
	}
COULAUD Olivier's avatar
COULAUD Olivier committed
276 277 278 279 280 281 282 283 284 285 286 287
	/**
	 * The box width from the simulation file opened by the loader
	 * @return the number of data per record (Particle)
	 */
	unsigned int getNbRecordPerline(){
		return typeData[1]; }
	/**
	 * To know if the data are in float or in double type
	 * @return the type of the values float (4) or double (8)
	 */
	unsigned int getDataType(){
		return typeData[0]; }
288 289

	/**
COULAUD Olivier's avatar
COULAUD Olivier committed
290 291 292 293 294
	 * Fill a particle form the current position in the file
	 *
	 * @param outParticlePositions the position of particle to fill (FPoint class)
	 * @param outPhysicalValue the physical value of particle to fill (FReal)
	 *
295 296 297 298 299
	 */
	void fillParticle(FPoint*const outParticlePositions, FReal*const outPhysicalValue){
		if(binaryFile){
			file->read((char*)(outParticlePositions), sizeof(FReal)*3);
			file->read((char*)(outPhysicalValue), sizeof(FReal));
COULAUD Olivier's avatar
COULAUD Olivier committed
300 301 302
			if(otherDataToRead> 0){
				file->read((char*)(this->tmpVal), sizeof(FReal)*otherDataToRead);
			}
303 304
		}
		else{
COULAUD Olivier's avatar
COULAUD Olivier committed
305
			FReal x,y,z;
COULAUD Olivier's avatar
COULAUD Olivier committed
306
			//			(*this->file)  >> &outParticlePositions>> &outPhysicalValue;
COULAUD Olivier's avatar
COULAUD Olivier committed
307
			(*this->file)  >> x >> y >> z >> (*outPhysicalValue);
308
			outParticlePositions->setPosition(x,y,z);
COULAUD Olivier's avatar
COULAUD Olivier committed
309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329

			if(otherDataToRead> 0){
				for (int 	i = 0 ; i <otherDataToRead; ++i){
					(*this->file) >> x ;
				}
			}
		}
		//		std::cout <<  "X Y Z Q " << *outParticlePositions<<" "<<*outPhysicalValue <<std::endl;

	}
	/**
	 * Fill a particle form the current position in the file
	 * @param dataToRead is an array of type FReal. It contains all the values of a particles (for instance X,Y,Z,Q, ..
	 * @param nbDataToRead number of value to read (I.e. size of the array)
	 */
	void fillParticle(FReal* dataToRead, const int nbDataToRead){
		if(binaryFile){
			file->read((char*)(dataToRead), sizeof(FReal)*nbDataToRead);
			if(nbDataToRead< typeData[1]){
				file->read((char*)(this->tmpVal), sizeof(FReal)*(typeData[1]-nbDataToRead));
			}
330
		}
COULAUD Olivier's avatar
COULAUD Olivier committed
331
		else{
332

COULAUD Olivier's avatar
COULAUD Olivier committed
333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351
			for (int i = 0 ; i <nbDataToRead; ++i){
				(*this->file)  >>dataToRead[i];
			}
			if(nbDataToRead< typeData[1]){
				FReal x;
				for (int 	i = 0 ; i <typeData[1]-nbDataToRead; ++i){
					(*this->file) >> x ;
				}
			}
		}
	}

	/**
	 * Fill a particle form the current position in the file
	 * @param dataToRead is the particle  class. If the class is different from FmaBasicParticle, FmaRParticle or FmaParticle the we have to implement the following method
	 *   getPtrFirstData(), getReadDataNumber() see FmaBasicParticle class for more details.
	 */
	template <class dataPart>
	void fillParticle(dataPart &dataToRead){
352 353
		int otherDataRead = typeData[1] - dataToRead.getReadDataNumber() ;
		if(binaryFile){
COULAUD Olivier's avatar
COULAUD Olivier committed
354
			file->read((char*)(dataToRead.getPtrFirstData()), sizeof(FReal)*(dataToRead.getReadDataNumber()));
355 356
			if( otherDataRead > 0){
				file->read((char*)(this->tmpVal), sizeof(FReal)*(otherDataRead));
COULAUD Olivier's avatar
COULAUD Olivier committed
357 358 359 360 361 362 363 364
			}
		}
		else{
			FReal * val = dataToRead.getPtrFirstData();
			for (int i = 0 ; i <dataToRead.getReadDataNumber(); ++i){
				(*this->file)  >>*val;
				++val;
			}
365
			if( otherDataRead > 0){
COULAUD Olivier's avatar
COULAUD Olivier committed
366
				FReal x;
367
				for (int i = 0 ; i <otherDataRead ;++i){
COULAUD Olivier's avatar
COULAUD Olivier committed
368 369 370 371 372 373 374 375 376 377 378 379 380
					(*this->file)  >>x;
				}
			}
		}
	}
	/**
	 * Fill a set of particles form the current position in the file.
	 *  If the file is a binary file and we read all record per particle then we read and fill the array in one instruction
	 * @param dataToRead is an array of the particle  class. If the class is different from FmaBasicParticle, FmaRParticle or FmaParticle the we have to implement the following method
	 *   getPtrFirstData(), getReadDataNumber() see FmaBasicParticle class for more details.
	 */
	template <class dataPart>
	void fillParticle(dataPart *dataToRead, const int N){
381 382
		int otherDataRead = typeData[1] - (*dataToRead).getReadDataNumber() ;
		if(binaryFile && otherDataRead == 0 ){
COULAUD Olivier's avatar
COULAUD Olivier committed
383 384 385 386 387 388 389
			file->read((char*)((*dataToRead).getPtrFirstData()), sizeof(FReal)*(N*(*dataToRead).getReadDataNumber()));
		}
		else {
			for (int i = 0 ; i <N; ++i){
				this->fillParticle(dataToRead[i]) ;
			}
		}
390
	}
COULAUD Olivier's avatar
COULAUD Olivier committed
391
private:
COULAUD Olivier's avatar
COULAUD Olivier committed
392 393 394 395 396 397 398 399 400 401 402 403
	void readHeader() {
		if(this->binaryFile){
			this->readBinaryHeader();
		}
		else {
			this->readAscciHeader();
		}

		std::cout << "   nbParticles: " <<this->nbParticles << std::endl
				<< "   Box width:   " <<this->boxWidth << std::endl
				<< "   Center:        " << this->centerOfBox << std::endl;
	}
COULAUD Olivier's avatar
COULAUD Olivier committed
404
	void readAscciHeader() {
COULAUD Olivier's avatar
COULAUD Olivier committed
405
		std::cout << " File open in ASCII mode "<< std::endl ;
COULAUD Olivier's avatar
COULAUD Olivier committed
406
		FReal x,y,z;
COULAUD Olivier's avatar
COULAUD Olivier committed
407 408
		(*this->file) >> typeData[0]>> typeData[1];
		std::cout << "   Datatype "<< typeData[0] << " "<< typeData[1] << std::endl;
COULAUD Olivier's avatar
COULAUD Olivier committed
409 410 411
		(*this->file) >> this->nbParticles >> this->boxWidth >> x >> y >> z;
		this->centerOfBox.setPosition(x,y,z);
		this->boxWidth *= 2;
COULAUD Olivier's avatar
COULAUD Olivier committed
412 413 414 415 416 417 418 419 420 421 422 423 424 425 426 427 428 429 430 431 432 433 434 435 436 437 438 439 440 441 442 443 444 445 446 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 469 470 471 472 473 474 475 476 477 478 479 480 481 482 483 484 485 486 487 488 489 490 491 492 493 494 495 496 497 498 499 500 501 502 503 504 505 506 507 508 509 510 511 512 513 514 515 516 517 518 519 520 521 522 523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552 553 554 555 556 557 558 559 560 561 562 563 564 565 566
		otherDataToRead = typeData[1] -4;
	};
	void readBinaryHeader(){
		std::cout << " File open in binary mode "<< std::endl;
		file->seekg (std::ios::beg);
		file->read((char*)&typeData,2*sizeof(unsigned int));
		std::cout << "   Datatype "<< typeData[0] << " "<< typeData[1] << std::endl;
		if(typeData[0] != sizeof(FReal)){
			std::cerr << "Size of elements in part file " << typeData[0] << " is different from size of FReal " << sizeof(FReal)<<std::endl;
			exit(-1);
		}
		else{
			file->read( (char*)&(this->nbParticles), sizeof(FSize) );
			file->read( (char*)&(this->boxWidth) ,sizeof(this->boxWidth) );
			this->boxWidth *= 2;

			FReal x[3];
			file->read( (char*)x,sizeof(FReal)*3);
			this->centerOfBox.setPosition(x[0],x[1],x[2]);
		}
		otherDataToRead = typeData[1] -4;
		if(otherDataToRead>0){
			tmpVal = new FReal[otherDataToRead];
		}
	}

};
//
//! \class  FFmaGenericWriter
//!
//! \brief Write a set of particles in FMA format (ASCII or binary mode)
//!
//! The FMA format is a simplest format to store the particles in a file.
//!
//!  It is organized as follow<br>
//!    DatatypeSise Number_of_record_per_line
//!    NB_particles half_Box_width Center_X Center_Y Center_Z <br>
//!    X Y Z PhysicalValue // one particle per line in ascii mode
//!
//!  DatatypeSise = 4 float ; 8 double
//!
//!  Number_of_record_per_line 4 (X,Y,Z,Q)  ; 8 (X,Y,Z,Q,P,FX,FY,FZ  )
//!  \code
//!     FFmaGenericWriter writer ("data.bfma");    // open binary fma file (extension .bfma)
//!
//!      // write the header of the file
//!     writer.writeHeader( loader.getCenterOfBox(), loader.getBoxWidth() , NbPoints, sizeof(FReal), nbData) ;
//!
//!      // write the data here particles is an array of float or double and a particle has nbData value
//!     writer.writeArrayOfReal(particles, nbData, NbPoints);
//!
//! \endcode
//! \warning It works only in shared memory (doesn't work with MPI)
//!

class FFmaGenericWriter {

protected:
	std::fstream *file;                   ///< the stream used to read the file
	bool binaryFile  ;                    ///< if true the file to read is in binary mode

public:
	/**
	 * The constructor need the file name
	 * @param filename the name of the file to open
	 * @param binary   true if the file to open is in binary mode
	 *
	 *  This function also read the header of the fma file (First line)
	 *  This means that we can obtain after the call the number of particles, ...
	 * you can test if file is successfully open by calling hasNotFinished()
	 */
	FFmaGenericWriter(const std::string & filename): binaryFile(false){
		std::string ext(".bfma");
		// open particle file
		if(filename.find(ext) !=std::string::npos) {
			binaryFile = true;
			this->file = new std::fstream (filename.c_str(),std::ifstream::out| std::ios::binary);
		}
		else if(filename.find(".fma")!=std::string::npos ) {
			this->file = new std::fstream(filename.c_str(),std::ifstream::out) ;
			this->file->precision(10);
		}
		else  {
			std::cout << "Input file not allowed only .fma or .bfma extensions" <<std::endl;
			exit (-1) ;
		}
		// test if open
		if(! this->file->is_open()){
			std::cerr << "File "<< filename<<" not opened! " <<std::endl;
			exit(-1);
		}
	}
	/**
	 * The constructor need the file name
	 * @param filename the name of the file to open
	 * @param binary   true if the file to open is in binary mode
	 *
	 *  This function also read the header of the fma file (First line)
	 *  This means that we can obtain after the call the number of particles, ...
	 * you can test if file is successfully open by calling hasNotFinished()
	 */
	FFmaGenericWriter(const std::string & filename, const bool binary ) : file(nullptr), binaryFile(binary)
	{
		if(binary) {
			this->file = new std::fstream (filename.c_str(),std::ifstream::out| std::ios::binary);
		}
		else {
			this->file = new std::fstream(filename.c_str(),std::ifstream::out) ;
			this->file->precision(10);
		}
		// test if open
		if(! this->file->is_open()){
			std::cerr << "File "<< filename<<" not opened! " <<std::endl;
			exit(-1);
		}
	}
	/**
	 * Default destructor, simply close the file
	 */
	virtual ~FFmaGenericWriter(){
		file->close();
		delete file ;
	}

	/**
	 * To know if file is open and ready to read
	 * @return true if loader can work
	 */
	bool isOpen() const{
		return this->file->is_open() && !this->file->eof();
	}
	//!
	//! Write the header of FMA file
	//! \warning All values inside typePart should be of the same type {float or double}
	//!
	//!@param centerOfBox  The centre of the Box (FPoint class)
	//!@param boxWidth       The width of the box
	//!@param nbParticles     Number of particles in the box (or to save)
	//!@param data               Data type of the particle class (FmaBasicParticle, FmaRParticle or FmaParticle)
	//!
	template <class typePart>
	void writeHeader(const FPoint &centerOfBox,const FReal &boxWidth, const FSize &nbParticles, const typePart  data) {
		unsigned int typeFReal[2]  = {sizeof(FReal) , sizeof(typePart) / sizeof(FReal) };
		const unsigned int ndata  = data.getWriteDataNumber();
		std::cout <<"    WriteHeader: typeFReal: " << typeFReal[0]  << "  nb Elts: " << typeFReal[1]  <<"   NData to write "<< ndata<< "\n";
		if (ndata != typeFReal[1]){
			typeFReal[1] = ndata;
		}
		FReal x = boxWidth *0.5;
		if(this->binaryFile) {
			this->writerBinaryHeader(centerOfBox,x,nbParticles,typeFReal);
		}
		else {
			this->writerAscciHeader(centerOfBox,x,nbParticles,typeFReal);
		}
COULAUD Olivier's avatar
COULAUD Olivier committed
567
	}
COULAUD Olivier's avatar
COULAUD Olivier committed
568 569 570 571 572 573 574 575 576 577 578 579 580 581 582 583 584 585 586 587 588 589 590 591 592 593 594 595 596 597 598 599 600 601 602 603 604 605 606 607 608 609 610 611 612 613 614 615 616 617 618 619 620 621 622 623 624 625 626 627 628 629 630 631 632 633 634 635 636 637 638 639 640 641 642 643 644
	//!
	//! Write the header of FMA file. Should be used if we write the particles with writeArrayOfReal method
	//!
	//!
	//!@param centerOfBox      The center of the Box (FPoint class)
	//!@param boxWidth            The width of the box
	//!@param nbParticles          Number of particles in the box (or to save)
	//!@param dataType             Data type of the values in particle
	//!@param nbDataPerRecord  Number of record/value per particle
	//!
	void writeHeader(const FPoint &centerOfBox,const FReal &boxWidth, const FSize &nbParticles,
			const unsigned int  dataType, const unsigned int  nbDataPerRecord) {
		unsigned int typeFReal[2]  = {dataType , nbDataPerRecord };
		FReal x = boxWidth *0.5;
		if(this->binaryFile) {
			this->writerBinaryHeader(centerOfBox,x,nbParticles,typeFReal);
		}
		else {
			this->writerAscciHeader(centerOfBox,x,nbParticles,typeFReal);
		}
	}

	//!
	//!@warning the type dataPart should be FmaBasicParticle, FmaRParticle or FmaParticle  class or should implement all methods inside FmaBasicParticle.
	//!
	//!@param dataToWrite array of particles of type dataPart
	//!@param N number of element in the array
	//!
	//!  example 1
	//!\code
	//!	FmaRParticle *  particles = new FmaRParticle[nbParticles];
	//!     memset(particles, 0, sizeof(FmaRParticle) * nbParticles) ;
	//!    ...
	//! 	FFmaGenericWriter writer(filenameOut) ;
	//! 	Fwriter.writeHeader(Centre,BoxWith, nbParticles,*particles) ;
	//! 	Fwriter.writeArrayOfParticles(particles, nbParticles);
	//! \endcode
	//!
	//!  example2
	//!\code
	//!	FReal *  particles = new FReal[4*NbPoints] ; // store 4 data per particle
	//!     memset(particles, 0, sizeof(FmaRParticle) * nbParticles) ;
	//!    ...
	//!    FmaBasicParticle *ppart = (FmaBasicParticle*)(&particles[0]);
	//! 	FFmaGenericWriter writer(filenameOut) ;
	//! 	Fwriter.writeHeader(Centre,BoxWith, nbParticles,*particles) ;
	//! 	Fwriter.writeArrayOfParticles(particles, nbParticles);
	//! \endcode

	template <class dataPart>
	void writeArrayOfParticles(const dataPart *dataToWrite, const FSize N){
		//		std::cout << "NB points to write: "<< N <<std::endl;
		if(binaryFile){
			unsigned int recordSize=  dataToWrite[0].getWriteDataSize() ;
			unsigned int typeFReal[2]      = {sizeof(FReal) , sizeof(dataPart) / sizeof(FReal) };
			//			std::cout << "typeData "<< typeFReal[0] << " "<< typeFReal[1] <<"  "<< std::endl;
			//
			if (sizeof(dataPart) == recordSize){
				//				std::cout << "Size to write:  "<<N*dataToWrite[0].getWriteDataSize() <<std::endl;
				file->write((char*)(dataToWrite[0].getPtrFirstData()), N*recordSize);
			}
			else {
				file->write((char* )&typeFReal[0],2*sizeof(unsigned int));
				//				std::cout << "Size to write:   N* "<<typeFReal[0] *typeFReal[1]  <<std::endl;
				for (int i = 0 ; i <N ; ++i){
					file->write((char*)(dataToWrite[i].getPtrFirstData()), recordSize );
					//					const FReal * val = dataToWrite[i].getPtrFirstData() ;
					//					std::cout << i <<"   ";
					//					for( int j=0; j<typeFReal[1] ; ++j){
					//						std::cout << *val << "   ";++val;
					//					}
					//					std::cout <<std::endl;
				}
			}
		}
		else{ // ASCII part
			const int ndata = dataToWrite[0].getWriteDataNumber();
645 646
//			std::cout << "typeData "<< sizeof(FReal) << " "<<ndata << std::endl;
			this->file->precision(10);
COULAUD Olivier's avatar
COULAUD Olivier committed
647 648 649 650 651 652 653 654 655 656 657 658 659 660 661 662 663 664 665 666 667 668 669 670 671 672 673 674 675 676 677 678 679 680

			for (int i = 0 ; i <N ; ++i){
				const FReal * val = dataToWrite[i].getPtrFirstData() ;
				for (int j= 0 ; j <ndata ; ++j){
					(*this->file)  << *val << "    "; ++val;
				}
				(*this->file)  <<std::endl;
			}
		}
	}
	//!
	//! Write an array of data in a file Fill
	//!
	//!@param dataToWrite array of particles of type FReal
	//!@param nbData number of data per particle
	//!@param N number of particles
	//!
	//!  The size of the array is N*nbData
	//!
	/*	//!  example
	//!\code
	//!	FmaRParticle * const particles = new FmaRParticle[nbParticles];
	//!     memset(particles, 0, sizeof(FmaRParticle) * nbParticles) ;
	//!    ...
	//! 	FFmaGenericWriter writer(filenameOut) ;
	//! 	Fwriter.writeHeader(Centre,BoxWith, nbParticles,*particles) ;
	//! 	Fwriter.writeArrayOfReal(particles, nbParticles);
	//! \endcode
	 */
	void writeArrayOfReal(const FReal *dataToWrite, const unsigned int nbData, const FSize N){
		if(binaryFile){
			file->write((char*)(dataToWrite), N*nbData*sizeof(FReal));
		}
		else{
681
			this->file->precision(10);
COULAUD Olivier's avatar
COULAUD Olivier committed
682 683 684 685 686 687 688 689 690 691 692 693 694 695 696 697 698 699
			//			std::cout << "N "<< N << " nbData "<< nbData<<std::endl;
			//			exit(-1);
			int k = 0;
			for (int i = 0 ; i <N ; ++i){
				//				std::cout << "i "<< i << "  ";
				for (int jj= 0 ; jj<nbData ; ++jj, ++k){
					(*this->file)  << dataToWrite[k] << "    ";
					//					std::cout      << dataToWrite[k]<< "  ";
				}
				(*this->file)  <<std::endl;
				//				std::cout <<std::endl;
			}
			//			std::cout << "END"<<std::endl;
		}
	}
private:
	void writerAscciHeader( const FPoint &centerOfBox,const FReal &boxWidth,
			const FSize &nbParticles, const unsigned int *typeFReal) {
700
		this->file->precision(10);
COULAUD Olivier's avatar
COULAUD Olivier committed
701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718 719 720
		(*this->file) << typeFReal[0] <<"   "<<typeFReal[1]<<std::endl;
		(*this->file) << nbParticles << "   "<<  boxWidth << "   "
				<<  centerOfBox.getX()  << "  " << centerOfBox.getY() << " "<<centerOfBox.getZ()
				<< std::endl;
	}
	void writerBinaryHeader(const FPoint &centerOfBox,const FReal &boxWidth,
			const FSize &nbParticles, const unsigned int *typeFReal) {
		int sizeOfElement;
		file->seekg (std::ios::beg);
		file->write((char*)typeFReal,2*sizeof(unsigned int));
		if(typeFReal[0]  != sizeof(FReal)){
			std::cout << "Size of elements in part file " << typeFReal[0] << " is different from size of FReal " << sizeof(FReal)<<std::endl;
			exit(-1);
		}
		else{
			file->write( (char*)&(nbParticles), sizeof(FSize) );
			//			std::cout << "nbParticles "<< nbParticles<<std::endl;
			file->write( (char*)&(boxWidth) ,sizeof(boxWidth) );
			file->write( (char*)(centerOfBox.getDataValue()),sizeof(FReal)*3);
		}
COULAUD Olivier's avatar
COULAUD Olivier committed
721 722
	}

723 724 725 726 727 728
};


#endif //FFmaGenericLoader_HPP