TTK
Loading...
Searching...
No Matches
TopologicalCompression.h
Go to the documentation of this file.
1
13
14#pragma once
15
16// base code includes
17#include <FTMTreePP.h>
19#include <Triangulation.h>
20
21// std
22#include <cstring>
23#include <stack>
24#include <type_traits>
25
26namespace ttk {
27
28 enum class CompressionType { PersistenceDiagram = 0, Other = 1 };
29
30 class TopologicalCompression : virtual public Debug {
31
32 public:
33 // Base code methods.
35
44 template <class dataType,
45 typename triangulationType = AbstractTriangulation>
46 int execute(const dataType *const inputData,
47 const SimplexId *const inputOffsets,
48 dataType *outputData,
49 const triangulationType &triangulation);
50
51 // Persistence compression methods.
52 template <class dataType, typename triangulationType>
54 std::vector<std::tuple<SimplexId, SimplexId, dataType>> &JTPairs,
55 std::vector<std::tuple<SimplexId, SimplexId, dataType>> &STPairs,
56 const dataType *const inputScalars_,
57 const SimplexId *const inputOffsets,
58 const triangulationType &triangulation);
59 template <typename dataType, typename triangulationType>
60 int compressForPersistenceDiagram(int vertexNumber,
61 const dataType *const inputData,
62 const SimplexId *const inputOffset,
63 dataType *outputData,
64 const double &tol,
65 const triangulationType &triangulation);
66
67 // Other compression methods.
68 int computeOther() const;
69
70 template <typename dataType>
71 int compressForOther(int vertexNumber,
72 const dataType *const inputData,
73 const SimplexId *const inputOffsets,
74 dataType *outputData,
75 const double &tol) const;
76
77 // Getters and setters.
78 inline void setCompressionType(int compressionType) {
79 compressionType_ = compressionType;
80 }
81 inline void setSQ(const std::string &sqMethod) {
82 SQMethod = sqMethod;
83 }
84 inline void setZFPOnly(bool z) {
85 ZFPOnly = z;
86 }
87 inline void setSubdivide(bool b) {
88 Subdivide = b;
89 }
90 inline void setMaximumError(double maximumError) {
91 MaximumError = maximumError;
92 }
93 inline void setTolerance(const double data) {
94 Tolerance = data;
95 }
96 inline void
97 setUseTopologicalSimplification(bool useTopologicalSimplification) {
98 UseTopologicalSimplification = useTopologicalSimplification;
99 }
100 inline void setFileName(char *fn) {
101 fileName = fn;
102 }
103 inline void
105 if(triangulation != nullptr) {
106 triangulation->preconditionVertexNeighbors();
108 ftmTreePP.preconditionTriangulation(triangulation, false);
109 }
110 }
111
112 inline int getNbVertices() const {
113 return NbVertices;
114 }
115 inline int getNbSegments() const {
116 return NbSegments;
117 }
118 inline std::vector<int> &getSegmentation() {
119 return segmentation_;
120 }
121 inline std::vector<std::tuple<double, int>> &getMapping() {
122 return mapping_;
123 }
124 inline std::vector<std::tuple<int, double, int>> &getCriticalConstraints() {
126 }
127
128 inline int getCompressionType() const {
129 return compressionType_;
130 }
131 inline int getSQMethod() const {
132 return SQMethodInt;
133 }
134 inline int getDataScalarType() const {
135 return dataScalarType_;
136 }
137 inline int *getDataExtent() {
138 return dataExtent_;
139 }
140 inline double *getDataSpacing() {
141 return dataSpacing_;
142 }
143 inline double *getDataOrigin() {
144 return dataOrigin_;
145 }
146 inline double getTolerance() const {
147 return Tolerance;
148 }
149 inline double getZFPTolerance() const {
150 return ZFPTolerance;
151 }
156 inline void relToAbsZFPTolerance(const double zfpRelTol,
157 const std::array<double, 2> &sfRange) {
158 this->ZFPTolerance = zfpRelTol * (sfRange[1] - sfRange[0]) / 100.0;
159 }
160 inline bool getZFPOnly() const {
161 return ZFPOnly;
162 }
163
164 inline const std::vector<char> &getDataArrayName() const {
165 return dataArrayName_;
166 }
167 inline std::vector<double> &getDecompressedData() {
168 return decompressedData_;
169 }
170 inline std::vector<SimplexId> &getDecompressedOffsets() {
172 }
173 inline std::vector<SimplexId> &getCompressedOffsets() {
174 return compressedOffsets_;
175 }
176
177 // IO management.
178 static unsigned int log2(int val);
179 inline static bool cmp(const std::tuple<double, int> &a,
180 const std::tuple<double, int> &b) {
181 return std::get<1>(a) > std::get<1>(b);
182 }
183 inline static bool cmp2(const std::tuple<double, int> &a,
184 const std::tuple<double, int> &b) {
185 return std::get<1>(a) < std::get<1>(b);
186 }
187
188 template <typename T>
189 T Read(FILE *fm) const {
190 static_assert(std::is_same<T, uint8_t>() || std::is_same<T, int32_t>()
191 || std::is_same<T, uint64_t>()
192 || std::is_same<T, double>(),
193 "Function should have only those types");
194 T ret;
195 const auto status = std::fread(&ret, sizeof(T), 1, fm);
196 if(status == 0) {
197 this->printErr("Error reading " + std::string(typeid(T).name()) + "!");
198 }
199 return ret;
200 }
201 template <typename T>
202 void ReadByteArray(FILE *fm, T *buffer, size_t length) const {
203 const auto status = std::fread(buffer, sizeof(T), length, fm);
204 if(status == 0) {
205 this->printErr("Error reading " + std::string(typeid(T).name())
206 + "array!");
207 }
208 }
209 template <typename T>
210 void Write(FILE *fm, T data) const {
211 static_assert(std::is_same<T, uint8_t>() || std::is_same<T, int32_t>()
212 || std::is_same<T, uint64_t>()
213 || std::is_same<T, double>(),
214 "Function should have only those types");
215 auto status = std::fwrite(&data, sizeof(T), 1, fm);
216 if(status == 0) {
217 this->printErr("Error writing " + std::string(typeid(T).name()) + "!");
218 }
219 }
220 template <typename T>
221 void WriteByteArray(FILE *fm, const T *buffer, size_t length) const {
222 const auto status = std::fwrite(buffer, sizeof(T), length, fm);
223 if(status == 0) {
224 this->printErr("Error writing " + std::string(typeid(T).name())
225 + "array!");
226 }
227 }
228
229 int ReadCompactSegmentation(FILE *fm,
230 std::vector<int> &segmentation,
231 int &numberOfVertices,
232 int &numberOfSegments) const;
234 FILE *fm,
235 std::vector<std::tuple<double, int>> &mappings,
236 std::vector<std::tuple<double, int>> &mappingsSortedPerValue,
237 std::vector<std::tuple<int, double, int>> &constraints,
238 double &min,
239 double &max,
240 int &nbConstraints) const;
241
242 int ReadMetaData(FILE *fm);
243 template <typename triangulationType>
244 int ReadFromFile(FILE *fm, const triangulationType &triangulation);
245
246 int WriteCompactSegmentation(FILE *fm,
247 const std::vector<int> &segmentation,
248 int numberOfVertices,
249 int numberOfSegments) const;
251 FILE *fm,
252 std::vector<std::tuple<double, int>> &mapping,
253 std::vector<std::tuple<int, double, int>> &constraints) const;
254
255 int WriteMetaData(FILE *fp,
256 int compressionType,
257 bool zfpOnly,
258 const char *sqMethod,
259 int dataType,
260 int *dataExtent,
261 double *dataSpacing,
262 double *dataOrigin,
263 double tolerance,
264 double zfpTolerance,
265 const std::string &dataArrayName);
266
267 int WriteToFile(FILE *fp,
268 int compressionType,
269 bool zfpOnly,
270 const char *sqMethod,
271 int dataType,
272 int *dataExtent,
273 double *dataSpacing,
274 double *dataOrigin,
275 double *data,
276 double tolerance,
277 double zfpTolerance,
278 const std::string &dataArrayName);
279
280 template <typename dataType>
281 void CropIntervals(
282 std::vector<std::tuple<dataType, int>> &mappings,
283 std::vector<std::tuple<dataType, int>> &mappingsSortedPerValue,
284 double min,
285 double max,
286 int vertexNumber,
287 double *array,
288 std::vector<int> &Seg) const;
289
290 // API management.
291
292#ifdef TTK_ENABLE_ZFP
293 int CompressWithZFP(FILE *file,
294 const bool decompress,
295 std::vector<double> &array,
296 const int nx,
297 const int ny,
298 const int nz,
299 const double zfpTolerance) const;
300#endif
301
302#ifdef TTK_ENABLE_ZLIB
303 unsigned long GetZlibDestLen(const unsigned long sourceLen) const;
304 void CompressWithZlib(bool decompress,
305 unsigned char *dest,
306 unsigned long &destLen,
307 const unsigned char *const source,
308 const unsigned long sourceLen) const;
309#endif
310
311 private:
312 // Internal read/write.
313
314 int ComputeTotalSizeForOther() const;
315
316 int ComputeTotalSizeForPersistenceDiagram(
317 std::vector<std::tuple<double, int>> &mapping,
318 std::vector<std::tuple<int, double, int>> &criticalConstraints,
319 bool zfpOnly,
320 int nbSegments,
321 int nbVertices,
322 double zfpTolerance) const;
323
324 int ReadPersistenceTopology(FILE *fm);
325 int ReadOtherTopology(FILE *fm) const;
326 template <typename triangulationType>
327 int ReadPersistenceGeometry(FILE *fm,
328 const triangulationType &triangulation);
329 int ReadOtherGeometry(FILE *fm) const;
330
331 int WritePersistenceTopology(FILE *fm);
332 int WriteOtherTopology(FILE *fm) const;
333 int WritePersistenceGeometry(FILE *fm,
334 int *dataExtent,
335 bool zfpOnly,
336 double zfpTolerance,
337 double *toCompress);
338 int WriteOtherGeometry(FILE *fm) const;
339
340 template <typename dataType, typename triangulationType>
341 int PerformSimplification(
342 const std::vector<std::tuple<int, double, int>> &constraints,
343 int nbConstraints,
344 int vertexNumber,
345 double *array,
346 const triangulationType &triangulation);
347
348 // Numeric management.
349
350 template <typename type>
351 static type abs(const type var) {
352 return (var >= 0) ? var : -var;
353 }
354
355 template <typename type>
356 static type abs_diff(const type var1, const type var2) {
357 return (var1 > var2) ? var1 - var2 : var2 - var1;
358 }
359
360 protected:
361 // General.
364
365 // Parameters
367 bool ZFPOnly{false};
368 double ZFPTolerance{50};
370
371 double Tolerance{10};
372 double MaximumError{10};
374 std::string SQMethod{};
375 bool Subdivide{false};
377
380 double dataSpacing_[3];
381 double dataOrigin_[3];
382
383 std::vector<char> dataArrayName_{};
384
385 // Persistence compression.
386 std::vector<int> segmentation_{};
387 std::vector<std::tuple<double, int>> mapping_{};
388 std::vector<std::tuple<int, double, int>> criticalConstraints_{};
389
390 // IO.
394 std::vector<double> decompressedData_{};
395 std::vector<SimplexId> decompressedOffsets_{};
396 std::vector<SimplexId> compressedOffsets_{};
398 char *fileName{};
399
400 // Char array that identifies the file format.
401 const char *magicBytes_{"TTKCompressedFileFormat"};
402 // Current version of the file format. To be incremented at every
403 // breaking change to keep backward compatibility.
404 const unsigned long formatVersion_{2};
405 };
406
407} // namespace ttk
408
409#include <OtherCompression.h>
411
412template <class dataType, typename triangulationType>
414 const dataType *const inputData,
415 const SimplexId *const inputOffsets,
416 dataType *outputData,
417 const triangulationType &triangulation) {
418 this->printMsg("Starting compression...");
419
420// check the consistency of the variables -- to adapt
421#ifndef TTK_ENABLE_KAMIKAZE
422 if(inputData == nullptr)
423 return -2;
424 if(outputData == nullptr)
425 return -3;
426// if (tol < 0 || tol > 100) return -4;
427#endif
428
429 int const vertexNumber = triangulation.getNumberOfVertices();
430
431 int const res = 0;
433 compressForPersistenceDiagram(vertexNumber, inputData, inputOffsets,
434 outputData, Tolerance, triangulation);
437 vertexNumber, inputData, inputOffsets, outputData, Tolerance);
438
439 return res;
440}
441
442template <typename triangulationType>
444 FILE *fp, const triangulationType &triangulation) {
445 // [fp->] Read headers.
446
447 this->printMsg("Successfully read metadata.");
448
449 if(ZFPOnly && ZFPTolerance < 0.0) {
450 this->printMsg("Wrong ZFP absolute error tolerance for ZFP-only use.");
451 return -4;
452 }
453
454 bool const useZlib = Read<uint8_t>(fp);
455 unsigned char *dest;
456 std::vector<unsigned char> ddest;
457 unsigned long destLen;
458
459#ifdef TTK_ENABLE_ZLIB
460 if(useZlib) {
461 // [fp->ff] Read compressed data.
462 auto sl = Read<uint64_t>(fp); // Compressed size...
463 auto dl = Read<uint64_t>(fp); // Uncompressed size...
464
465 destLen = dl;
466 std::vector<unsigned char> ssource(sl);
467 ReadByteArray(fp, ssource.data(), sl);
468 this->printMsg("Successfully read compressed data.");
469
470 // [ff->fm] Decompress data.
471 ddest.resize(destLen);
472 dest = ddest.data();
473 CompressWithZlib(true, dest, destLen, ssource.data(), sl);
474 this->printMsg("Successfully uncompressed data.");
475
476 } else {
477 this->printMsg("File was not compressed with ZLIB.");
478
479 Read<uint64_t>(fp); // Compressed size...
480 const auto dl = Read<uint64_t>(fp); // Uncompressed size...
481
482 destLen = dl;
483 ddest.resize(destLen);
484 dest = ddest.data();
485 ReadByteArray(fp, dest, destLen);
486 }
487#else
488 if(useZlib) {
489 this->printMsg(" File compressed but ZLIB not installed! Aborting.");
490 return -4;
491
492 } else {
493 this->printMsg(" ZLIB not installed, but file was not compressed anyways.");
494
495 Read<uint64_t>(fp); // Compressed size...
496 const auto dl = Read<uint64_t>(fp); // Uncompressed size...
497
498 destLen = dl;
499 ddest.resize(destLen);
500 dest = ddest.data();
501 ReadByteArray(fp, dest, destLen);
502 }
503#endif
504
505 // [fm->] Read data.
506 char *buf = reinterpret_cast<char *>(dest);
507
508 //#ifndef _MSC_VER
509 // FILE *fm = fmemopen(buf, destLen, "r+");
510 //#else
511 const std::string s = fileName + std::string(".temp");
512 const char *ffn = s.c_str();
513 FILE *ftemp = fopen(ffn, "wb");
514 fwrite(buf, destLen, sizeof(char), ftemp);
515 fclose(ftemp);
516 FILE *fm = fopen(ffn, "rb");
517 //#endif
518
519 // Do read topology.
520 if(!(ZFPOnly)) {
521 if(compressionType_ == (int)ttk::CompressionType::PersistenceDiagram)
522 ReadPersistenceTopology(fm);
523 else if(compressionType_ == (int)ttk::CompressionType::Other)
524 ReadOtherTopology(fm);
525 }
526
527 this->printMsg("Successfully read topology.");
528
529 // Get altered geometry.
530 // Rebuild topologically consistent geometry.
531 int status = 0;
532 if(compressionType_ == (int)ttk::CompressionType::PersistenceDiagram)
533 status = ReadPersistenceGeometry(fm, triangulation);
534 else if(compressionType_ == (int)ttk::CompressionType::Other)
535 status = ReadOtherGeometry(fm);
536
537 fclose(fm);
538 // #ifdef _MSC_VER
539 remove(ffn);
540 // #endif
541 fclose(fp);
542
543 if(status == 0) {
544 this->printMsg("Successfully read geometry.");
545 this->printMsg("Successfully read file.");
546 } else {
547 this->printMsg("Failed to write (possibly ZFP)!");
548 this->printMsg("File may be corrupted!");
549 }
550
551 return status;
552}
AbstractTriangulation is an interface class that defines an interface for efficient traversal methods...
Minimalist debugging class.
Definition Debug.h:88
int printErr(const std::string &msg, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cerr) const
Definition Debug.h:149
TTK processing package for the topological simplification of scalar data.
int preconditionTriangulation(AbstractTriangulation *triangulation)
TTK processing package for the computation of persistence diagrams.
TTK topologicalCompression processing package.
static bool cmp2(const std::tuple< double, int > &a, const std::tuple< double, int > &b)
int ReadFromFile(FILE *fm, const triangulationType &triangulation)
int ReadPersistenceIndex(FILE *fm, std::vector< std::tuple< double, int > > &mappings, std::vector< std::tuple< double, int > > &mappingsSortedPerValue, std::vector< std::tuple< int, double, int > > &constraints, double &min, double &max, int &nbConstraints) const
void Write(FILE *fm, T data) const
int WriteCompactSegmentation(FILE *fm, const std::vector< int > &segmentation, int numberOfVertices, int numberOfSegments) const
std::vector< int > & getSegmentation()
std::vector< SimplexId > compressedOffsets_
std::vector< std::tuple< double, int > > & getMapping()
void WriteByteArray(FILE *fm, const T *buffer, size_t length) const
int compressForOther(int vertexNumber, const dataType *const inputData, const SimplexId *const inputOffsets, dataType *outputData, const double &tol) const
void relToAbsZFPTolerance(const double zfpRelTol, const std::array< double, 2 > &sfRange)
Switch from a relative (% of scalar field range) to an absolute ZFP tolerance.
int WritePersistenceIndex(FILE *fm, std::vector< std::tuple< double, int > > &mapping, std::vector< std::tuple< int, double, int > > &constraints) const
int execute(const dataType *const inputData, const SimplexId *const inputOffsets, dataType *outputData, const triangulationType &triangulation)
void setUseTopologicalSimplification(bool useTopologicalSimplification)
int compressForPersistenceDiagram(int vertexNumber, const dataType *const inputData, const SimplexId *const inputOffset, dataType *outputData, const double &tol, const triangulationType &triangulation)
int computePersistencePairs(std::vector< std::tuple< SimplexId, SimplexId, dataType > > &JTPairs, std::vector< std::tuple< SimplexId, SimplexId, dataType > > &STPairs, const dataType *const inputScalars_, const SimplexId *const inputOffsets, const triangulationType &triangulation)
void ReadByteArray(FILE *fm, T *buffer, size_t length) const
int WriteToFile(FILE *fp, int compressionType, bool zfpOnly, const char *sqMethod, int dataType, int *dataExtent, double *dataSpacing, double *dataOrigin, double *data, double tolerance, double zfpTolerance, const std::string &dataArrayName)
std::vector< std::tuple< int, double, int > > & getCriticalConstraints()
void setTolerance(const double data)
const std::vector< char > & getDataArrayName() const
static bool cmp(const std::tuple< double, int > &a, const std::tuple< double, int > &b)
std::vector< double > decompressedData_
std::vector< std::tuple< int, double, int > > criticalConstraints_
std::vector< SimplexId > & getDecompressedOffsets()
void CropIntervals(std::vector< std::tuple< dataType, int > > &mappings, std::vector< std::tuple< dataType, int > > &mappingsSortedPerValue, double min, double max, int vertexNumber, double *array, std::vector< int > &Seg) const
void setCompressionType(int compressionType)
int ReadCompactSegmentation(FILE *fm, std::vector< int > &segmentation, int &numberOfVertices, int &numberOfSegments) const
std::vector< std::tuple< double, int > > mapping_
void preconditionTriangulation(AbstractTriangulation *const triangulation)
LegacyTopologicalSimplification topologicalSimplification
int WriteMetaData(FILE *fp, int compressionType, bool zfpOnly, const char *sqMethod, int dataType, int *dataExtent, double *dataSpacing, double *dataOrigin, double tolerance, double zfpTolerance, const std::string &dataArrayName)
std::vector< SimplexId > & getCompressedOffsets()
static unsigned int log2(int val)
void setSQ(const std::string &sqMethod)
void setMaximumError(double maximumError)
std::vector< SimplexId > decompressedOffsets_
std::vector< double > & getDecompressedData()
void preconditionTriangulation(AbstractTriangulation *tri, const bool preproc=true)
Definition FTMTree_CT.h:72
The Topology ToolKit.
int SimplexId
Identifier type for simplices of any dimension.
Definition DataTypes.h:22
printMsg(debug::output::BOLD+" | | | | | . \\ | | (__| | / __/| |_| / __/|__ _|"+debug::output::ENDCOLOR, debug::Priority::PERFORMANCE, debug::LineMode::NEW, stream)