1 /*
2 MDAL - Mesh Data Abstraction Library (MIT License)
3 Copyright (C) 2019 Peter Petrik (zilolv at gmail dot com)
4 */
5
6 #include <stdio.h>
7 #include <string.h>
8 #include <stddef.h>
9 #include <iosfwd>
10 #include <iostream>
11 #include <fstream>
12 #include <sstream>
13 #include <string>
14 #include <vector>
15 #include <map>
16 #include <cassert>
17 #include <memory>
18 #include <algorithm>
19
20 #include "mdal_selafin.hpp"
21 #include "mdal.h"
22 #include "mdal_utils.hpp"
23 #include <math.h>
24 #include "mdal_logger.hpp"
25
26 #define BUFFER_SIZE 2000
27
28 // //////////////////////////////
29 // SelafinFile
30 // //////////////////////////////
31
SelafinFile(const std::string & fileName)32 MDAL::SelafinFile::SelafinFile( const std::string &fileName ):
33 mFileName( fileName )
34 {}
35
initialize()36 void MDAL::SelafinFile::initialize()
37 {
38 if ( !MDAL::fileExists( mFileName ) )
39 {
40 throw MDAL::Error( MDAL_Status::Err_FileNotFound, "Did not find file " + mFileName );
41 }
42 mIn = MDAL::openInputFile( mFileName, std::ios_base::in | std::ios_base::binary );
43 if ( !mIn )
44 throw MDAL::Error( MDAL_Status::Err_FileNotFound, "File " + mFileName + " could not be open" ); // Couldn't open the file
45
46 // get length of file:
47 mIn.seekg( 0, mIn.end );
48 mFileSize = mIn.tellg();
49 mIn.seekg( 0, mIn.beg );
50
51 mChangeEndianness = MDAL::isNativeLittleEndian();
52
53 //Check if need to change the endianness
54 // read first size_t that has to be 80
55 size_t firstInt = readSizeT();
56 mIn.seekg( 0, mIn.beg );
57 if ( firstInt != 80 )
58 {
59 mChangeEndianness = !mChangeEndianness;
60 //Retry
61 firstInt = readSizeT();
62 if ( firstInt != 80 )
63 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "File " + mFileName + " could not be open" );
64 mIn.seekg( 0, mIn.beg );
65 }
66
67 mParsed = false;
68 }
69
parseMeshFrame()70 void MDAL::SelafinFile::parseMeshFrame()
71 {
72 /* 1 record containing the title of the study (72 characters) and a 8 characters
73 string indicating the type of format (SERAFIN or SERAFIND)
74 */
75 readHeader();
76
77 /* 1 record containing the two integers NBV(1) and NBV(2) (number of linear
78 and quadratic variables, NBV(2) with the value of 0 for Telemac, as
79 quadratic values are not saved so far), numbered from 1 in docs
80 */
81 std::vector<int> nbv = readIntArr( 2 );
82
83 /* NBV(1) records containing the names and units of each variable (over 32 characters)
84 */
85 mVariableNames.clear();
86 for ( int i = 0; i < nbv[0]; ++i )
87 {
88 mVariableNames.push_back( readString( 32 ) );
89 }
90
91 /* 1 record containing the integers table IPARAM (10 integers, of which only
92 the 6 are currently being used), numbered from 1 in docs
93
94 - if IPARAM (3) != 0: the value corresponds to the x-coordinate of the
95 origin of the mesh,
96
97 - if IPARAM (4) != 0: the value corresponds to the y-coordinate of the
98 origin of the mesh,
99
100 - if IPARAM (7) != 0: the value corresponds to the number of planes
101 on the vertical (3D computation),
102
103 - if IPARAM (8) != 0: the value corresponds to the number of
104 boundary points (in parallel),
105
106 - if IPARAM (9) != 0: the value corresponds to the number of interface
107 points (in parallel),
108
109 - if IPARAM(8) or IPARAM(9) != 0: the array IPOBO below is replaced
110 by the array KNOLG (total initial number of points). All the other
111 numbers are local to the sub-domain, including IKLE
112 */
113 mParameters = readIntArr( 10 );
114 mXOrigin = static_cast<double>( mParameters[2] );
115 mYOrigin = static_cast<double>( mParameters[3] );
116
117 if ( mParameters[6] != 0 && mParameters[6] != 1 ) //some tools set this value to one for 2D mesh
118 {
119 // would need additional parsing
120 throw MDAL::Error( MDAL_Status::Err_MissingDriver, "File " + mFileName + " would need additional parsing" );
121 }
122
123 /*
124 if IPARAM (10) = 1: a record containing the computation starting date
125 */
126
127 if ( mParameters[9] == 1 )
128 {
129 std::vector<int> datetime = readIntArr( 6 );
130 mReferenceTime = DateTime( datetime[0], datetime[1], datetime[2], datetime[3], datetime[4], double( datetime[5] ) );
131 }
132
133 /* 1 record containing the integers NELEM,NPOIN,NDP,1 (number of
134 elements, number of points, number of points per element and the value 1)
135 */
136 std::vector<int> numbers = readIntArr( 4 );
137 mFacesCount = numbers[0];
138 mVerticesCount = numbers[1];
139 mVerticesPerFace = numbers[2];
140
141 /* 1 record containing table IKLE (integer array
142 of dimension (NDP,NELEM)
143 which is the connectivity table.
144 */
145 size_t size = mFacesCount * mVerticesPerFace;
146 if ( ! checkIntArraySize( size ) )
147 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "File format problem while reading connectivity table" );
148 mConnectivityStreamPosition = passThroughIntArray( size );
149
150 /* 1 record containing table IPOBO (integer array of dimension NPOIN); the
151 value of one element is 0 for an internal point, and
152 gives the numbering of boundary points for the others
153 */
154 size = mVerticesCount;
155 if ( ! checkIntArraySize( size ) )
156 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "File format problem while reading IPOBO table" );
157 mIPOBOStreamPosition = passThroughIntArray( size );
158
159 /* 1 record containing table X (real array of dimension NPOIN containing the
160 abscisse of the points)
161 AND here, we can know if float of this file is simple or double precision:
162 result of size of record divided by number of vertices gives the byte size of the float:
163 -> 4 : simple precision -> 8 : double precision
164 */
165 size = mVerticesCount;
166 size_t recordSize = readSizeT();
167 mStreamInFloatPrecision = recordSize / size == 4;
168 if ( !mStreamInFloatPrecision && recordSize / size != 8 )
169 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "File format problem: could not determine if simple or double precision" );
170 mXStreamPosition = passThroughDoubleArray( size );
171
172 /* 1 record containing table Y (real array of dimension NPOIN containing the
173 ordinate of the points)
174 */
175 size = mVerticesCount;
176 if ( ! checkDoubleArraySize( size ) )
177 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "File format problem while reading abscisse values" );
178 mYStreamPosition = passThroughDoubleArray( size );
179 }
180
parseFile()181 void MDAL::SelafinFile::parseFile()
182 {
183 parseMeshFrame();
184
185 /* Next, for each time step, the following are found:
186 - 1 record containing time T (real),
187 - NBV(1)+NBV(2) records containing the results tables for each variable at time
188 */
189
190 size_t realSize = mStreamInFloatPrecision ? 4 : 8;
191 size_t nTimesteps = remainingBytes() / ( 8 + realSize + ( 4 + ( mVerticesCount ) * realSize + 4 ) * mVariableNames.size() );
192 mVariableStreamPosition.resize( mVariableNames.size(), std::vector<std::streampos>( nTimesteps ) );
193 mTimeSteps.resize( nTimesteps );
194 for ( size_t nT = 0; nT < nTimesteps; ++nT )
195 {
196 std::vector<double> times = readDoubleArr( 1 );
197 mTimeSteps[nT] = RelativeTimestamp( times[0], RelativeTimestamp::seconds );
198 for ( size_t i = 0; i < mVariableNames.size(); ++i )
199 {
200 if ( ! checkDoubleArraySize( mVerticesCount ) )
201 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "File format problem while reading dataset values" );
202 mVariableStreamPosition[i][nT] = passThroughDoubleArray( mVerticesCount );
203 }
204 }
205
206 mParsed = true;
207 }
208
readHeader()209 std::string MDAL::SelafinFile::readHeader()
210 {
211 initialize();
212 std::string header = readString( 80 );
213
214 std::string title = header.substr( 0, 72 );
215 title = trim( title );
216
217 if ( header.size() < 80 ) // IF "SERAFIN", the readString method remove the last character that is a space
218 header.append( " " );
219 return header;
220 }
221
connectivityIndex(size_t offset,size_t count)222 std::vector<int> MDAL::SelafinFile::connectivityIndex( size_t offset, size_t count )
223 {
224 return readIntArr( mConnectivityStreamPosition, offset, count );
225 }
226
vertices(size_t offset,size_t count)227 std::vector<double> MDAL::SelafinFile::vertices( size_t offset, size_t count )
228 {
229 std::vector<double> xValues = readDoubleArr( mXStreamPosition, offset, count );
230 std::vector<double> yValues = readDoubleArr( mYStreamPosition, offset, count );
231
232 if ( xValues.size() != count || yValues.size() != count )
233 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "File format problem while reading vertices" );
234
235 std::vector<double> coordinates( count * 3 );
236 for ( size_t i = 0; i < count; ++i )
237 {
238 coordinates[i * 3] = xValues.at( i ) + mXOrigin;
239 coordinates[i * 3 + 1] = yValues.at( i ) + mYOrigin;
240 coordinates[i * 3 + 2] = 0;
241 }
242 return coordinates;
243 }
244
createMesh(const std::string & fileName)245 std::unique_ptr<MDAL::Mesh> MDAL::SelafinFile::createMesh( const std::string &fileName )
246 {
247 std::shared_ptr<SelafinFile> reader = std::make_shared<SelafinFile>( fileName );
248 reader->initialize();
249 reader->parseFile();
250
251 std::unique_ptr<Mesh> mesh( new MeshSelafin( fileName, reader ) );
252 populateDataset( mesh.get(), reader );
253
254 return mesh;
255 }
256
populateDataset(MDAL::Mesh * mesh,const std::string & fileName)257 void MDAL::SelafinFile::populateDataset( MDAL::Mesh *mesh, const std::string &fileName )
258 {
259 std::shared_ptr<SelafinFile> reader = std::make_shared<SelafinFile>( fileName );
260 reader->initialize();
261 reader->parseFile();
262
263 if ( mesh->verticesCount() != reader->verticesCount() || mesh->facesCount() != reader->facesCount() )
264 throw MDAL::Error( MDAL_Status::Err_IncompatibleDataset, "Faces or vertices counts in the file are not the same" );
265
266 populateDataset( mesh, reader );
267 }
268
facesCount()269 size_t MDAL::SelafinFile::facesCount()
270 {
271 if ( !mParsed )
272 parseFile();
273 return mFacesCount;
274 }
275
verticesCount()276 size_t MDAL::SelafinFile::verticesCount()
277 {
278 if ( !mParsed )
279 parseFile();
280 return mVerticesCount;
281 }
282
verticesPerFace()283 size_t MDAL::SelafinFile::verticesPerFace()
284 {
285 if ( !mParsed )
286 parseFile();
287 return mVerticesPerFace;
288 }
289
datasetValues(size_t timeStepIndex,size_t variableIndex,size_t offset,size_t count)290 std::vector<double> MDAL::SelafinFile::datasetValues( size_t timeStepIndex, size_t variableIndex, size_t offset, size_t count )
291 {
292 if ( !mParsed )
293 parseFile();
294 if ( variableIndex < mVariableStreamPosition.size() && timeStepIndex < mVariableStreamPosition[variableIndex].size() )
295 return readDoubleArr( mVariableStreamPosition[variableIndex][timeStepIndex], offset, count );
296 else
297 return std::vector<double>();
298 }
299
populateDataset(MDAL::Mesh * mesh,std::shared_ptr<MDAL::SelafinFile> reader)300 void MDAL::SelafinFile::populateDataset( MDAL::Mesh *mesh, std::shared_ptr<MDAL::SelafinFile> reader )
301 {
302 std::map<std::string, std::shared_ptr<DatasetGroup>> groupsByName;
303 std::vector< std::shared_ptr<DatasetGroup>> groupsInOrder;
304
305 for ( size_t nName = 0; nName < reader->mVariableNames.size(); ++nName )
306 {
307 std::string var_name( reader->mVariableNames[nName] );
308 var_name = MDAL::toLower( MDAL::trim( var_name ) );
309 var_name = MDAL::replace( var_name, "/", "" ); // slash is represented as sub-dataset group but we do not have such subdatasets groups in Selafin
310
311 bool is_vector = false;
312 bool is_x = true;
313
314 if ( MDAL::contains( var_name, "velocity u" ) || MDAL::contains( var_name, "along x" ) )
315 {
316 is_vector = true;
317 var_name = MDAL::replace( var_name, "velocity u", "velocity" );
318 var_name = MDAL::replace( var_name, "along x", "" );
319 }
320 else if ( MDAL::contains( var_name, "velocity v" ) || MDAL::contains( var_name, "along y" ) )
321 {
322 is_vector = true;
323 is_x = false;
324 var_name = MDAL::replace( var_name, "velocity v", "velocity" );
325 var_name = MDAL::replace( var_name, "along y", "" );
326 }
327
328 if ( MDAL::contains( var_name, "vitesse u" ) || MDAL::contains( var_name, "suivant x" ) )
329 {
330 is_vector = true;
331 var_name = MDAL::replace( var_name, "vitesse u", "vitesse" );
332 var_name = MDAL::replace( var_name, "suivant x", "" );
333 }
334 else if ( MDAL::contains( var_name, "vitesse v" ) || MDAL::contains( var_name, "suivant y" ) )
335 {
336 is_vector = true;
337 is_x = false;
338 var_name = MDAL::replace( var_name, "vitesse v", "vitesse" );
339 var_name = MDAL::replace( var_name, "suivant y", "" );
340 }
341
342 std::shared_ptr<DatasetGroup> group;
343 auto it = groupsByName.find( var_name );
344 if ( it == groupsByName.end() )
345 {
346 group = std::make_shared< DatasetGroup >(
347 mesh->driverName(),
348 mesh,
349 mesh->uri(),
350 var_name
351 );
352 group->setDataLocation( MDAL_DataLocation::DataOnVertices );
353 group->setIsScalar( !is_vector );
354
355 groupsInOrder.push_back( group );
356 groupsByName[var_name] = group;
357 }
358 else
359 group = it->second;
360
361 group->setReferenceTime( reader->mReferenceTime );
362
363 for ( size_t nT = 0; nT < reader->mTimeSteps.size(); nT++ )
364 {
365 std::shared_ptr<MDAL::DatasetSelafin> dataset;
366 if ( group->datasets.size() > nT )
367 {
368 dataset = std::dynamic_pointer_cast<DatasetSelafin>( group->datasets[nT] );
369 }
370 else
371 {
372 dataset = std::make_shared< DatasetSelafin >( group.get(), reader, nT );
373 dataset->setTime( reader->mTimeSteps.at( nT ) );
374 group->datasets.push_back( dataset );
375 }
376 if ( is_vector )
377 {
378 if ( is_x )
379 {
380 dataset->setXVariableIndex( nName );
381 }
382 else
383 {
384 dataset->setYVariableIndex( nName );
385 }
386 }
387 else
388 {
389 dataset->setXVariableIndex( nName );
390 }
391
392 }
393 }
394
395 // now calculate statistics
396 for ( const std::shared_ptr<DatasetGroup> &group : groupsInOrder )
397 {
398 for ( const std::shared_ptr<Dataset> &dataset : group->datasets )
399 {
400 MDAL::Statistics stats = MDAL::calculateStatistics( dataset );
401 dataset->setStatistics( stats );
402 }
403
404 MDAL::Statistics stats = MDAL::calculateStatistics( group );
405 group->setStatistics( stats );
406 }
407
408 // As everything seems to be ok (no exception thrown), push the groups in the mesh
409 for ( const std::shared_ptr<DatasetGroup> &group : groupsInOrder )
410 mesh->datasetGroups.push_back( group );
411 }
412
readString(size_t len)413 std::string MDAL::SelafinFile::readString( size_t len )
414 {
415 size_t length = readSizeT();
416 if ( length != len ) throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Unable to read string" );
417 std::string ret = readStringWithoutLength( len );
418 ignoreArrayLength();
419 return ret;
420 }
421
readDoubleArr(size_t len)422 std::vector<double> MDAL::SelafinFile::readDoubleArr( size_t len )
423 {
424 size_t length = readSizeT();
425 if ( mStreamInFloatPrecision )
426 {
427 if ( length != len * 4 )
428 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "File format problem while reading double array" );
429 }
430 else
431 {
432 if ( length != len * 8 )
433 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "File format problem while reading double array" );
434 }
435 std::vector<double> ret( len );
436 for ( size_t i = 0; i < len; ++i )
437 {
438 ret[i] = readDouble();
439 }
440 ignoreArrayLength();
441 return ret;
442 }
443
readDoubleArr(const std::streampos & position,size_t offset,size_t len)444 std::vector<double> MDAL::SelafinFile::readDoubleArr( const std::streampos &position, size_t offset, size_t len )
445 {
446 std::vector<double> ret( len );
447 std::streamoff off;
448 if ( mStreamInFloatPrecision )
449 off = offset * 4;
450 else
451 off = offset * 8;
452
453 mIn.seekg( position + off );
454 for ( size_t i = 0; i < len; ++i )
455 ret[i] = readDouble();
456
457 return ret;
458 }
459
readIntArr(size_t len)460 std::vector<int> MDAL::SelafinFile::readIntArr( size_t len )
461 {
462 size_t length = readSizeT();
463 if ( length != len * 4 ) throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "File format problem while reading int array" );
464 std::vector<int> ret( len );
465 for ( size_t i = 0; i < len; ++i )
466 {
467 ret[i] = readInt();
468 }
469 ignoreArrayLength();
470 return ret;
471 }
472
readIntArr(const std::streampos & position,size_t offset,size_t len)473 std::vector<int> MDAL::SelafinFile::readIntArr( const std::streampos &position, size_t offset, size_t len )
474 {
475 std::vector<int> ret( len );
476 std::streamoff off = offset * 4;
477
478 mIn.seekg( position + off );
479 for ( size_t i = 0; i < len; ++i )
480 ret[i] = readInt();
481
482 return ret;
483 }
484
readStringWithoutLength(size_t len)485 std::string MDAL::SelafinFile::readStringWithoutLength( size_t len )
486 {
487 std::vector<char> ptr( len );
488 mIn.read( ptr.data(), static_cast<int>( len ) );
489 if ( !mIn )
490 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Unable to open stream for reading string without length" );
491
492 size_t str_length = 0;
493 for ( size_t i = len; i > 0; --i )
494 {
495 if ( ptr[i - 1] != 0x20 )
496 {
497 str_length = i;
498 break;
499 }
500 }
501 std::string ret( ptr.data(), str_length );
502 return ret;
503 }
504
readDouble()505 double MDAL::SelafinFile::readDouble( )
506 {
507 double ret;
508
509 if ( mStreamInFloatPrecision )
510 {
511 float ret_f;
512 if ( !readValue( ret_f, mIn, mChangeEndianness ) )
513 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Reading double failed" );
514 ret = static_cast<double>( ret_f );
515 }
516 else
517 {
518 if ( !readValue( ret, mIn, mChangeEndianness ) )
519 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Reading double failed" );
520 }
521 return ret;
522 }
523
524
readInt()525 int MDAL::SelafinFile::readInt( )
526 {
527 unsigned char data[4];
528
529 if ( mIn.read( reinterpret_cast< char * >( &data ), 4 ) )
530 if ( !mIn )
531 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Unable to open stream for reading int" );
532 if ( mChangeEndianness )
533 {
534 std::reverse( reinterpret_cast< char * >( &data ), reinterpret_cast< char * >( &data ) + 4 );
535 }
536 int var;
537 memcpy( reinterpret_cast< char * >( &var ),
538 reinterpret_cast< char * >( &data ),
539 4 );
540 return var;
541 }
542
readSizeT()543 size_t MDAL::SelafinFile::readSizeT()
544 {
545 int var = readInt( );
546 return static_cast<size_t>( var );
547 }
548
checkIntArraySize(size_t len)549 bool MDAL::SelafinFile::checkIntArraySize( size_t len )
550 {
551 return ( len * 4 == readSizeT() );
552 }
553
checkDoubleArraySize(size_t len)554 bool MDAL::SelafinFile::checkDoubleArraySize( size_t len )
555 {
556 if ( mStreamInFloatPrecision )
557 {
558 return ( len * 4 ) == readSizeT();
559 }
560 else
561 {
562 return ( len * 8 ) == readSizeT();
563 }
564 }
565
remainingBytes()566 size_t MDAL::SelafinFile::remainingBytes()
567 {
568 if ( mIn.eof() )
569 return 0;
570 return static_cast<size_t>( mFileSize - mIn.tellg() );
571 }
572
passThroughIntArray(size_t size)573 std::streampos MDAL::SelafinFile::passThroughIntArray( size_t size )
574 {
575 std::streampos pos = mIn.tellg();
576 mIn.seekg( size * 4, std::ios_base::cur );
577 ignoreArrayLength();
578 return pos;
579 }
580
passThroughDoubleArray(size_t size)581 std::streampos MDAL::SelafinFile::passThroughDoubleArray( size_t size )
582 {
583 std::streampos pos = mIn.tellg();
584 if ( mStreamInFloatPrecision )
585 size *= 4;
586 else
587 size *= 8;
588
589 mIn.seekg( size, std::ios_base::cur );
590 ignoreArrayLength();
591 return pos;
592 }
593
ignore(int len)594 void MDAL::SelafinFile::ignore( int len )
595 {
596 mIn.ignore( len );
597 if ( !mIn )
598 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Unable to ignore characters (invalid stream)" );
599 }
600
ignoreArrayLength()601 void MDAL::SelafinFile::ignoreArrayLength( )
602 {
603 ignore( 4 );
604 }
605
606 // //////////////////////////////
607 // DRIVER
608 // //////////////////////////////
609
DriverSelafin()610 MDAL::DriverSelafin::DriverSelafin():
611 Driver( "SELAFIN",
612 "Selafin File",
613 "*.slf;;*.ser;;*.geo;;*.res",
614 Capability::ReadMesh | Capability::SaveMesh | Capability::WriteDatasetsOnVertices | Capability::ReadDatasets
615 )
616 {
617 }
618
619 MDAL::DriverSelafin::~DriverSelafin() = default;
620
create()621 MDAL::DriverSelafin *MDAL::DriverSelafin::create()
622 {
623 return new DriverSelafin();
624 }
625
canReadMesh(const std::string & uri)626 bool MDAL::DriverSelafin::canReadMesh( const std::string &uri )
627 {
628 if ( !MDAL::fileExists( uri ) ) return false;
629
630 try
631 {
632 SelafinFile file( uri );
633 file.parseMeshFrame();
634 return true;
635 }
636 catch ( ... )
637 {
638 return false;
639 }
640 }
641
canReadDatasets(const std::string & uri)642 bool MDAL::DriverSelafin::canReadDatasets( const std::string &uri )
643 {
644 if ( !MDAL::fileExists( uri ) ) return false;
645
646 try
647 {
648 SelafinFile file( uri );
649 file.parseMeshFrame();
650 return true;
651 }
652 catch ( ... )
653 {
654 return false;
655 }
656 }
657
load(const std::string & meshFile,const std::string &)658 std::unique_ptr<MDAL::Mesh> MDAL::DriverSelafin::load( const std::string &meshFile, const std::string & )
659 {
660 MDAL::Log::resetLastStatus();
661 std::unique_ptr<Mesh> mesh;
662
663 try
664 {
665 mesh = SelafinFile::createMesh( meshFile );
666 }
667 catch ( MDAL_Status error )
668 {
669 MDAL::Log::error( error, name(), "Error while loading file " + meshFile );
670 mesh.reset();
671 }
672 catch ( MDAL::Error err )
673 {
674 MDAL::Log::error( err, name() );
675 mesh.reset();
676 }
677 return mesh;
678 }
679
load(const std::string & datFile,MDAL::Mesh * mesh)680 void MDAL::DriverSelafin::load( const std::string &datFile, MDAL::Mesh *mesh )
681 {
682 MDAL::Log::resetLastStatus();
683
684 try
685 {
686 SelafinFile::populateDataset( mesh, datFile );
687 }
688 catch ( MDAL_Status error )
689 {
690 MDAL::Log::error( error, name(), "Error while loading dataset in file " + datFile );
691 }
692 catch ( MDAL::Error err )
693 {
694 MDAL::Log::error( err, name() );
695 }
696 }
697
persist(MDAL::DatasetGroup * group)698 bool MDAL::DriverSelafin::persist( MDAL::DatasetGroup *group )
699 {
700 if ( !group || ( group->dataLocation() != MDAL_DataLocation::DataOnVertices ) )
701 {
702 MDAL::Log::error( MDAL_Status::Err_IncompatibleDataset, name(), "Selafin can store only 2D vertices datasets" );
703 return true;
704 }
705
706 try
707 {
708 saveDatasetGroupOnFile( group );
709 return false;
710 }
711 catch ( MDAL::Error err )
712 {
713 MDAL::Log::error( err, name() );
714 return true;
715 }
716 }
717
saveDatasetGroupOnFile(MDAL::DatasetGroup * datasetGroup)718 bool MDAL::DriverSelafin::saveDatasetGroupOnFile( MDAL::DatasetGroup *datasetGroup )
719 {
720 const std::string &fileName = datasetGroup->uri();
721
722 if ( ! MDAL::fileExists( fileName ) )
723 {
724 //create a new mesh file
725 save( fileName, "", datasetGroup->mesh() );
726
727 if ( ! MDAL::fileExists( fileName ) )
728 throw MDAL::Error( MDAL_Status::Err_FailToWriteToDisk, "Unable to create new file" );
729 }
730
731 SelafinFile file( fileName );
732 return file.addDatasetGroup( datasetGroup );
733 }
734
735 // //////////////////////////////
736 // MeshSelafin
737 // //////////////////////////////
738
MeshSelafin(const std::string & uri,std::shared_ptr<MDAL::SelafinFile> reader)739 MDAL::MeshSelafin::MeshSelafin( const std::string &uri,
740 std::shared_ptr<MDAL::SelafinFile> reader ):
741 Mesh( "SELAFIN", reader->verticesPerFace(), uri )
742 , mReader( reader )
743 {}
744
readVertices()745 std::unique_ptr<MDAL::MeshVertexIterator> MDAL::MeshSelafin::readVertices()
746 {
747 return std::unique_ptr<MDAL::MeshVertexIterator>( new MeshSelafinVertexIterator(
748 mReader ) );
749 }
750
readEdges()751 std::unique_ptr<MDAL::MeshEdgeIterator> MDAL::MeshSelafin::readEdges()
752 {
753 return std::unique_ptr<MDAL::MeshEdgeIterator>();
754 }
755
readFaces()756 std::unique_ptr<MDAL::MeshFaceIterator> MDAL::MeshSelafin::readFaces()
757 {
758 return std::unique_ptr<MDAL::MeshFaceIterator>(
759 new MeshSelafinFaceIterator( mReader ) );
760 }
761
extent() const762 MDAL::BBox MDAL::MeshSelafin::extent() const
763 {
764 if ( mIsExtentUpToDate )
765 return mExtent;
766 calculateExtent();
767 return mExtent;
768 }
769
closeSource()770 void MDAL::MeshSelafin::closeSource()
771 {
772 if ( mReader )
773 {
774 mReader->mIn.close();
775 mReader->mParsed = false;
776 }
777 }
778
calculateExtent() const779 void MDAL::MeshSelafin::calculateExtent() const
780 {
781 size_t bufferSize = 1000;
782 std::unique_ptr<MeshSelafinVertexIterator> vertexIt(
783 new MeshSelafinVertexIterator( mReader ) );
784 size_t count = 0;
785 BBox bbox;
786 std::vector<Vertex> vertices( mReader->verticesCount() );
787 size_t index = 0;
788 do
789 {
790 std::vector<double> verticesCoord( bufferSize * 3 );
791 count = vertexIt->next( bufferSize, verticesCoord.data() );
792
793 for ( size_t i = 0; i < count; ++i )
794 {
795 vertices[index + i].x = verticesCoord.at( i * 3 );
796 vertices[index + i].y = verticesCoord.at( i * 3 + 1 );
797 vertices[index + i].z = verticesCoord.at( i * 3 + 2 );
798 }
799 index += count;
800 }
801 while ( count != 0 );
802
803 mExtent = MDAL::computeExtent( vertices );
804 mIsExtentUpToDate = true;
805 }
806
MeshSelafinVertexIterator(std::shared_ptr<MDAL::SelafinFile> reader)807 MDAL::MeshSelafinVertexIterator::MeshSelafinVertexIterator( std::shared_ptr<MDAL::SelafinFile> reader ):
808 mReader( reader )
809 {}
810
next(size_t vertexCount,double * coordinates)811 size_t MDAL::MeshSelafinVertexIterator::next( size_t vertexCount, double *coordinates )
812 {
813 size_t count = std::min( vertexCount, mReader->verticesCount() - mPosition );
814
815 if ( count == 0 )
816 return 0;
817
818 std::vector<double> coord = mReader->vertices( mPosition, count );
819
820 memcpy( coordinates, coord.data(), count * 24 );
821
822 mPosition += count;
823
824 return count;
825 }
826
MeshSelafinFaceIterator(std::shared_ptr<MDAL::SelafinFile> reader)827 MDAL::MeshSelafinFaceIterator::MeshSelafinFaceIterator( std::shared_ptr<MDAL::SelafinFile> reader ):
828 mReader( reader )
829 {}
830
next(size_t faceOffsetsBufferLen,int * faceOffsetsBuffer,size_t vertexIndicesBufferLen,int * vertexIndicesBuffer)831 size_t MDAL::MeshSelafinFaceIterator::next( size_t faceOffsetsBufferLen, int *faceOffsetsBuffer, size_t vertexIndicesBufferLen, int *vertexIndicesBuffer )
832 {
833 assert( faceOffsetsBuffer );
834 assert( vertexIndicesBuffer );
835 assert( mReader->verticesPerFace() != 0 );
836
837 const size_t verticesPerFace = mReader->verticesPerFace();
838 size_t count = std::min( faceOffsetsBufferLen, mReader->facesCount() - mPosition );
839
840 count = std::min( count, vertexIndicesBufferLen / verticesPerFace );
841
842 if ( count == 0 )
843 return 0;
844
845 std::vector<int> indexes = mReader->connectivityIndex( mPosition * verticesPerFace, count * verticesPerFace );
846
847 if ( indexes.size() != count * verticesPerFace )
848 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "File format problem while reading faces" );
849
850 int vertexLocalIndex = 0;
851
852 for ( size_t i = 0; i < count; i++ )
853 {
854 for ( size_t j = 0; j < verticesPerFace; ++j )
855 {
856 if ( size_t( indexes[j + i * verticesPerFace] ) > mReader->verticesCount() )
857 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "File format problem while reading faces" );
858 vertexIndicesBuffer[vertexLocalIndex + j] = indexes[j + i * verticesPerFace] - 1;
859 }
860 vertexLocalIndex += MDAL::toInt( verticesPerFace );
861 faceOffsetsBuffer[i] = vertexLocalIndex;
862 }
863
864 mPosition += count;
865
866 return count;
867
868 }
869
DatasetSelafin(MDAL::DatasetGroup * parent,std::shared_ptr<MDAL::SelafinFile> reader,size_t timeStepIndex)870 MDAL::DatasetSelafin::DatasetSelafin( MDAL::DatasetGroup *parent,
871 std::shared_ptr<MDAL::SelafinFile> reader, size_t timeStepIndex ):
872 Dataset2D( parent )
873 , mReader( reader )
874 , mTimeStepIndex( timeStepIndex )
875 {
876 }
877
scalarData(size_t indexStart,size_t count,double * buffer)878 size_t MDAL::DatasetSelafin::scalarData( size_t indexStart, size_t count, double *buffer )
879 {
880 count = std::min( mReader->verticesCount() - indexStart, count );
881 std::vector<double> values = mReader->datasetValues( mTimeStepIndex, mXVariableIndex, indexStart, count );
882 if ( values.size() != count )
883 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "File format problem while reading dataset value" );
884
885 memcpy( buffer, values.data(), count * 8 );
886
887 return count;
888 }
889
vectorData(size_t indexStart,size_t count,double * buffer)890 size_t MDAL::DatasetSelafin::vectorData( size_t indexStart, size_t count, double *buffer )
891 {
892 count = std::min( mReader->verticesCount() - indexStart, count );
893 std::vector<double> xValues = mReader->datasetValues( mTimeStepIndex, mXVariableIndex, indexStart, count );
894 std::vector<double> yValues = mReader->datasetValues( mTimeStepIndex, mYVariableIndex, indexStart, count );
895
896 if ( xValues.size() != count || yValues.size() != count )
897 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "File format problem while reading dataset value" );
898
899 for ( size_t i = 0; i < count; ++i )
900 {
901 buffer[2 * i] = xValues[i];
902 buffer[2 * i + 1] = yValues[i];
903 }
904
905 return count;
906 }
907
setXVariableIndex(size_t index)908 void MDAL::DatasetSelafin::setXVariableIndex( size_t index )
909 {
910 mXVariableIndex = index;
911 }
912
setYVariableIndex(size_t index)913 void MDAL::DatasetSelafin::setYVariableIndex( size_t index )
914 {
915 mYVariableIndex = index;
916 }
917
918 template<typename T>
writeValue(std::ofstream & file,T value)919 static void writeValue( std::ofstream &file, T value )
920 {
921 bool isNativeLittleEndian = MDAL::isNativeLittleEndian();
922 MDAL::writeValue( value, file, isNativeLittleEndian );
923 }
924
writeInt(std::ofstream & file,int i)925 static void writeInt( std::ofstream &file, int i )
926 {
927 bool isNativeLittleEndian = MDAL::isNativeLittleEndian();
928 MDAL::writeValue( i, file, isNativeLittleEndian );
929 }
930
writeStringRecord(std::ofstream & file,const std::string & str)931 static void writeStringRecord( std::ofstream &file, const std::string &str )
932 {
933 writeInt( file, MDAL::toInt( str.size() ) );
934 file.write( str.data(), str.size() );
935 writeInt( file, MDAL::toInt( str.size() ) );
936 }
937
938 template<typename T>
writeValueArrayRecord(std::ofstream & file,const std::vector<T> & array)939 static void writeValueArrayRecord( std::ofstream &file, const std::vector<T> &array )
940 {
941 writeValue( file, int( array.size()*sizeof( T ) ) );
942 for ( const T value : array )
943 writeValue( file, value );
944 writeValue( file, int( array.size()*sizeof( T ) ) );
945 }
946
947 template<typename T>
writeValueArray(std::ofstream & file,const std::vector<T> & array)948 static void writeValueArray( std::ofstream &file, const std::vector<T> &array )
949 {
950 for ( const T value : array )
951 writeValue( file, value );
952 }
953
954 template<typename T>
writeVertices(std::ofstream & file,MDAL::Mesh * mesh)955 static void writeVertices( std::ofstream &file, MDAL::Mesh *mesh )
956 {
957 std::unique_ptr<MDAL::MeshVertexIterator> vertexIter = mesh->readVertices();
958 size_t verticesCount = mesh->verticesCount();
959 std::vector<T> xValues( verticesCount );
960 std::vector<T> yValues( verticesCount );
961 size_t count = 0;
962 size_t vertexIndex = 0;
963 size_t bufferSize = BUFFER_SIZE;
964 do
965 {
966 std::vector<double> coordinates( bufferSize * 3 );
967 count = vertexIter->next( bufferSize, coordinates.data() );
968 for ( size_t i = 0; i < count; ++i )
969 {
970 xValues[vertexIndex + i] = coordinates[i * 3];
971 yValues[vertexIndex + i] = coordinates[i * 3 + 1];
972 }
973 vertexIndex += count;
974 }
975 while ( count != 0 );
976 writeValueArrayRecord( file, xValues );
977 writeValueArrayRecord( file, yValues );
978 }
979
save(const std::string & fileName,const std::string &,MDAL::Mesh * mesh)980 void MDAL::DriverSelafin::save( const std::string &fileName, const std::string &, MDAL::Mesh *mesh )
981 {
982 std::ofstream file = MDAL::openOutputFile( fileName.c_str(), std::ofstream::binary );
983
984 std::string header( "Selafin file created by MDAL library" );
985 std::string remainingStr( " ", 72 - header.size() );
986 header.append( remainingStr );
987 header.append( "SERAFIND" );
988 assert( header.size() == 80 );
989 writeStringRecord( file, header );
990
991 // NBV(1) NBV(2) size
992 std::vector<int> nbvSize( 2 );
993 nbvSize[0] = 0;
994 nbvSize[1] = 0;
995 writeValueArrayRecord( file, nbvSize );
996
997 //don't write variable name
998
999 //parameter table, all values are 0
1000 std::vector<int> param( 10, 0 );
1001 writeValueArrayRecord( file, param );
1002
1003 //NELEM,NPOIN,NDP,1
1004 size_t verticesPerFace = mesh->faceVerticesMaximumCount();
1005 size_t verticesCount = mesh->verticesCount();
1006 size_t facesCount = mesh->facesCount();
1007 std::vector<int> elem( 4 );
1008 elem[0] = MDAL::toInt( facesCount );
1009 elem[1] = MDAL::toInt( verticesCount );
1010 elem[2] = MDAL::toInt( verticesPerFace );
1011 elem[3] = 1;
1012 writeValueArrayRecord( file, elem );
1013
1014 //connectivity table
1015 int bufferSize = BUFFER_SIZE;
1016 std::vector<int> faceOffsetBuffer( bufferSize );
1017 std::unique_ptr<MeshFaceIterator> faceIter = mesh->readFaces();
1018 size_t count = 0;
1019 writeInt( file, MDAL::toInt( facesCount * verticesPerFace * 4 ) );
1020 if ( facesCount > 0 )
1021 {
1022 do
1023 {
1024 std::vector<int> inkle( bufferSize * verticesPerFace );
1025 count = faceIter->next( bufferSize, faceOffsetBuffer.data(), bufferSize * verticesPerFace, inkle.data() );
1026 inkle.resize( count * verticesPerFace );
1027 for ( size_t i = 0; i < inkle.size(); ++i )
1028 inkle[i]++;
1029
1030 writeValueArray( file, inkle );
1031 }
1032 while ( count != 0 );
1033 }
1034 writeInt( file, MDAL::toInt( facesCount * verticesPerFace * 4 ) );
1035
1036 // IPOBO filled with 0
1037 writeValueArrayRecord( file, std::vector<int>( verticesCount, 0 ) );
1038
1039 //Vertices
1040 writeVertices<double>( file, mesh );
1041
1042 file.close();
1043 }
1044
writeDatasetOnFileSuffix() const1045 std::string MDAL::DriverSelafin::writeDatasetOnFileSuffix() const
1046 {
1047 return "slf";
1048 }
1049
saveMeshOnFileSuffix() const1050 std::string MDAL::DriverSelafin::saveMeshOnFileSuffix() const
1051 {
1052 return "slf";
1053 }
1054
1055 // return false if fails
streamToStream(std::ostream & destination,std::ifstream & source,std::streampos sourceStartPosition,std::streamoff len,std::streamoff maxBufferSize)1056 static void streamToStream( std::ostream &destination,
1057 std::ifstream &source,
1058 std::streampos sourceStartPosition,
1059 std::streamoff len,
1060 std::streamoff maxBufferSize )
1061 {
1062 assert( maxBufferSize != 0 );
1063 std::streampos position = sourceStartPosition;
1064 source.seekg( 0, source.end );
1065 std::streampos end = std::min( source.tellg(), sourceStartPosition + len );
1066 source.seekg( sourceStartPosition );
1067
1068 while ( position < end )
1069 {
1070 size_t bufferSize = std::min( maxBufferSize, end - position );
1071 std::vector<char> buffer( bufferSize );
1072 source.read( &buffer[0], bufferSize );
1073 destination.write( &buffer[0], bufferSize );
1074 position += bufferSize;
1075 }
1076 }
1077
writeScalarDataset(std::ofstream & file,MDAL::Dataset * dataset,bool isFloat)1078 static void writeScalarDataset( std::ofstream &file, MDAL::Dataset *dataset, bool isFloat )
1079 {
1080 size_t valuesCount = dataset->valuesCount();
1081 size_t bufferSize = BUFFER_SIZE;
1082 size_t count = 0;
1083 size_t indexStart = 0;
1084 writeInt( file, MDAL::toInt( valuesCount * ( isFloat ? 4 : 8 ) ) );
1085 do
1086 {
1087 std::vector<double> values( bufferSize );
1088 count = dataset->scalarData( indexStart, bufferSize, values.data() );
1089 values.resize( count );
1090 if ( isFloat )
1091 {
1092 std::vector<float> floatValues( count );
1093 for ( int i = 0; i < MDAL::toInt( count ); ++i )
1094 floatValues[i] = static_cast<float>( values[i] );
1095 writeValueArray( file, floatValues );
1096 }
1097 else
1098 writeValueArray( file, values );
1099
1100 indexStart += count;
1101 }
1102 while ( count != 0 );
1103 writeInt( file, MDAL::toInt( valuesCount * ( isFloat ? 4 : 8 ) ) );
1104 }
1105
1106 template<typename T>
writeVectorDataset(std::ofstream & file,MDAL::Dataset * dataset)1107 static void writeVectorDataset( std::ofstream &file, MDAL::Dataset *dataset )
1108 {
1109 size_t valuesCount = dataset->valuesCount();
1110
1111 std::vector<T> valuesX( valuesCount );
1112 std::vector<T> valuesY( valuesCount );
1113 size_t bufferSize = BUFFER_SIZE;
1114 size_t count = 0;
1115 size_t indexStart = 0;
1116 do
1117 {
1118 std::vector<double> values( bufferSize * 2 );
1119 count = dataset->vectorData( indexStart, bufferSize, values.data() );
1120 values.resize( count * 2 );
1121 for ( int i = 0; i < MDAL::toInt( count ); ++i )
1122 {
1123 valuesX[indexStart + i] = static_cast< T >( values[i * 2] );
1124 valuesY[indexStart + i] = static_cast< T >( values[i * 2 + 1] );
1125 }
1126 indexStart += count;
1127 }
1128 while ( count != 0 );
1129 writeValueArrayRecord( file, valuesX );
1130 writeValueArrayRecord( file, valuesY );
1131 }
1132
addDatasetGroup(MDAL::DatasetGroup * datasetGroup)1133 bool MDAL::SelafinFile::addDatasetGroup( MDAL::DatasetGroup *datasetGroup )
1134 {
1135 // Create a new file with same data but with another datasetGroup
1136 initialize();
1137
1138 if ( !mIn.is_open() )
1139 return false;
1140
1141 parseFile();
1142
1143 size_t realSize;
1144 if ( mStreamInFloatPrecision )
1145 realSize = 4;
1146 else
1147 realSize = 8;
1148
1149 std::string tempFileName = mFileName;
1150 tempFileName.append( ".tmp" );
1151
1152 mIn.seekg( 0, mIn.beg );
1153 std::ofstream out = MDAL::openOutputFile( tempFileName, std::ios_base::binary );
1154 if ( ! out.is_open() )
1155 throw MDAL::Error( MDAL_Status::Err_FailToWriteToDisk, "Unable to add dataset in file" );
1156
1157 //write the same header
1158 writeStringRecord( out, readHeader() );
1159
1160 //Read the NBV1//NBV2 size, and add 1 to NBV1
1161 std::vector<int> nbv = readIntArr( 2 );
1162
1163 int addedVariable = 1;
1164 if ( !datasetGroup->isScalar() )
1165 addedVariable = 2;
1166
1167 nbv[0] = nbv[0] + addedVariable;
1168 writeValueArrayRecord( out, nbv );
1169
1170 // write pre-existing dataset name
1171 for ( size_t i = 0; i < mVariableNames.size(); ++i )
1172 {
1173 std::string variableName = mVariableNames.at( i );
1174 variableName.resize( 32, ' ' );
1175 writeStringRecord( out, variableName );
1176 }
1177
1178 // write new(s) variable name
1179 std::string datasetGroupName = datasetGroup->name();
1180 if ( datasetGroup->isScalar() )
1181 {
1182 datasetGroupName.resize( 32, ' ' );
1183 writeStringRecord( out, datasetGroupName );
1184 }
1185 else
1186 {
1187 if ( datasetGroupName.size() > 27 )
1188 datasetGroupName.resize( 27 );
1189 std::string xName = datasetGroupName + " along x";
1190 std::string yName = datasetGroupName + " along y";
1191 xName.resize( 32 );
1192 yName.resize( 32 );
1193 writeStringRecord( out, xName );
1194 writeStringRecord( out, yName );
1195 }
1196
1197 //check if valid reference time
1198 if ( !mReferenceTime.isValid() && datasetGroup->referenceTime().isValid() )
1199 mReferenceTime = datasetGroup->referenceTime();
1200
1201 //handlle time step
1202 if ( mVariableNames.empty() ) //if no variable name, no valid time step before
1203 {
1204 mTimeSteps = std::vector<RelativeTimestamp>( datasetGroup->datasets.size() );
1205 for ( size_t i = 0; i < datasetGroup->datasets.size(); ++i )
1206 {
1207 mTimeSteps[i] = RelativeTimestamp( datasetGroup->datasets.at( i )->timestamp() );
1208 }
1209 }
1210 else
1211 {
1212 if ( datasetGroup->datasets.size() != mTimeSteps.size() )
1213 throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Incompatible dataset group : time steps count are not the same" );
1214 }
1215
1216 //parameters table
1217 mParameters[9] = mReferenceTime.isValid() ? 1 : 0;
1218 writeValueArrayRecord( out, mParameters );
1219
1220 if ( mReferenceTime.isValid() )
1221 {
1222 writeValueArrayRecord( out, mReferenceTime.expandToCalendarArray() );
1223 }
1224
1225 //elems count
1226 writeValueArrayRecord( out, std::vector<int> {int( mFacesCount ), int( mVerticesCount ), int( mVerticesPerFace ), 1} );
1227
1228 //IKLE
1229 writeInt( out, MDAL::toInt( mFacesCount * mVerticesPerFace * 4 ) );
1230 streamToStream( out, mIn, mConnectivityStreamPosition, mFacesCount * mVerticesPerFace * 4, BUFFER_SIZE );
1231 writeInt( out, MDAL::toInt( mFacesCount * mVerticesPerFace * 4 ) );
1232 //vertices
1233
1234 //IPOBO
1235 writeInt( out, MDAL::toInt( mVerticesCount * 4 ) );
1236 streamToStream( out, mIn, mIPOBOStreamPosition, mVerticesCount * 4, BUFFER_SIZE );
1237 writeInt( out, MDAL::toInt( mVerticesCount * 4 ) );
1238
1239 //X Vertices
1240 writeInt( out, MDAL::toInt( mVerticesCount * realSize ) );
1241 streamToStream( out, mIn, mXStreamPosition, mVerticesCount * realSize, BUFFER_SIZE );
1242 writeInt( out, MDAL::toInt( mVerticesCount * realSize ) );
1243 //Y Vertices
1244 writeInt( out, MDAL::toInt( mVerticesCount * realSize ) );
1245 streamToStream( out, mIn, mYStreamPosition, mVerticesCount * realSize, BUFFER_SIZE );
1246 writeInt( out, MDAL::toInt( mVerticesCount * realSize ) );
1247
1248 // Write datasets
1249 for ( size_t nT = 0; nT < mTimeSteps.size(); nT++ )
1250 {
1251 // Time step
1252 if ( mStreamInFloatPrecision )
1253 {
1254 std::vector<float> time( 1, static_cast<float>( mTimeSteps.at( nT ).value( RelativeTimestamp::seconds ) ) );
1255 writeValueArrayRecord( out, time );
1256 }
1257 else
1258 {
1259 std::vector<double> time( 1, mTimeSteps.at( nT ).value( RelativeTimestamp::seconds ) );
1260 writeValueArrayRecord( out, time );
1261 }
1262
1263 // First, prexisting datasets from the original file
1264 for ( int i = 0; i < nbv[0] - addedVariable; ++i )
1265 {
1266 writeInt( out, MDAL::toInt( mVerticesCount * realSize ) );
1267 streamToStream( out, mIn, mVariableStreamPosition[i][nT], realSize * mVerticesCount, BUFFER_SIZE );
1268 writeInt( out, MDAL::toInt( mVerticesCount * realSize ) );
1269 }
1270
1271 // Then, new datasets from the new dataset group
1272 Dataset *dataset = datasetGroup->datasets[nT].get();
1273 if ( datasetGroup->isScalar() )
1274 {
1275 writeScalarDataset( out, dataset, mStreamInFloatPrecision );
1276 }
1277 else
1278 {
1279 if ( mStreamInFloatPrecision )
1280 writeVectorDataset<float>( out, dataset );
1281 else
1282 writeVectorDataset<double>( out, dataset );
1283 }
1284 }
1285
1286 out.close();
1287 mIn.close();
1288
1289 // if the uri of the dataset group is the same than the file name, be sure to close it before replace it
1290 if ( datasetGroup->uri() == mFileName )
1291 datasetGroup->mesh()->closeSource();
1292
1293 if ( std::remove( mFileName.c_str() ) != 0 )
1294 {
1295 std::remove( tempFileName.c_str() );
1296 throw MDAL::Error( MDAL_Status::Err_FailToWriteToDisk, "Unable to write dataset in file" );
1297 }
1298
1299 std::rename( tempFileName.c_str(), mFileName.c_str() );
1300
1301 parseFile();
1302
1303 return true;
1304 }
1305