diff --git a/.clang-format b/.clang-format new file mode 100644 index 0000000..293b66b --- /dev/null +++ b/.clang-format @@ -0,0 +1,134 @@ +# SPDX-FileCopyrightText: 2021 SeisSol Group +# +# SPDX-License-Identifier: BSD-3-Clause +# SPDX-LicenseComments: Full text under /LICENSE and /LICENSES/ +# +# SPDX-FileContributor: Author lists in /AUTHORS and /CITATION.cff + +Language: Cpp +AccessModifierOffset: -4 +AlignAfterOpenBracket: Align +AlignConsecutiveMacros: false +AlignConsecutiveAssignments: false +AlignConsecutiveDeclarations: false +AlignEscapedNewlines: Right +AlignOperands: true +AlignTrailingComments: true +AllowAllArgumentsOnNextLine: true +AllowAllConstructorInitializersOnNextLine: true +AllowAllParametersOfDeclarationOnNextLine: true +AllowShortBlocksOnASingleLine: false +AllowShortCaseLabelsOnASingleLine: false +AllowShortFunctionsOnASingleLine: All +AllowShortLambdasOnASingleLine: All +AllowShortIfStatementsOnASingleLine: Never +AllowShortLoopsOnASingleLine: false +AlwaysBreakAfterDefinitionReturnType: None +AlwaysBreakAfterReturnType: None +AlwaysBreakBeforeMultilineStrings: false +AlwaysBreakTemplateDeclarations: Yes +BinPackArguments: false +BinPackParameters: false +BraceWrapping: + AfterCaseLabel: false + AfterClass: false + AfterControlStatement: false + AfterEnum: false + AfterFunction: false + AfterNamespace: false + AfterObjCDeclaration: false + AfterStruct: false + AfterUnion: false + AfterExternBlock: false + BeforeCatch: false + BeforeElse: false + IndentBraces: false + SplitEmptyFunction: true + SplitEmptyRecord: true + SplitEmptyNamespace: true +BreakBeforeBinaryOperators: None +BreakBeforeBraces: Attach +BreakBeforeInheritanceComma: false +BreakInheritanceList: BeforeColon +BreakBeforeTernaryOperators: true +BreakConstructorInitializersBeforeComma: false +BreakConstructorInitializers: BeforeColon +BreakAfterJavaFieldAnnotations: false +BreakStringLiterals: true +ColumnLimit: 100 +CommentPragmas: '^ IWYU pragma:' +CompactNamespaces: false +ConstructorInitializerAllOnOneLineOrOnePerLine: false +ConstructorInitializerIndentWidth: 4 +ContinuationIndentWidth: 4 +Cpp11BracedListStyle: true +DerivePointerAlignment: false +DisableFormat: false +ExperimentalAutoDetectBinPacking: false +FixNamespaceComments: true +ForEachMacros: + - foreach + - Q_FOREACH + - BOOST_FOREACH +IncludeBlocks: Regroup +IncludeCategories: + # keep the doctest headers in front + - Regex: '^(<|")doctest' + Priority: 1 + - Regex: '^"(llvm|llvm-c|clang|clang-c)/' + Priority: 3 + - Regex: '^(<|"(gtest|gmock|isl|json)/)' + Priority: 4 + - Regex: '.*' + Priority: 2 +IncludeIsMainRegex: '(Test)?$' +IndentCaseLabels: false +IndentPPDirectives: None +IndentWidth: 2 +IndentWrappedFunctionNames: true +JavaScriptQuotes: Leave +JavaScriptWrapImports: true +KeepEmptyLinesAtTheStartOfBlocks: true +MacroBlockBegin: '' +MacroBlockEnd: '' +MaxEmptyLinesToKeep: 1 +NamespaceIndentation: None +ObjCBinPackProtocolList: Auto +ObjCBlockIndentWidth: 2 +ObjCSpaceAfterProperty: false +ObjCSpaceBeforeProtocolList: true +PenaltyBreakAssignment: 2 +PenaltyBreakBeforeFirstCallParameter: 19 +PenaltyBreakComment: 300 +PenaltyBreakFirstLessLess: 120 +PenaltyBreakString: 1000 +PenaltyBreakTemplateDeclaration: 10 +PenaltyExcessCharacter: 1000000 +PenaltyReturnTypeOnItsOwnLine: 60 +PointerAlignment: Left +QualifierAlignment: Left +ReflowComments: true +SortIncludes: true +SortUsingDeclarations: true +SpaceAfterCStyleCast: false +SpaceAfterLogicalNot: false +SpaceAfterTemplateKeyword: true +SpaceBeforeAssignmentOperators: true +SpaceBeforeCpp11BracedList: false +SpaceBeforeCtorInitializerColon: true +SpaceBeforeInheritanceColon: true +SpaceBeforeParens: ControlStatements +SpaceBeforeRangeBasedForLoopColon: true +SpaceInEmptyParentheses: false +SpacesBeforeTrailingComments: 1 +SpacesInAngles: false +SpacesInContainerLiterals: true +SpacesInCStyleCastParentheses: false +SpacesInParentheses: false +SpacesInSquareBrackets: false +Standard: Cpp11 +StatementMacros: + - Q_UNUSED + - QT_REQUIRE_VERSION +TabWidth: 8 +UseTab: Never diff --git a/.git-blame-ignore-revs b/.git-blame-ignore-revs new file mode 100644 index 0000000..ed592c1 --- /dev/null +++ b/.git-blame-ignore-revs @@ -0,0 +1,9 @@ +# SPDX-FileCopyrightText: 2026 SeisSol Group +# +# SPDX-License-Identifier: BSD-3-Clause +# SPDX-LicenseComments: Full text under /LICENSE and /LICENSES/ +# +# SPDX-FileContributor: Author lists in /AUTHORS and /CITATION.cff + +# apply clang-format and pre-commit +a557b85c5e4b925fc739aa1ded1a76ccde4ee6fa diff --git a/.github/workflows/pre-commit.yml b/.github/workflows/pre-commit.yml new file mode 100644 index 0000000..6f56a26 --- /dev/null +++ b/.github/workflows/pre-commit.yml @@ -0,0 +1,21 @@ +# SPDX-FileCopyrightText: 2025 SeisSol Group +# +# SPDX-License-Identifier: BSD-3-Clause +# SPDX-LicenseComments: Full text under /LICENSE and /LICENSES/ +# +# SPDX-FileContributor: Author lists in /AUTHORS and /CITATION.cff + +name: pre-commit +on: + - push + +jobs: + pre-commit: + name: pre-commit + runs-on: ubuntu-24.04 + steps: + - uses: actions/checkout@v6 + with: + submodules: recursive + + - uses: pre-commit/action@v3.0.1 diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml new file mode 100644 index 0000000..fc7c532 --- /dev/null +++ b/.pre-commit-config.yaml @@ -0,0 +1,92 @@ +# SPDX-FileCopyrightText: 2025 SeisSol Group +# +# SPDX-License-Identifier: BSD-3-Clause +# SPDX-LicenseComments: Full text under /LICENSE and /LICENSES/ +# +# SPDX-FileContributor: Author lists in /AUTHORS and /CITATION.cff + +--- + +repos: +- repo: https://github.com/pre-commit/pre-commit-hooks + rev: v6.0.0 + hooks: + - id: check-merge-conflict + name: '[GENERIC] merge conflict check' + - id: check-symlinks + name: '[GENERIC] symlink check' + - id: destroyed-symlinks + name: '[GENERIC] detect broken symlinks' + - id: detect-private-key + name: '[GENERIC] detect private keys uploaded by accident' + - id: check-case-conflict + name: '[GENERIC] detect OS file naming case conflicts' + - id: check-executables-have-shebangs + name: '[GENERIC] check for shebangs in executable files' + - id: check-illegal-windows-names + name: '[GENERIC] detect illegal Windows file names' + - id: check-json + name: '[JSON] check' + - id: check-xml + name: '[XML] check' + - id: check-shebang-scripts-are-executable + name: '[GENERIC] check that shebang-containing files are executable' + + +- repo: https://github.com/DavidAnson/markdownlint-cli2 + rev: v0.18.1 + hooks: + - id: markdownlint-cli2 + name: '[MARKDOWN] lint' + +#- repo: https://github.com/fsfe/reuse-tool +# rev: v6.0.0 +# hooks: +# - id: reuse +# name: '[GENERIC] REUSE compatibiltiy' + +#- repo: https://github.com/psf/black-pre-commit-mirror +# rev: 25.1.0 +# hooks: +# - id: black +# files: ^(?!preprocessing|postprocessing) +# name: '[PYTHON] black' +#- repo: https://github.com/pycqa/isort +# rev: 6.0.1 +# hooks: +# - id: isort +# files: ^(?!preprocessing|postprocessing) +# args: ["--profile", "black"] +# name: '[PYTHON] isort' +- repo: https://github.com/pycqa/bandit + rev: 1.8.6 + hooks: + - id: bandit + args: ["--confidence-level", "high", "--severity-level", "high"] + name: '[PYTHON] bandit' +#- repo: https://github.com/pycqa/flake8 +# rev: '7.3.0' +# hooks: +# - id: flake8 +# files: ^(?!preprocessing|postprocessing) +# name: '[PYTHON] Flake8' + +- repo: https://github.com/sphinx-contrib/sphinx-lint + rev: 'v1.0.0' + hooks: + - id: sphinx-lint + name: '[SPHINX/RST] sphinx lint' + +- repo: https://github.com/pre-commit/mirrors-clang-format + rev: 'v22.1.0' + hooks: + - id: clang-format + name: '[C++] clang-format' + +- repo: https://github.com/pre-commit/pre-commit-hooks + rev: v6.0.0 + hooks: + - id: end-of-file-fixer + name: '[GENERIC] newline eof' + - id: trailing-whitespace + name: '[GENERIC] remove trailing whitespace' diff --git a/README.md b/README.md index c92b2d7..735a2bb 100644 --- a/README.md +++ b/README.md @@ -1 +1,3 @@ -A collection of useful tools for visualization of [SeisSol](http://www.seissol.org)'s results. \ No newline at end of file +# Postprocessing Tools for SeisSol + +Useful tools for visualizing [SeisSol](http://www.seissol.org) results. diff --git a/compress_xmdf/SConstruct b/compress_xmdf/SConstruct index ac62c2c..f845ed9 100644 --- a/compress_xmdf/SConstruct +++ b/compress_xmdf/SConstruct @@ -1,4 +1,3 @@ -#!/usr/bin/python ## # @file # This file is part of SeisSol. @@ -76,4 +75,4 @@ Export('env') SConscript('src/SConscript', variant_dir='#/build/'+buildDir, src_dir='#/', duplicate=0) Import('env') -env.Program('#/build/compress-xdmf', env.sourceFiles) \ No newline at end of file +env.Program('#/build/compress-xdmf', env.sourceFiles) diff --git a/compress_xmdf/src/SConscript b/compress_xmdf/src/SConscript index 0a1a958..dce56aa 100644 --- a/compress_xmdf/src/SConscript +++ b/compress_xmdf/src/SConscript @@ -1,4 +1,3 @@ -#!/usr/bin/python ## # @file # This file is part of SeisSol. @@ -44,4 +43,4 @@ sourceFiles = ['main.cpp'] for f in sourceFiles: env.sourceFiles.append(env.Object(f)) -Export('env') \ No newline at end of file +Export('env') diff --git a/compress_xmdf/src/hdf5_helper.h b/compress_xmdf/src/hdf5_helper.h index d64db9a..a879897 100644 --- a/compress_xmdf/src/hdf5_helper.h +++ b/compress_xmdf/src/hdf5_helper.h @@ -2,7 +2,8 @@ * @file * This file is part of SeisSol. * - * @author Sebastian Rettenberger (sebastian.rettenberger AT tum.de, http://www5.in.tum.de/wiki/index.php/Sebastian_Rettenberger) + * @author Sebastian Rettenberger (sebastian.rettenberger AT tum.de, + * http://www5.in.tum.de/wiki/index.php/Sebastian_Rettenberger) * * @section LICENSE * Copyright (c) 2017, SeisSol Group @@ -40,21 +41,16 @@ #ifndef HDF5_HELPER_H #define HDF5_HELPER_H -#include - #include "utils/logger.h" -template -static void checkH5Err(T status) -{ - if (status < 0) - logError() << "An HDF5 error occurred"; -} +#include -inline -static hsize_t sizeofType(hid_t type) -{ - return H5Tget_size(type); +template +static void checkH5Err(T status) { + if (status < 0) + logError() << "An HDF5 error occurred"; } -#endif // HDF5_HELPER_H \ No newline at end of file +inline static hsize_t sizeofType(hid_t type) { return H5Tget_size(type); } + +#endif // HDF5_HELPER_H diff --git a/compress_xmdf/src/input.h b/compress_xmdf/src/input.h index 969b26a..fa066f3 100644 --- a/compress_xmdf/src/input.h +++ b/compress_xmdf/src/input.h @@ -2,7 +2,8 @@ * @file * This file is part of SeisSol. * - * @author Sebastian Rettenberger (sebastian.rettenberger AT tum.de, http://www5.in.tum.de/wiki/index.php/Sebastian_Rettenberger) + * @author Sebastian Rettenberger (sebastian.rettenberger AT tum.de, + * http://www5.in.tum.de/wiki/index.php/Sebastian_Rettenberger) * * @section LICENSE * Copyright (c) 2017, SeisSol Group @@ -40,398 +41,360 @@ #ifndef INPUT_H #define INPUT_H +#include "hdf5_helper.h" +#include "output.h" +#include "utils/stringutils.h" + #include #include +#include #include +#include #include #include -#include -#include +struct Variable { + std::string name; + unsigned int timesteps; +}; -#include "utils/stringutils.h" +class Input { + protected: + size_t m_numElements; + size_t m_numVertices; -#include "hdf5_helper.h" -#include "output.h" + unsigned int m_verticesPerElement; -struct Variable -{ - std::string name; - unsigned int timesteps; -}; + char* m_buffer; -class Input -{ -protected: - size_t m_numElements; - size_t m_numVertices; - - unsigned int m_verticesPerElement; - - char* m_buffer; - -public: - Input() - : m_numElements(0), m_numVertices(0), m_verticesPerElement(0) - { - m_buffer = new char[OutputVar::CHUNK_SIZE]; - } - - virtual ~Input() - { - delete [] m_buffer; - } - - size_t numElements() const - { - return m_numElements; - } - - size_t numVertices() const - { - return m_numVertices; - } - - unsigned int verticesPerElement() const - { - return m_verticesPerElement; - } - - virtual std::vector getVarList() = 0; - - virtual void writeVariable(const Variable &variable, hid_t nativeType, bool isVertex, OutputVar &writer) = 0; + public: + Input() : m_numElements(0), m_numVertices(0), m_verticesPerElement(0) { + m_buffer = new char[OutputVar::CHUNK_SIZE]; + } + + virtual ~Input() { delete[] m_buffer; } + + size_t numElements() const { return m_numElements; } + + size_t numVertices() const { return m_numVertices; } + + unsigned int verticesPerElement() const { return m_verticesPerElement; } + + virtual std::vector getVarList() = 0; + + virtual void writeVariable(const Variable& variable, + hid_t nativeType, + bool isVertex, + OutputVar& writer) = 0; }; +class HDF5Input : public Input { + private: + hid_t m_file; + + public: + HDF5Input(const std::string& xdmfFile) { + std::string h5File(xdmfFile); + utils::StringUtils::replaceLast(h5File, ".xdmf", ".h5"); + + logInfo() << "Reading heavy data from" << h5File; + + m_file = H5Fopen(h5File.c_str(), H5F_ACC_RDONLY, H5P_DEFAULT); + checkH5Err(m_file); + } + + virtual ~HDF5Input() { checkH5Err(H5Fclose(m_file)); } + + std::vector getVarList() { + std::vector variables; + + hsize_t idx = 0; + checkH5Err(H5Literate_by_name( + m_file, "/", H5_INDEX_NAME, H5_ITER_NATIVE, &idx, extractVars, &variables, H5P_DEFAULT)); + + for (std::vector::iterator it = variables.begin(); it != variables.end(); ++it) { + + hid_t var = H5Dopen(m_file, ("/" + it->name).c_str(), H5P_DEFAULT); + checkH5Err(var); + + hid_t space = H5Dget_space(var); + checkH5Err(space); + + hsize_t size[2]; + int ndims = H5Sget_simple_extent_dims(space, size, 0L); + checkH5Err(ndims); + + if (it->name == "connect") { + m_numElements = size[0]; + m_verticesPerElement = size[1]; + } else if (it->name == "geometry") { + m_numVertices = size[0]; + } else { + if (ndims > 1) + it->timesteps = size[0]; + else + it->timesteps = 0; + } + + checkH5Err(H5Sclose(space)); + checkH5Err(H5Dclose(var)); + } + + return variables; + } + + void writeVariable(const Variable& variable, hid_t nativeType, bool isVertex, OutputVar& writer) { + unsigned int timesteps = std::max(variable.timesteps, 1u); -class HDF5Input : public Input -{ -private: - hid_t m_file; - -public: - HDF5Input(const std::string &xdmfFile) - { - std::string h5File(xdmfFile); - utils::StringUtils::replaceLast(h5File, ".xdmf", ".h5"); - - logInfo() << "Reading heavy data from" << h5File; - - m_file = H5Fopen(h5File.c_str(), H5F_ACC_RDONLY, H5P_DEFAULT); - checkH5Err(m_file); - } - - virtual ~HDF5Input() - { - checkH5Err(H5Fclose(m_file)); - } - - std::vector getVarList() - { - std::vector variables; - - hsize_t idx = 0; - checkH5Err(H5Literate_by_name(m_file, "/", H5_INDEX_NAME, H5_ITER_NATIVE, &idx, extractVars, &variables, H5P_DEFAULT)); - - for (std::vector::iterator it = variables.begin(); - it != variables.end(); ++it) { - - hid_t var = H5Dopen(m_file, ("/"+it->name).c_str(), H5P_DEFAULT); - checkH5Err(var); - - hid_t space = H5Dget_space(var); - checkH5Err(space); - - hsize_t size[2]; - int ndims = H5Sget_simple_extent_dims(space, size, 0L); - checkH5Err(ndims); - - if (it->name == "connect") { - m_numElements = size[0]; - m_verticesPerElement = size[1]; - } else if (it->name == "geometry") { - m_numVertices = size[0]; - } else { - if (ndims > 1) - it->timesteps = size[0]; - else - it->timesteps = 0; - } - - checkH5Err(H5Sclose(space)); - checkH5Err(H5Dclose(var)); - } - - return variables; - } - - void writeVariable(const Variable &variable, hid_t nativeType, bool isVertex, OutputVar &writer) - { - unsigned int timesteps = std::max(variable.timesteps, 1u); - - hid_t var = H5Dopen(m_file, variable.name.c_str(), H5P_DEFAULT); - checkH5Err(var); - - hid_t space = H5Dget_space(var); - checkH5Err(space); - - hsize_t extent[2]; - int ndims = H5Sget_simple_extent_dims(space, extent, 0L); - checkH5Err(ndims); - if (ndims > 2) - logError() << "Dimension > 2 are not supported"; - - unsigned int dim2 = 1; - hsize_t nElements; - - if (variable.timesteps == 0) { - nElements = extent[0]; - if (ndims > 1) - dim2 = extent[1]; - } else { - nElements = extent[1]; - } - - hsize_t chunkSize = OutputVar::CHUNK_SIZE / dim2 / sizeofType(nativeType); - - // Write data - for (unsigned int t = 0; t < timesteps; t++) { - unsigned long pos = 0; - while (pos < nElements) { - const unsigned long left = nElements - pos; - hsize_t tmpChunkSize = chunkSize; - - if (left < tmpChunkSize) - tmpChunkSize = left; - - hsize_t offset[2]; - hsize_t size[2]; - if (variable.timesteps > 0) { - offset[0] = t; - offset[1] = pos; - size[0] = 1; - size[1] = tmpChunkSize; - } else { - offset[0] = pos; - offset[1] = 0; - size[0] = tmpChunkSize; - size[1] = dim2; - } - - hid_t memspace = H5Screate_simple(ndims, size, 0L); - checkH5Err(memspace); - - checkH5Err(H5Sselect_hyperslab(space, H5S_SELECT_SET, offset, 0L, size, 0L)); - checkH5Err(H5Dread(var, nativeType, memspace, space, H5P_DEFAULT, m_buffer)); - - checkH5Err(H5Sclose(memspace)); - - writer.write(m_buffer, nativeType, offset, size); - - pos += chunkSize; - } - } - - } - -private: - static herr_t extractVars(hid_t group, const char* name, const H5L_info_t* info, void* op_data) - { - std::vector *variables = static_cast*>(op_data); - - Variable var; - var.name = name; - var.timesteps = false; - - variables->push_back(var); - - return 0; - } + hid_t var = H5Dopen(m_file, variable.name.c_str(), H5P_DEFAULT); + checkH5Err(var); + + hid_t space = H5Dget_space(var); + checkH5Err(space); + + hsize_t extent[2]; + int ndims = H5Sget_simple_extent_dims(space, extent, 0L); + checkH5Err(ndims); + if (ndims > 2) + logError() << "Dimension > 2 are not supported"; + + unsigned int dim2 = 1; + hsize_t nElements; + + if (variable.timesteps == 0) { + nElements = extent[0]; + if (ndims > 1) + dim2 = extent[1]; + } else { + nElements = extent[1]; + } + + hsize_t chunkSize = OutputVar::CHUNK_SIZE / dim2 / sizeofType(nativeType); + + // Write data + for (unsigned int t = 0; t < timesteps; t++) { + unsigned long pos = 0; + while (pos < nElements) { + const unsigned long left = nElements - pos; + hsize_t tmpChunkSize = chunkSize; + + if (left < tmpChunkSize) + tmpChunkSize = left; + + hsize_t offset[2]; + hsize_t size[2]; + if (variable.timesteps > 0) { + offset[0] = t; + offset[1] = pos; + size[0] = 1; + size[1] = tmpChunkSize; + } else { + offset[0] = pos; + offset[1] = 0; + size[0] = tmpChunkSize; + size[1] = dim2; + } + + hid_t memspace = H5Screate_simple(ndims, size, 0L); + checkH5Err(memspace); + + checkH5Err(H5Sselect_hyperslab(space, H5S_SELECT_SET, offset, 0L, size, 0L)); + checkH5Err(H5Dread(var, nativeType, memspace, space, H5P_DEFAULT, m_buffer)); + + checkH5Err(H5Sclose(memspace)); + + writer.write(m_buffer, nativeType, offset, size); + + pos += chunkSize; + } + } + } + + private: + static herr_t extractVars(hid_t group, const char* name, const H5L_info_t* info, void* op_data) { + std::vector* variables = static_cast*>(op_data); + + Variable var; + var.name = name; + var.timesteps = false; + + variables->push_back(var); + + return 0; + } }; -class BinaryInput : public Input -{ -private: - std::string m_fileBase; - std::string m_fileRegex; - - char* m_readBuffer; - -public: - BinaryInput(const std::string &xdmfFile) - : m_fileBase(xdmfFile) - { - utils::StringUtils::replaceLast(m_fileBase, ".xdmf", "_"); - m_fileRegex = m_fileBase + "*.bin"; - - m_readBuffer = new char[2*OutputVar::CHUNK_SIZE]; - - logInfo() << "Reading heavy data from" << m_fileRegex.c_str(); - } - - virtual ~BinaryInput() - { - delete [] m_readBuffer; - } - - std::vector getVarList() - { - std::vector variables; - - glob_t globResult; - glob(m_fileRegex.c_str(), 0, 0L, &globResult); - - for (size_t i = 0; i < globResult.gl_pathc; i++) { - std::string file(globResult.gl_pathv[i]); - - Variable var; - var.name = file.substr(m_fileRegex.size()-5, file.size()-m_fileRegex.size()+1); - - if (var.name == "geometry") { - int fd = openByVar(var.name); - m_numVertices = getFileSize(fd) / (3 * sizeof(double)); - close(fd); - } else if (var.name == "partition") { - int fd = openByVar(var.name); - m_numElements = getFileSize(fd) / sizeof(int); - close(fd); - } - - var.timesteps = 0; - - variables.push_back(var); - } - - globfree(&globResult); - - // Get the timesteps and the vertices per element - for (std::vector::iterator it = variables.begin(); - it != variables.end(); ++it) { - - if (it->name == "connect") { - int fd = openByVar(it->name); - m_verticesPerElement = getFileSize(fd) / (m_numElements * sizeof(unsigned long)); - close(fd); - } else if (it->name == "geometry" || it->name == "partition") { - // Do nothing - } else { - int fd = openByVar(it->name); - it->timesteps = getFileSize(fd) / (m_numElements * sizeof(double)); - close(fd); - - if (it->timesteps == 1) - it->timesteps = 0; - } - } - - return variables; - } - - void writeVariable(const Variable &variable, hid_t nativeType, bool isVertex, OutputVar &writer) - { - - int fd = openByVar(variable.name); - - size_t typeSize = sizeofType(nativeType); - if (doCompression(nativeType)) - typeSize *= 2; - - size_t fileSize = getFileSize(fd); - lseek(fd, 0, SEEK_SET); - - unsigned int dim2 = 1; - hsize_t nElements = isVertex ? m_numVertices : m_numElements; - - if (variable.timesteps == 0) { - dim2 = fileSize / (nElements * typeSize); - } - - unsigned int timesteps = std::max(variable.timesteps, 1u); - - hsize_t chunkSize = OutputVar::CHUNK_SIZE / dim2 / sizeofType(nativeType); - - // Write data - for (unsigned int t = 0; t < timesteps; t++) { - unsigned long pos = 0; - while (pos < nElements) { - const unsigned long left = nElements - pos; - hsize_t tmpChunkSize = chunkSize; - - if (left < tmpChunkSize) - tmpChunkSize = left; - - hsize_t offset[2]; - hsize_t size[2]; - if (variable.timesteps > 0) { - offset[0] = t; - offset[1] = pos; - size[0] = 1; - size[1] = tmpChunkSize; - } else { - offset[0] = pos; - offset[1] = 0; - size[0] = tmpChunkSize; - size[1] = dim2; - } - - size_t readSize = tmpChunkSize * dim2 * sizeofType(nativeType); - if (doCompression(nativeType)) - readSize *= 2; - - if (read(fd, m_readBuffer, readSize) != readSize) - logError() << "Could not read data" << readSize; - - void* buffer; - if (doCompression(nativeType)) { - compress(); - buffer = m_buffer; - } else { - buffer = m_readBuffer; - } - - writer.write(buffer, nativeType, offset, size); - - pos += chunkSize; - } - } - - if (variable.timesteps > 0) { - } else { - } - - close(fd); - } - -private: - int openByVar(const std::string &var) - { - int fd = open((m_fileBase+var+".bin").c_str(), 0); - if (fd < 0) - logError() << "Could not open file" << (m_fileBase+var+".bin").c_str(); - - return fd; - } - - void compress() - { - for (unsigned int i = 0; i < OutputVar::CHUNK_SIZE / sizeof(float); i++) { - reinterpret_cast(m_buffer)[i] = reinterpret_cast(m_readBuffer)[i]; - } - } - -private: - static size_t getFileSize(int fd) - { - size_t len = lseek(fd, 0, SEEK_END); - return len; - } - - static bool doCompression(hid_t type) - { - if (type == H5T_NATIVE_FLOAT) - return true; - - return false; - } +class BinaryInput : public Input { + private: + std::string m_fileBase; + std::string m_fileRegex; + + char* m_readBuffer; + + public: + BinaryInput(const std::string& xdmfFile) : m_fileBase(xdmfFile) { + utils::StringUtils::replaceLast(m_fileBase, ".xdmf", "_"); + m_fileRegex = m_fileBase + "*.bin"; + + m_readBuffer = new char[2 * OutputVar::CHUNK_SIZE]; + + logInfo() << "Reading heavy data from" << m_fileRegex.c_str(); + } + + virtual ~BinaryInput() { delete[] m_readBuffer; } + + std::vector getVarList() { + std::vector variables; + + glob_t globResult; + glob(m_fileRegex.c_str(), 0, 0L, &globResult); + + for (size_t i = 0; i < globResult.gl_pathc; i++) { + std::string file(globResult.gl_pathv[i]); + + Variable var; + var.name = file.substr(m_fileRegex.size() - 5, file.size() - m_fileRegex.size() + 1); + + if (var.name == "geometry") { + int fd = openByVar(var.name); + m_numVertices = getFileSize(fd) / (3 * sizeof(double)); + close(fd); + } else if (var.name == "partition") { + int fd = openByVar(var.name); + m_numElements = getFileSize(fd) / sizeof(int); + close(fd); + } + + var.timesteps = 0; + + variables.push_back(var); + } + + globfree(&globResult); + + // Get the timesteps and the vertices per element + for (std::vector::iterator it = variables.begin(); it != variables.end(); ++it) { + + if (it->name == "connect") { + int fd = openByVar(it->name); + m_verticesPerElement = getFileSize(fd) / (m_numElements * sizeof(unsigned long)); + close(fd); + } else if (it->name == "geometry" || it->name == "partition") { + // Do nothing + } else { + int fd = openByVar(it->name); + it->timesteps = getFileSize(fd) / (m_numElements * sizeof(double)); + close(fd); + + if (it->timesteps == 1) + it->timesteps = 0; + } + } + + return variables; + } + + void writeVariable(const Variable& variable, hid_t nativeType, bool isVertex, OutputVar& writer) { + + int fd = openByVar(variable.name); + + size_t typeSize = sizeofType(nativeType); + if (doCompression(nativeType)) + typeSize *= 2; + + size_t fileSize = getFileSize(fd); + lseek(fd, 0, SEEK_SET); + + unsigned int dim2 = 1; + hsize_t nElements = isVertex ? m_numVertices : m_numElements; + + if (variable.timesteps == 0) { + dim2 = fileSize / (nElements * typeSize); + } + + unsigned int timesteps = std::max(variable.timesteps, 1u); + + hsize_t chunkSize = OutputVar::CHUNK_SIZE / dim2 / sizeofType(nativeType); + + // Write data + for (unsigned int t = 0; t < timesteps; t++) { + unsigned long pos = 0; + while (pos < nElements) { + const unsigned long left = nElements - pos; + hsize_t tmpChunkSize = chunkSize; + + if (left < tmpChunkSize) + tmpChunkSize = left; + + hsize_t offset[2]; + hsize_t size[2]; + if (variable.timesteps > 0) { + offset[0] = t; + offset[1] = pos; + size[0] = 1; + size[1] = tmpChunkSize; + } else { + offset[0] = pos; + offset[1] = 0; + size[0] = tmpChunkSize; + size[1] = dim2; + } + + size_t readSize = tmpChunkSize * dim2 * sizeofType(nativeType); + if (doCompression(nativeType)) + readSize *= 2; + + if (read(fd, m_readBuffer, readSize) != readSize) + logError() << "Could not read data" << readSize; + + void* buffer; + if (doCompression(nativeType)) { + compress(); + buffer = m_buffer; + } else { + buffer = m_readBuffer; + } + + writer.write(buffer, nativeType, offset, size); + + pos += chunkSize; + } + } + + if (variable.timesteps > 0) { + } else { + } + + close(fd); + } + + private: + int openByVar(const std::string& var) { + int fd = open((m_fileBase + var + ".bin").c_str(), 0); + if (fd < 0) + logError() << "Could not open file" << (m_fileBase + var + ".bin").c_str(); + + return fd; + } + + void compress() { + for (unsigned int i = 0; i < OutputVar::CHUNK_SIZE / sizeof(float); i++) { + reinterpret_cast(m_buffer)[i] = reinterpret_cast(m_readBuffer)[i]; + } + } + + private: + static size_t getFileSize(int fd) { + size_t len = lseek(fd, 0, SEEK_END); + return len; + } + + static bool doCompression(hid_t type) { + if (type == H5T_NATIVE_FLOAT) + return true; + + return false; + } }; -#endif // INPUT_H \ No newline at end of file +#endif // INPUT_H diff --git a/compress_xmdf/src/main.cpp b/compress_xmdf/src/main.cpp index 76c884c..b9da232 100644 --- a/compress_xmdf/src/main.cpp +++ b/compress_xmdf/src/main.cpp @@ -2,7 +2,8 @@ * @file * This file is part of SeisSol. * - * @author Sebastian Rettenberger (sebastian.rettenberger AT tum.de, http://www5.in.tum.de/wiki/index.php/Sebastian_Rettenberger) + * @author Sebastian Rettenberger (sebastian.rettenberger AT tum.de, + * http://www5.in.tum.de/wiki/index.php/Sebastian_Rettenberger) * * @section LICENSE * Copyright (c) 2016-2017, SeisSol Group @@ -37,250 +38,255 @@ * @section DESCRIPTION */ -#include -#include - -#include - +#include "hdf5_helper.h" +#include "input.h" #include "utils/args.h" #include "utils/logger.h" #include "utils/stringutils.h" -#include "hdf5_helper.h" -#include "input.h" - -static void compressData(unsigned int ndims, hsize_t offset[], hsize_t size[], - hid_t h5invar, hid_t h5inspace, hid_t h5outvar, hid_t h5outspace, - hid_t h5native_type, void* buffer) -{ - hid_t h5memspace = H5Screate_simple(ndims, size, 0L); - checkH5Err(h5memspace); - - checkH5Err(H5Sselect_hyperslab(h5inspace, H5S_SELECT_SET, offset, 0L, size, 0L)); - checkH5Err(H5Dread(h5invar, h5native_type, h5memspace, h5inspace, H5P_DEFAULT, buffer)); - - checkH5Err(H5Sselect_hyperslab(h5outspace, H5S_SELECT_SET, offset, 0L, size, 0L)); - checkH5Err(H5Dwrite(h5outvar, h5native_type, h5memspace, h5outspace, H5P_DEFAULT, buffer)); +#include +#include +#include - checkH5Err(H5Sclose(h5memspace)); +static void compressData(unsigned int ndims, + hsize_t offset[], + hsize_t size[], + hid_t h5invar, + hid_t h5inspace, + hid_t h5outvar, + hid_t h5outspace, + hid_t h5native_type, + void* buffer) { + hid_t h5memspace = H5Screate_simple(ndims, size, 0L); + checkH5Err(h5memspace); + + checkH5Err(H5Sselect_hyperslab(h5inspace, H5S_SELECT_SET, offset, 0L, size, 0L)); + checkH5Err(H5Dread(h5invar, h5native_type, h5memspace, h5inspace, H5P_DEFAULT, buffer)); + + checkH5Err(H5Sselect_hyperslab(h5outspace, H5S_SELECT_SET, offset, 0L, size, 0L)); + checkH5Err(H5Dwrite(h5outvar, h5native_type, h5memspace, h5outspace, H5P_DEFAULT, buffer)); + + checkH5Err(H5Sclose(h5memspace)); } -template -static void compressDataset(hid_t h5ifile, hid_t h5ofile, - const char* varname, hid_t h5native_type, hid_t h5type, - unsigned int compressionLevel, void* buffer, unsigned long bufferSize) -{ - // Open original var - hid_t h5ivar = H5Dopen(h5ifile, varname, H5P_DEFAULT); - checkH5Err(h5ivar); - - // Get dimension information - hid_t h5ispace = H5Dget_space(h5ivar); - checkH5Err(h5ispace); - - int ndims = H5Sget_simple_extent_ndims(h5ispace); - checkH5Err(ndims); - if (ndims > 2) - logError() << "Dimension > 2 are not supported"; - - hsize_t extent[2]; - checkH5Err(H5Sget_simple_extent_dims(h5ispace, extent, 0L)); - hsize_t nelements = extent[0]; - unsigned int dim2 = 1; - if (ndims > 1) - dim2 = extent[1]; - - // Create new dataset - hsize_t chunkSize = bufferSize / dim2 / sizeof(T); - - hsize_t dims[2] = {nelements, dim2}; // Change this for other elements - hid_t h5ospace = H5Screate_simple(ndims, dims, 0L); - checkH5Err(h5ospace); - hid_t h5opcreate = H5Pcreate(H5P_DATASET_CREATE); - checkH5Err(h5opcreate); - hsize_t chunkDims[2] = {std::min(chunkSize, nelements), dim2}; - checkH5Err(H5Pset_chunk(h5opcreate, ndims, chunkDims)); -// checkH5Err(H5Pset_szip(h5opcreate, H5_SZIP_NN_OPTION_MASK, 4)); - checkH5Err(H5Pset_deflate(h5opcreate, compressionLevel)); - hid_t h5ovar = H5Dcreate(h5ofile, varname, h5type, h5ospace, - H5P_DEFAULT, h5opcreate, H5P_DEFAULT); - checkH5Err(h5ovar); - checkH5Err(H5Pclose(h5opcreate)); - - // Transfer data - unsigned long pos = 0; - while (pos < nelements) { - const unsigned long left = nelements - pos; - if (left < chunkSize) - chunkSize = left; - - hsize_t offset[2] = {pos, 0}; - hsize_t size[2] = {chunkSize, static_cast(dim2)}; - - compressData(ndims, offset, size, - h5ivar, h5ispace, h5ovar, h5ospace, - h5native_type, buffer); - - pos += chunkSize; - } - - checkH5Err(H5Sclose(h5ispace)); - checkH5Err(H5Dclose(h5ivar)); - checkH5Err(H5Dclose(h5ovar)); - checkH5Err(H5Sclose(h5ospace)); +template +static void compressDataset(hid_t h5ifile, + hid_t h5ofile, + const char* varname, + hid_t h5native_type, + hid_t h5type, + unsigned int compressionLevel, + void* buffer, + unsigned long bufferSize) { + // Open original var + hid_t h5ivar = H5Dopen(h5ifile, varname, H5P_DEFAULT); + checkH5Err(h5ivar); + + // Get dimension information + hid_t h5ispace = H5Dget_space(h5ivar); + checkH5Err(h5ispace); + + int ndims = H5Sget_simple_extent_ndims(h5ispace); + checkH5Err(ndims); + if (ndims > 2) + logError() << "Dimension > 2 are not supported"; + + hsize_t extent[2]; + checkH5Err(H5Sget_simple_extent_dims(h5ispace, extent, 0L)); + hsize_t nelements = extent[0]; + unsigned int dim2 = 1; + if (ndims > 1) + dim2 = extent[1]; + + // Create new dataset + hsize_t chunkSize = bufferSize / dim2 / sizeof(T); + + hsize_t dims[2] = {nelements, dim2}; // Change this for other elements + hid_t h5ospace = H5Screate_simple(ndims, dims, 0L); + checkH5Err(h5ospace); + hid_t h5opcreate = H5Pcreate(H5P_DATASET_CREATE); + checkH5Err(h5opcreate); + hsize_t chunkDims[2] = {std::min(chunkSize, nelements), dim2}; + checkH5Err(H5Pset_chunk(h5opcreate, ndims, chunkDims)); + // checkH5Err(H5Pset_szip(h5opcreate, H5_SZIP_NN_OPTION_MASK, 4)); + checkH5Err(H5Pset_deflate(h5opcreate, compressionLevel)); + hid_t h5ovar = + H5Dcreate(h5ofile, varname, h5type, h5ospace, H5P_DEFAULT, h5opcreate, H5P_DEFAULT); + checkH5Err(h5ovar); + checkH5Err(H5Pclose(h5opcreate)); + + // Transfer data + unsigned long pos = 0; + while (pos < nelements) { + const unsigned long left = nelements - pos; + if (left < chunkSize) + chunkSize = left; + + hsize_t offset[2] = {pos, 0}; + hsize_t size[2] = {chunkSize, static_cast(dim2)}; + + compressData(ndims, offset, size, h5ivar, h5ispace, h5ovar, h5ospace, h5native_type, buffer); + + pos += chunkSize; + } + + checkH5Err(H5Sclose(h5ispace)); + checkH5Err(H5Dclose(h5ivar)); + checkH5Err(H5Dclose(h5ovar)); + checkH5Err(H5Sclose(h5ospace)); } -template -static void compressTimeDataset(hid_t h5ifile, hid_t h5ofile, - const char* varname, hid_t h5native_type, hid_t h5type, - unsigned int compressionLevel, void* buffer, unsigned long bufferSize) -{ - // Open original var - hid_t h5ivar = H5Dopen(h5ifile, varname, H5P_DEFAULT); - checkH5Err(h5ivar); - - // Get dimension information - hid_t h5ispace = H5Dget_space(h5ivar); - checkH5Err(h5ispace); - - int ndims = H5Sget_simple_extent_ndims(h5ispace); - checkH5Err(ndims); - if (ndims > 2) - logError() << "Dimension > 1 are not supported for time datasets"; - - hsize_t extent[2]; - checkH5Err(H5Sget_simple_extent_dims(h5ispace, extent, 0L)); - hsize_t timesteps = extent[0]; - hsize_t nelements = extent[1]; - - // Create new dataset - const hsize_t chunkSize = bufferSize / sizeof(T); - - hsize_t dims[2] = {timesteps, nelements}; // Change this for other elements - hid_t h5ospace = H5Screate_simple(ndims, dims, 0L); - checkH5Err(h5ospace); - hid_t h5opcreate = H5Pcreate(H5P_DATASET_CREATE); - checkH5Err(h5opcreate); - hsize_t chunkDims[2] = {1, std::min(chunkSize, nelements)}; - checkH5Err(H5Pset_chunk(h5opcreate, ndims, chunkDims)); -// checkH5Err(H5Pset_szip(h5opcreate, H5_SZIP_NN_OPTION_MASK, 4)); - checkH5Err(H5Pset_deflate(h5opcreate, compressionLevel)); - hid_t h5ovar = H5Dcreate(h5ofile, varname, h5type, h5ospace, - H5P_DEFAULT, h5opcreate, H5P_DEFAULT); - checkH5Err(h5ovar); - checkH5Err(H5Pclose(h5opcreate)); - - // Transfer data - for (unsigned long t = 0; t < timesteps; t++) { - unsigned long pos = 0; - while (pos < nelements) { - const unsigned long left = nelements - pos; - hsize_t tmpChunkSize = chunkSize; - - if (left < tmpChunkSize) - tmpChunkSize = left; - - hsize_t offset[2] = {t, pos}; - hsize_t size[2] = {1, tmpChunkSize}; - - compressData(ndims, offset, size, - h5ivar, h5ispace, h5ovar, h5ospace, - h5native_type, buffer); - - pos += chunkSize; - } - } - - checkH5Err(H5Sclose(h5ispace)); - checkH5Err(H5Dclose(h5ivar)); - checkH5Err(H5Dclose(h5ovar)); - checkH5Err(H5Sclose(h5ospace)); +template +static void compressTimeDataset(hid_t h5ifile, + hid_t h5ofile, + const char* varname, + hid_t h5native_type, + hid_t h5type, + unsigned int compressionLevel, + void* buffer, + unsigned long bufferSize) { + // Open original var + hid_t h5ivar = H5Dopen(h5ifile, varname, H5P_DEFAULT); + checkH5Err(h5ivar); + + // Get dimension information + hid_t h5ispace = H5Dget_space(h5ivar); + checkH5Err(h5ispace); + + int ndims = H5Sget_simple_extent_ndims(h5ispace); + checkH5Err(ndims); + if (ndims > 2) + logError() << "Dimension > 1 are not supported for time datasets"; + + hsize_t extent[2]; + checkH5Err(H5Sget_simple_extent_dims(h5ispace, extent, 0L)); + hsize_t timesteps = extent[0]; + hsize_t nelements = extent[1]; + + // Create new dataset + const hsize_t chunkSize = bufferSize / sizeof(T); + + hsize_t dims[2] = {timesteps, nelements}; // Change this for other elements + hid_t h5ospace = H5Screate_simple(ndims, dims, 0L); + checkH5Err(h5ospace); + hid_t h5opcreate = H5Pcreate(H5P_DATASET_CREATE); + checkH5Err(h5opcreate); + hsize_t chunkDims[2] = {1, std::min(chunkSize, nelements)}; + checkH5Err(H5Pset_chunk(h5opcreate, ndims, chunkDims)); + // checkH5Err(H5Pset_szip(h5opcreate, H5_SZIP_NN_OPTION_MASK, 4)); + checkH5Err(H5Pset_deflate(h5opcreate, compressionLevel)); + hid_t h5ovar = + H5Dcreate(h5ofile, varname, h5type, h5ospace, H5P_DEFAULT, h5opcreate, H5P_DEFAULT); + checkH5Err(h5ovar); + checkH5Err(H5Pclose(h5opcreate)); + + // Transfer data + for (unsigned long t = 0; t < timesteps; t++) { + unsigned long pos = 0; + while (pos < nelements) { + const unsigned long left = nelements - pos; + hsize_t tmpChunkSize = chunkSize; + + if (left < tmpChunkSize) + tmpChunkSize = left; + + hsize_t offset[2] = {t, pos}; + hsize_t size[2] = {1, tmpChunkSize}; + + compressData(ndims, offset, size, h5ivar, h5ispace, h5ovar, h5ospace, h5native_type, buffer); + + pos += chunkSize; + } + } + + checkH5Err(H5Sclose(h5ispace)); + checkH5Err(H5Dclose(h5ivar)); + checkH5Err(H5Dclose(h5ovar)); + checkH5Err(H5Sclose(h5ospace)); } -int main(int argc, char* argv[]) -{ - utils::Args args; - args.addOption("binary", 'b', "assume binary XDMF file", utils::Args::No, false); - args.addOption("level", 'l', "gzip compressen level [0-9]", utils::Args::Required, false); - args.addAdditionalOption("input", "input file"); - args.addAdditionalOption("output", "output file", false); - - switch (args.parse(argc, argv)) { - case utils::Args::Help: - return 0; - case utils::Args::Error: - return 1; - } - - std::string input = args.getAdditionalArgument("input"); - - bool binary = args.isSet("binary"); - if (!binary && utils::StringUtils::endsWith(input, ".h5")) - utils::StringUtils::replaceLast(input, ".h5", ".xdmf"); - - std::string output; - if (args.isSetAdditional("output")) - output = args.getAdditionalArgument("output"); - else { - output = input; - utils::StringUtils::replaceLast(output, ".xdmf", ""); - output += "_compressed.h5"; - } - - unsigned int compressionLevel = args.getArgument("level", 5); - - Input* inputHandle = 0L; - if (binary) - inputHandle = new BinaryInput(input); - else - inputHandle = new HDF5Input(input); - - std::vector variables = inputHandle->getVarList(); - - logInfo() << "Found dataset with" << inputHandle->numElements() << "elements and" - << inputHandle->numVertices() << "vertices"; - - hid_t h5ofile = H5Fcreate(output.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT); - checkH5Err(h5ofile); - - for (std::vector::const_iterator it = variables.begin(); - it != variables.end(); ++it) { - size_t numElements = inputHandle->numElements(); - hid_t type; - hid_t nativeType; - unsigned dim2 = 0; - bool isVertex = false; - - if (it->name == "connect") { - logInfo() << "Compressing connectivity..."; - type = H5T_STD_U64LE; - nativeType = H5T_NATIVE_UINT64; - dim2 = inputHandle->verticesPerElement(); - } else if (it->name == "geometry") { - logInfo() << "Compressing geometry..."; - numElements = inputHandle->numVertices(); - type = H5T_IEEE_F32LE; - nativeType = H5T_NATIVE_FLOAT; - dim2 = 3; - isVertex = true; - } else if (it->name == "partition") { - logInfo() << "Compressing partition..."; - type = H5T_STD_U32LE; - nativeType = H5T_NATIVE_UINT32; - } else { - logInfo() << "Compressing" << utils::nospace << it->name.c_str() << "..."; - type = H5T_IEEE_F32LE; - nativeType = H5T_NATIVE_FLOAT; - } - - OutputVar writer(h5ofile, it->name.c_str(), type, it->timesteps, numElements, dim2, - compressionLevel); - inputHandle->writeVariable(*it, nativeType, isVertex, writer); - } - - checkH5Err(H5Fclose(h5ofile)); - - delete inputHandle; - - return 0; +int main(int argc, char* argv[]) { + utils::Args args; + args.addOption("binary", 'b', "assume binary XDMF file", utils::Args::No, false); + args.addOption("level", 'l', "gzip compressen level [0-9]", utils::Args::Required, false); + args.addAdditionalOption("input", "input file"); + args.addAdditionalOption("output", "output file", false); + + switch (args.parse(argc, argv)) { + case utils::Args::Help: + return 0; + case utils::Args::Error: + return 1; + } + + std::string input = args.getAdditionalArgument("input"); + + bool binary = args.isSet("binary"); + if (!binary && utils::StringUtils::endsWith(input, ".h5")) + utils::StringUtils::replaceLast(input, ".h5", ".xdmf"); + + std::string output; + if (args.isSetAdditional("output")) + output = args.getAdditionalArgument("output"); + else { + output = input; + utils::StringUtils::replaceLast(output, ".xdmf", ""); + output += "_compressed.h5"; + } + + unsigned int compressionLevel = args.getArgument("level", 5); + + Input* inputHandle = 0L; + if (binary) + inputHandle = new BinaryInput(input); + else + inputHandle = new HDF5Input(input); + + std::vector variables = inputHandle->getVarList(); + + logInfo() << "Found dataset with" << inputHandle->numElements() << "elements and" + << inputHandle->numVertices() << "vertices"; + + hid_t h5ofile = H5Fcreate(output.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT); + checkH5Err(h5ofile); + + for (std::vector::const_iterator it = variables.begin(); it != variables.end(); ++it) { + size_t numElements = inputHandle->numElements(); + hid_t type; + hid_t nativeType; + unsigned dim2 = 0; + bool isVertex = false; + + if (it->name == "connect") { + logInfo() << "Compressing connectivity..."; + type = H5T_STD_U64LE; + nativeType = H5T_NATIVE_UINT64; + dim2 = inputHandle->verticesPerElement(); + } else if (it->name == "geometry") { + logInfo() << "Compressing geometry..."; + numElements = inputHandle->numVertices(); + type = H5T_IEEE_F32LE; + nativeType = H5T_NATIVE_FLOAT; + dim2 = 3; + isVertex = true; + } else if (it->name == "partition") { + logInfo() << "Compressing partition..."; + type = H5T_STD_U32LE; + nativeType = H5T_NATIVE_UINT32; + } else { + logInfo() << "Compressing" << utils::nospace << it->name.c_str() << "..."; + type = H5T_IEEE_F32LE; + nativeType = H5T_NATIVE_FLOAT; + } + + OutputVar writer( + h5ofile, it->name.c_str(), type, it->timesteps, numElements, dim2, compressionLevel); + inputHandle->writeVariable(*it, nativeType, isVertex, writer); + } + + checkH5Err(H5Fclose(h5ofile)); + + delete inputHandle; + + return 0; } diff --git a/compress_xmdf/src/output.h b/compress_xmdf/src/output.h index e91b251..80d0b39 100644 --- a/compress_xmdf/src/output.h +++ b/compress_xmdf/src/output.h @@ -2,7 +2,8 @@ * @file * This file is part of SeisSol. * - * @author Sebastian Rettenberger (sebastian.rettenberger AT tum.de, http://www5.in.tum.de/wiki/index.php/Sebastian_Rettenberger) + * @author Sebastian Rettenberger (sebastian.rettenberger AT tum.de, + * http://www5.in.tum.de/wiki/index.php/Sebastian_Rettenberger) * * @section LICENSE * Copyright (c) 2017, SeisSol Group @@ -40,91 +41,90 @@ #ifndef OUTPUT_H #define OUTPUT_H +#include "hdf5_helper.h" +#include "utils/logger.h" + #include -#include "utils/logger.h" +class OutputVar { + private: + hid_t m_var; + hid_t m_space; -#include "hdf5_helper.h" + unsigned int m_ndims; + + public: + OutputVar(hid_t file, + const char* name, + hid_t type, + unsigned int timesteps, + size_t nElements, + unsigned int dim2, + unsigned int compressionLevel) { + if (timesteps > 0 && dim2 > 0) + logError() << "Time data set with multiple dimensions is not supported."; + + m_ndims = 2; + if (timesteps == 0 && nElements == 0) + m_ndims = 1; + + bool hasTime = true; + if (timesteps == 0) { + timesteps = 1; + hasTime = false; + } + if (dim2 == 0) + dim2 = 1; + + // Create new dataset + hsize_t chunkSize = CHUNK_SIZE / dim2 / sizeofType(type); + + hsize_t dims[2]; + if (hasTime) { + dims[0] = timesteps; + dims[1] = nElements; + } else { + dims[0] = nElements; + dims[1] = dim2; + } + + m_space = H5Screate_simple(m_ndims, dims, 0L); + checkH5Err(m_space); + hid_t h5pcreate = H5Pcreate(H5P_DATASET_CREATE); + checkH5Err(h5pcreate); + hsize_t chunkDims[2]; + if (hasTime) { + chunkDims[0] = 1; + chunkDims[1] = std::min(chunkSize, static_cast(nElements)); + } else { + chunkDims[0] = std::min(chunkSize, static_cast(nElements)); + chunkDims[1] = dim2; + } + checkH5Err(H5Pset_chunk(h5pcreate, m_ndims, chunkDims)); + // checkH5Err(H5Pset_szip(h5opcreate, H5_SZIP_NN_OPTION_MASK, 4)); + checkH5Err(H5Pset_deflate(h5pcreate, compressionLevel)); + m_var = H5Dcreate(file, name, type, m_space, H5P_DEFAULT, h5pcreate, H5P_DEFAULT); + checkH5Err(m_var); + checkH5Err(H5Pclose(h5pcreate)); + } + + ~OutputVar() { + checkH5Err(H5Dclose(m_var)); + checkH5Err(H5Sclose(m_space)); + } + + void write(void* buffer, hid_t nativeType, hsize_t offset[2], hsize_t size[2]) { + hid_t memspace = H5Screate_simple(m_ndims, size, 0L); + checkH5Err(memspace); + + checkH5Err(H5Sselect_hyperslab(m_space, H5S_SELECT_SET, offset, 0L, size, 0L)); + checkH5Err(H5Dwrite(m_var, nativeType, memspace, m_space, H5P_DEFAULT, buffer)); + + checkH5Err(H5Sclose(memspace)); + } -class OutputVar -{ -private: - hid_t m_var; - hid_t m_space; - - unsigned int m_ndims; - -public: - OutputVar(hid_t file, const char* name, hid_t type, unsigned int timesteps, size_t nElements, unsigned int dim2, - unsigned int compressionLevel) - { - if (timesteps > 0 && dim2 > 0) - logError() << "Time data set with multiple dimensions is not supported."; - - m_ndims = 2; - if (timesteps == 0 && nElements == 0) - m_ndims = 1; - - bool hasTime = true; - if (timesteps == 0) { - timesteps = 1; - hasTime = false; - } - if (dim2 == 0) - dim2 = 1; - - // Create new dataset - hsize_t chunkSize = CHUNK_SIZE / dim2 / sizeofType(type); - - hsize_t dims[2]; - if (hasTime) { - dims[0] = timesteps; - dims[1] = nElements; - } else { - dims[0] = nElements; - dims[1] = dim2; - } - - m_space = H5Screate_simple(m_ndims, dims, 0L); - checkH5Err(m_space); - hid_t h5pcreate = H5Pcreate(H5P_DATASET_CREATE); - checkH5Err(h5pcreate); - hsize_t chunkDims[2]; - if (hasTime) { - chunkDims[0] = 1; - chunkDims[1] = std::min(chunkSize, static_cast(nElements)); - } else { - chunkDims[0] = std::min(chunkSize, static_cast(nElements)); - chunkDims[1] = dim2; - } - checkH5Err(H5Pset_chunk(h5pcreate, m_ndims, chunkDims)); -// checkH5Err(H5Pset_szip(h5opcreate, H5_SZIP_NN_OPTION_MASK, 4)); - checkH5Err(H5Pset_deflate(h5pcreate, compressionLevel)); - m_var = H5Dcreate(file, name, type, m_space, - H5P_DEFAULT, h5pcreate, H5P_DEFAULT); - checkH5Err(m_var); - checkH5Err(H5Pclose(h5pcreate)); - } - - ~OutputVar() - { - checkH5Err(H5Dclose(m_var)); - checkH5Err(H5Sclose(m_space)); - } - - void write(void* buffer, hid_t nativeType, hsize_t offset[2], hsize_t size[2]) - { - hid_t memspace = H5Screate_simple(m_ndims, size, 0L); - checkH5Err(memspace); - - checkH5Err(H5Sselect_hyperslab(m_space, H5S_SELECT_SET, offset, 0L, size, 0L)); - checkH5Err(H5Dwrite(m_var, nativeType, memspace, m_space, H5P_DEFAULT, buffer)); - - checkH5Err(H5Sclose(memspace)); - } - -public: - static const size_t CHUNK_SIZE = 256*1024*1024; + public: + static const size_t CHUNK_SIZE = 256 * 1024 * 1024; }; -#endif // OUTPUT_H \ No newline at end of file +#endif // OUTPUT_H diff --git a/external/stb_image_write.h b/external/stb_image_write.h index 95943eb..36d64f8 100644 --- a/external/stb_image_write.h +++ b/external/stb_image_write.h @@ -26,11 +26,12 @@ You can #define STBIW_MALLOC(), STBIW_REALLOC(), and STBIW_FREE() to replace malloc,realloc,free. You can #define STBIW_MEMMOVE() to replace memmove() - You can #define STBIW_ZLIB_COMPRESS to use a custom zlib-style compress function - for PNG compression (instead of the builtin one), it must have the following signature: - unsigned char * my_compress(unsigned char *data, int data_len, int *out_len, int quality); - The returned data will be freed with STBIW_FREE() (free() by default), - so it must be heap allocated with STBIW_MALLOC() (malloc() by default), + You can #define STBIW_ZLIB_COMPRESS to use a custom zlib-style compress +function for PNG compression (instead of the builtin one), it must have the +following signature: unsigned char * my_compress(unsigned char *data, int +data_len, int *out_len, int quality); The returned data will be freed with +STBIW_FREE() (free() by default), so it must be heap allocated with +STBIW_MALLOC() (malloc() by default), UNICODE: @@ -44,30 +45,37 @@ There are five functions, one for each image file format: - int stbi_write_png(char const *filename, int w, int h, int comp, const void *data, int stride_in_bytes); - int stbi_write_bmp(char const *filename, int w, int h, int comp, const void *data); - int stbi_write_tga(char const *filename, int w, int h, int comp, const void *data); - int stbi_write_jpg(char const *filename, int w, int h, int comp, const void *data, int quality); - int stbi_write_hdr(char const *filename, int w, int h, int comp, const float *data); - - void stbi_flip_vertically_on_write(int flag); // flag is non-zero to flip data vertically - - There are also five equivalent functions that use an arbitrary write function. You are - expected to open/close your file-equivalent before and after calling these: - - int stbi_write_png_to_func(stbi_write_func *func, void *context, int w, int h, int comp, const void *data, int stride_in_bytes); - int stbi_write_bmp_to_func(stbi_write_func *func, void *context, int w, int h, int comp, const void *data); - int stbi_write_tga_to_func(stbi_write_func *func, void *context, int w, int h, int comp, const void *data); - int stbi_write_hdr_to_func(stbi_write_func *func, void *context, int w, int h, int comp, const float *data); - int stbi_write_jpg_to_func(stbi_write_func *func, void *context, int x, int y, int comp, const void *data, int quality); + int stbi_write_png(char const *filename, int w, int h, int comp, const void +*data, int stride_in_bytes); int stbi_write_bmp(char const *filename, int w, int +h, int comp, const void *data); int stbi_write_tga(char const *filename, int w, +int h, int comp, const void *data); int stbi_write_jpg(char const *filename, int +w, int h, int comp, const void *data, int quality); int stbi_write_hdr(char +const *filename, int w, int h, int comp, const float *data); + + void stbi_flip_vertically_on_write(int flag); // flag is non-zero to flip +data vertically + + There are also five equivalent functions that use an arbitrary write +function. You are expected to open/close your file-equivalent before and after +calling these: + + int stbi_write_png_to_func(stbi_write_func *func, void *context, int w, int +h, int comp, const void *data, int stride_in_bytes); int +stbi_write_bmp_to_func(stbi_write_func *func, void *context, int w, int h, int +comp, const void *data); int stbi_write_tga_to_func(stbi_write_func *func, void +*context, int w, int h, int comp, const void *data); int +stbi_write_hdr_to_func(stbi_write_func *func, void *context, int w, int h, int +comp, const float *data); int stbi_write_jpg_to_func(stbi_write_func *func, void +*context, int x, int y, int comp, const void *data, int quality); where the callback is: void stbi_write_func(void *context, void *data, int size); You can configure it with these global variables: - int stbi_write_tga_with_rle; // defaults to true; set to 0 to disable RLE - int stbi_write_png_compression_level; // defaults to 8; set to higher for more compression - int stbi_write_force_png_filter; // defaults to -1; set to 0..5 to force a filter mode + int stbi_write_tga_with_rle; // defaults to true; set to 0 to +disable RLE int stbi_write_png_compression_level; // defaults to 8; set to +higher for more compression int stbi_write_force_png_filter; // defaults +to -1; set to 0..5 to force a filter mode You can define STBI_WRITE_NO_STDIO to disable the file variant of these @@ -152,99 +160,113 @@ LICENSE #include -// if STB_IMAGE_WRITE_STATIC causes problems, try defining STBIWDEF to 'inline' or 'static inline' +// if STB_IMAGE_WRITE_STATIC causes problems, try defining STBIWDEF to 'inline' +// or 'static inline' #ifndef STBIWDEF #ifdef STB_IMAGE_WRITE_STATIC -#define STBIWDEF static +#define STBIWDEF static #else #ifdef __cplusplus -#define STBIWDEF extern "C" +#define STBIWDEF extern "C" #else -#define STBIWDEF extern +#define STBIWDEF extern #endif #endif #endif -#ifndef STB_IMAGE_WRITE_STATIC // C++ forbids static forward declarations +#ifndef STB_IMAGE_WRITE_STATIC // C++ forbids static forward declarations extern int stbi_write_tga_with_rle; extern int stbi_write_png_compression_level; extern int stbi_write_force_png_filter; #endif #ifndef STBI_WRITE_NO_STDIO -STBIWDEF int stbi_write_png(char const *filename, int w, int h, int comp, const void *data, int stride_in_bytes); -STBIWDEF int stbi_write_bmp(char const *filename, int w, int h, int comp, const void *data); -STBIWDEF int stbi_write_tga(char const *filename, int w, int h, int comp, const void *data); -STBIWDEF int stbi_write_hdr(char const *filename, int w, int h, int comp, const float *data); -STBIWDEF int stbi_write_jpg(char const *filename, int x, int y, int comp, const void *data, int quality); +STBIWDEF int stbi_write_png( + const char* filename, int w, int h, int comp, const void* data, int stride_in_bytes); +STBIWDEF int stbi_write_bmp(const char* filename, int w, int h, int comp, const void* data); +STBIWDEF int stbi_write_tga(const char* filename, int w, int h, int comp, const void* data); +STBIWDEF int stbi_write_hdr(const char* filename, int w, int h, int comp, const float* data); +STBIWDEF int + stbi_write_jpg(const char* filename, int x, int y, int comp, const void* data, int quality); #ifdef STBI_WINDOWS_UTF8 -STBIWDEF int stbiw_convert_wchar_to_utf8(char *buffer, size_t bufferlen, const wchar_t* input); +STBIWDEF int stbiw_convert_wchar_to_utf8(char* buffer, size_t bufferlen, const wchar_t* input); #endif #endif -typedef void stbi_write_func(void *context, void *data, int size); - -STBIWDEF int stbi_write_png_to_func(stbi_write_func *func, void *context, int w, int h, int comp, const void *data, int stride_in_bytes); -STBIWDEF int stbi_write_bmp_to_func(stbi_write_func *func, void *context, int w, int h, int comp, const void *data); -STBIWDEF int stbi_write_tga_to_func(stbi_write_func *func, void *context, int w, int h, int comp, const void *data); -STBIWDEF int stbi_write_hdr_to_func(stbi_write_func *func, void *context, int w, int h, int comp, const float *data); -STBIWDEF int stbi_write_jpg_to_func(stbi_write_func *func, void *context, int x, int y, int comp, const void *data, int quality); +typedef void stbi_write_func(void* context, void* data, int size); + +STBIWDEF int stbi_write_png_to_func(stbi_write_func* func, + void* context, + int w, + int h, + int comp, + const void* data, + int stride_in_bytes); +STBIWDEF int stbi_write_bmp_to_func( + stbi_write_func* func, void* context, int w, int h, int comp, const void* data); +STBIWDEF int stbi_write_tga_to_func( + stbi_write_func* func, void* context, int w, int h, int comp, const void* data); +STBIWDEF int stbi_write_hdr_to_func( + stbi_write_func* func, void* context, int w, int h, int comp, const float* data); +STBIWDEF int stbi_write_jpg_to_func( + stbi_write_func* func, void* context, int x, int y, int comp, const void* data, int quality); STBIWDEF void stbi_flip_vertically_on_write(int flip_boolean); -#endif//INCLUDE_STB_IMAGE_WRITE_H +#endif // INCLUDE_STB_IMAGE_WRITE_H #ifdef STB_IMAGE_WRITE_IMPLEMENTATION #ifdef _WIN32 - #ifndef _CRT_SECURE_NO_WARNINGS - #define _CRT_SECURE_NO_WARNINGS - #endif - #ifndef _CRT_NONSTDC_NO_DEPRECATE - #define _CRT_NONSTDC_NO_DEPRECATE - #endif +#ifndef _CRT_SECURE_NO_WARNINGS +#define _CRT_SECURE_NO_WARNINGS +#endif +#ifndef _CRT_NONSTDC_NO_DEPRECATE +#define _CRT_NONSTDC_NO_DEPRECATE +#endif #endif #ifndef STBI_WRITE_NO_STDIO #include #endif // STBI_WRITE_NO_STDIO +#include #include #include #include -#include -#if defined(STBIW_MALLOC) && defined(STBIW_FREE) && (defined(STBIW_REALLOC) || defined(STBIW_REALLOC_SIZED)) +#if defined(STBIW_MALLOC) && defined(STBIW_FREE) && \ + (defined(STBIW_REALLOC) || defined(STBIW_REALLOC_SIZED)) // ok -#elif !defined(STBIW_MALLOC) && !defined(STBIW_FREE) && !defined(STBIW_REALLOC) && !defined(STBIW_REALLOC_SIZED) +#elif !defined(STBIW_MALLOC) && !defined(STBIW_FREE) && !defined(STBIW_REALLOC) && \ + !defined(STBIW_REALLOC_SIZED) // ok #else -#error "Must define all or none of STBIW_MALLOC, STBIW_FREE, and STBIW_REALLOC (or STBIW_REALLOC_SIZED)." +#error \ + "Must define all or none of STBIW_MALLOC, STBIW_FREE, and STBIW_REALLOC (or STBIW_REALLOC_SIZED)." #endif #ifndef STBIW_MALLOC -#define STBIW_MALLOC(sz) malloc(sz) -#define STBIW_REALLOC(p,newsz) realloc(p,newsz) -#define STBIW_FREE(p) free(p) +#define STBIW_MALLOC(sz) malloc(sz) +#define STBIW_REALLOC(p, newsz) realloc(p, newsz) +#define STBIW_FREE(p) free(p) #endif #ifndef STBIW_REALLOC_SIZED -#define STBIW_REALLOC_SIZED(p,oldsz,newsz) STBIW_REALLOC(p,newsz) +#define STBIW_REALLOC_SIZED(p, oldsz, newsz) STBIW_REALLOC(p, newsz) #endif - #ifndef STBIW_MEMMOVE -#define STBIW_MEMMOVE(a,b,sz) memmove(a,b,sz) +#define STBIW_MEMMOVE(a, b, sz) memmove(a, b, sz) #endif - #ifndef STBIW_ASSERT #include #define STBIW_ASSERT(x) assert(x) #endif -#define STBIW_UCHAR(x) (unsigned char) ((x) & 0xff) +#define STBIW_UCHAR(x) (unsigned char)((x) & 0xff) #ifdef STB_IMAGE_WRITE_STATIC static int stbi_write_png_compression_level = 8; @@ -258,31 +280,25 @@ int stbi_write_force_png_filter = -1; static int stbi__flip_vertically_on_write = 0; -STBIWDEF void stbi_flip_vertically_on_write(int flag) -{ - stbi__flip_vertically_on_write = flag; -} +STBIWDEF void stbi_flip_vertically_on_write(int flag) { stbi__flip_vertically_on_write = flag; } -typedef struct -{ - stbi_write_func *func; - void *context; - unsigned char buffer[64]; - int buf_used; +typedef struct { + stbi_write_func* func; + void* context; + unsigned char buffer[64]; + int buf_used; } stbi__write_context; // initialize a callback-based context -static void stbi__start_write_callbacks(stbi__write_context *s, stbi_write_func *c, void *context) -{ - s->func = c; - s->context = context; +static void stbi__start_write_callbacks(stbi__write_context* s, stbi_write_func* c, void* context) { + s->func = c; + s->context = context; } #ifndef STBI_WRITE_NO_STDIO -static void stbi__stdio_write(void *context, void *data, int size) -{ - fwrite(data,1,size,(FILE*) context); +static void stbi__stdio_write(void* context, void* data, int size) { + fwrite(data, 1, size, (FILE*)context); } #if defined(_MSC_VER) && defined(STBI_WINDOWS_UTF8) @@ -291,328 +307,411 @@ static void stbi__stdio_write(void *context, void *data, int size) #else #define STBIW_EXTERN extern #endif -STBIW_EXTERN __declspec(dllimport) int __stdcall MultiByteToWideChar(unsigned int cp, unsigned long flags, const char *str, int cbmb, wchar_t *widestr, int cchwide); -STBIW_EXTERN __declspec(dllimport) int __stdcall WideCharToMultiByte(unsigned int cp, unsigned long flags, const wchar_t *widestr, int cchwide, char *str, int cbmb, const char *defchar, int *used_default); - -STBIWDEF int stbiw_convert_wchar_to_utf8(char *buffer, size_t bufferlen, const wchar_t* input) -{ - return WideCharToMultiByte(65001 /* UTF8 */, 0, input, -1, buffer, (int) bufferlen, NULL, NULL); +STBIW_EXTERN __declspec(dllimport) int __stdcall MultiByteToWideChar( + unsigned int cp, unsigned long flags, const char* str, int cbmb, wchar_t* widestr, int cchwide); +STBIW_EXTERN __declspec(dllimport) int __stdcall WideCharToMultiByte(unsigned int cp, + unsigned long flags, + const wchar_t* widestr, + int cchwide, + char* str, + int cbmb, + const char* defchar, + int* used_default); + +STBIWDEF int stbiw_convert_wchar_to_utf8(char* buffer, size_t bufferlen, const wchar_t* input) { + return WideCharToMultiByte(65001 /* UTF8 */, 0, input, -1, buffer, (int)bufferlen, NULL, NULL); } #endif -static FILE *stbiw__fopen(char const *filename, char const *mode) -{ - FILE *f; +static FILE* stbiw__fopen(const char* filename, const char* mode) { + FILE* f; #if defined(_MSC_VER) && defined(STBI_WINDOWS_UTF8) - wchar_t wMode[64]; - wchar_t wFilename[1024]; - if (0 == MultiByteToWideChar(65001 /* UTF8 */, 0, filename, -1, wFilename, sizeof(wFilename))) - return 0; + wchar_t wMode[64]; + wchar_t wFilename[1024]; + if (0 == MultiByteToWideChar(65001 /* UTF8 */, 0, filename, -1, wFilename, sizeof(wFilename))) + return 0; - if (0 == MultiByteToWideChar(65001 /* UTF8 */, 0, mode, -1, wMode, sizeof(wMode))) - return 0; + if (0 == MultiByteToWideChar(65001 /* UTF8 */, 0, mode, -1, wMode, sizeof(wMode))) + return 0; #if _MSC_VER >= 1400 - if (0 != _wfopen_s(&f, wFilename, wMode)) - f = 0; + if (0 != _wfopen_s(&f, wFilename, wMode)) + f = 0; #else - f = _wfopen(wFilename, wMode); + f = _wfopen(wFilename, wMode); #endif #elif defined(_MSC_VER) && _MSC_VER >= 1400 - if (0 != fopen_s(&f, filename, mode)) - f=0; + if (0 != fopen_s(&f, filename, mode)) + f = 0; #else - f = fopen(filename, mode); + f = fopen(filename, mode); #endif - return f; + return f; } -static int stbi__start_write_file(stbi__write_context *s, const char *filename) -{ - FILE *f = stbiw__fopen(filename, "wb"); - stbi__start_write_callbacks(s, stbi__stdio_write, (void *) f); - return f != NULL; +static int stbi__start_write_file(stbi__write_context* s, const char* filename) { + FILE* f = stbiw__fopen(filename, "wb"); + stbi__start_write_callbacks(s, stbi__stdio_write, (void*)f); + return f != NULL; } -static void stbi__end_write_file(stbi__write_context *s) -{ - fclose((FILE *)s->context); -} +static void stbi__end_write_file(stbi__write_context* s) { fclose((FILE*)s->context); } #endif // !STBI_WRITE_NO_STDIO typedef unsigned int stbiw_uint32; -typedef int stb_image_write_test[sizeof(stbiw_uint32)==4 ? 1 : -1]; - -static void stbiw__writefv(stbi__write_context *s, const char *fmt, va_list v) -{ - while (*fmt) { - switch (*fmt++) { - case ' ': break; - case '1': { unsigned char x = STBIW_UCHAR(va_arg(v, int)); - s->func(s->context,&x,1); - break; } - case '2': { int x = va_arg(v,int); - unsigned char b[2]; - b[0] = STBIW_UCHAR(x); - b[1] = STBIW_UCHAR(x>>8); - s->func(s->context,b,2); - break; } - case '4': { stbiw_uint32 x = va_arg(v,int); - unsigned char b[4]; - b[0]=STBIW_UCHAR(x); - b[1]=STBIW_UCHAR(x>>8); - b[2]=STBIW_UCHAR(x>>16); - b[3]=STBIW_UCHAR(x>>24); - s->func(s->context,b,4); - break; } - default: - STBIW_ASSERT(0); - return; - } - } +typedef int stb_image_write_test[sizeof(stbiw_uint32) == 4 ? 1 : -1]; + +static void stbiw__writefv(stbi__write_context* s, const char* fmt, va_list v) { + while (*fmt) { + switch (*fmt++) { + case ' ': + break; + case '1': { + unsigned char x = STBIW_UCHAR(va_arg(v, int)); + s->func(s->context, &x, 1); + break; + } + case '2': { + int x = va_arg(v, int); + unsigned char b[2]; + b[0] = STBIW_UCHAR(x); + b[1] = STBIW_UCHAR(x >> 8); + s->func(s->context, b, 2); + break; + } + case '4': { + stbiw_uint32 x = va_arg(v, int); + unsigned char b[4]; + b[0] = STBIW_UCHAR(x); + b[1] = STBIW_UCHAR(x >> 8); + b[2] = STBIW_UCHAR(x >> 16); + b[3] = STBIW_UCHAR(x >> 24); + s->func(s->context, b, 4); + break; + } + default: + STBIW_ASSERT(0); + return; + } + } } -static void stbiw__writef(stbi__write_context *s, const char *fmt, ...) -{ - va_list v; - va_start(v, fmt); - stbiw__writefv(s, fmt, v); - va_end(v); +static void stbiw__writef(stbi__write_context* s, const char* fmt, ...) { + va_list v; + va_start(v, fmt); + stbiw__writefv(s, fmt, v); + va_end(v); } -static void stbiw__write_flush(stbi__write_context *s) -{ - if (s->buf_used) { - s->func(s->context, &s->buffer, s->buf_used); - s->buf_used = 0; - } +static void stbiw__write_flush(stbi__write_context* s) { + if (s->buf_used) { + s->func(s->context, &s->buffer, s->buf_used); + s->buf_used = 0; + } } -static void stbiw__putc(stbi__write_context *s, unsigned char c) -{ - s->func(s->context, &c, 1); -} +static void stbiw__putc(stbi__write_context* s, unsigned char c) { s->func(s->context, &c, 1); } -static void stbiw__write1(stbi__write_context *s, unsigned char a) -{ - if (s->buf_used + 1 > sizeof(s->buffer)) - stbiw__write_flush(s); - s->buffer[s->buf_used++] = a; +static void stbiw__write1(stbi__write_context* s, unsigned char a) { + if (s->buf_used + 1 > sizeof(s->buffer)) + stbiw__write_flush(s); + s->buffer[s->buf_used++] = a; } -static void stbiw__write3(stbi__write_context *s, unsigned char a, unsigned char b, unsigned char c) -{ - int n; - if (s->buf_used + 3 > sizeof(s->buffer)) - stbiw__write_flush(s); - n = s->buf_used; - s->buf_used = n+3; - s->buffer[n+0] = a; - s->buffer[n+1] = b; - s->buffer[n+2] = c; +static void + stbiw__write3(stbi__write_context* s, unsigned char a, unsigned char b, unsigned char c) { + int n; + if (s->buf_used + 3 > sizeof(s->buffer)) + stbiw__write_flush(s); + n = s->buf_used; + s->buf_used = n + 3; + s->buffer[n + 0] = a; + s->buffer[n + 1] = b; + s->buffer[n + 2] = c; } -static void stbiw__write_pixel(stbi__write_context *s, int rgb_dir, int comp, int write_alpha, int expand_mono, unsigned char *d) -{ - unsigned char bg[3] = { 255, 0, 255}, px[3]; - int k; - - if (write_alpha < 0) - stbiw__write1(s, d[comp - 1]); - - switch (comp) { - case 2: // 2 pixels = mono + alpha, alpha is written separately, so same as 1-channel case - case 1: - if (expand_mono) - stbiw__write3(s, d[0], d[0], d[0]); // monochrome bmp - else - stbiw__write1(s, d[0]); // monochrome TGA - break; - case 4: - if (!write_alpha) { - // composite against pink background - for (k = 0; k < 3; ++k) - px[k] = bg[k] + ((d[k] - bg[k]) * d[3]) / 255; - stbiw__write3(s, px[1 - rgb_dir], px[1], px[1 + rgb_dir]); - break; - } - /* FALLTHROUGH */ - case 3: - stbiw__write3(s, d[1 - rgb_dir], d[1], d[1 + rgb_dir]); - break; - } - if (write_alpha > 0) - stbiw__write1(s, d[comp - 1]); +static void stbiw__write_pixel(stbi__write_context* s, + int rgb_dir, + int comp, + int write_alpha, + int expand_mono, + unsigned char* d) { + unsigned char bg[3] = {255, 0, 255}, px[3]; + int k; + + if (write_alpha < 0) + stbiw__write1(s, d[comp - 1]); + + switch (comp) { + case 2: // 2 pixels = mono + alpha, alpha is written separately, so same as + // 1-channel case + case 1: + if (expand_mono) + stbiw__write3(s, d[0], d[0], d[0]); // monochrome bmp + else + stbiw__write1(s, d[0]); // monochrome TGA + break; + case 4: + if (!write_alpha) { + // composite against pink background + for (k = 0; k < 3; ++k) + px[k] = bg[k] + ((d[k] - bg[k]) * d[3]) / 255; + stbiw__write3(s, px[1 - rgb_dir], px[1], px[1 + rgb_dir]); + break; + } + /* FALLTHROUGH */ + case 3: + stbiw__write3(s, d[1 - rgb_dir], d[1], d[1 + rgb_dir]); + break; + } + if (write_alpha > 0) + stbiw__write1(s, d[comp - 1]); } -static void stbiw__write_pixels(stbi__write_context *s, int rgb_dir, int vdir, int x, int y, int comp, void *data, int write_alpha, int scanline_pad, int expand_mono) -{ - stbiw_uint32 zero = 0; - int i,j, j_end; - - if (y <= 0) - return; - - if (stbi__flip_vertically_on_write) - vdir *= -1; - - if (vdir < 0) { - j_end = -1; j = y-1; - } else { - j_end = y; j = 0; - } - - for (; j != j_end; j += vdir) { - for (i=0; i < x; ++i) { - unsigned char *d = (unsigned char *) data + (j*x+i)*comp; - stbiw__write_pixel(s, rgb_dir, comp, write_alpha, expand_mono, d); - } - stbiw__write_flush(s); - s->func(s->context, &zero, scanline_pad); - } +static void stbiw__write_pixels(stbi__write_context* s, + int rgb_dir, + int vdir, + int x, + int y, + int comp, + void* data, + int write_alpha, + int scanline_pad, + int expand_mono) { + stbiw_uint32 zero = 0; + int i, j, j_end; + + if (y <= 0) + return; + + if (stbi__flip_vertically_on_write) + vdir *= -1; + + if (vdir < 0) { + j_end = -1; + j = y - 1; + } else { + j_end = y; + j = 0; + } + + for (; j != j_end; j += vdir) { + for (i = 0; i < x; ++i) { + unsigned char* d = (unsigned char*)data + (j * x + i) * comp; + stbiw__write_pixel(s, rgb_dir, comp, write_alpha, expand_mono, d); + } + stbiw__write_flush(s); + s->func(s->context, &zero, scanline_pad); + } } -static int stbiw__outfile(stbi__write_context *s, int rgb_dir, int vdir, int x, int y, int comp, int expand_mono, void *data, int alpha, int pad, const char *fmt, ...) -{ - if (y < 0 || x < 0) { - return 0; - } else { - va_list v; - va_start(v, fmt); - stbiw__writefv(s, fmt, v); - va_end(v); - stbiw__write_pixels(s,rgb_dir,vdir,x,y,comp,data,alpha,pad, expand_mono); - return 1; - } +static int stbiw__outfile(stbi__write_context* s, + int rgb_dir, + int vdir, + int x, + int y, + int comp, + int expand_mono, + void* data, + int alpha, + int pad, + const char* fmt, + ...) { + if (y < 0 || x < 0) { + return 0; + } else { + va_list v; + va_start(v, fmt); + stbiw__writefv(s, fmt, v); + va_end(v); + stbiw__write_pixels(s, rgb_dir, vdir, x, y, comp, data, alpha, pad, expand_mono); + return 1; + } } -static int stbi_write_bmp_core(stbi__write_context *s, int x, int y, int comp, const void *data) -{ - int pad = (-x*3) & 3; - return stbiw__outfile(s,-1,-1,x,y,comp,1,(void *) data,0,pad, - "11 4 22 4" "4 44 22 444444", - 'B', 'M', 14+40+(x*3+pad)*y, 0,0, 14+40, // file header - 40, x,y, 1,24, 0,0,0,0,0,0); // bitmap header +static int stbi_write_bmp_core(stbi__write_context* s, int x, int y, int comp, const void* data) { + int pad = (-x * 3) & 3; + return stbiw__outfile(s, + -1, + -1, + x, + y, + comp, + 1, + (void*)data, + 0, + pad, + "11 4 22 4" + "4 44 22 444444", + 'B', + 'M', + 14 + 40 + (x * 3 + pad) * y, + 0, + 0, + 14 + 40, // file header + 40, + x, + y, + 1, + 24, + 0, + 0, + 0, + 0, + 0, + 0); // bitmap header } -STBIWDEF int stbi_write_bmp_to_func(stbi_write_func *func, void *context, int x, int y, int comp, const void *data) -{ - stbi__write_context s = { 0 }; - stbi__start_write_callbacks(&s, func, context); - return stbi_write_bmp_core(&s, x, y, comp, data); +STBIWDEF int stbi_write_bmp_to_func( + stbi_write_func* func, void* context, int x, int y, int comp, const void* data) { + stbi__write_context s = {0}; + stbi__start_write_callbacks(&s, func, context); + return stbi_write_bmp_core(&s, x, y, comp, data); } #ifndef STBI_WRITE_NO_STDIO -STBIWDEF int stbi_write_bmp(char const *filename, int x, int y, int comp, const void *data) -{ - stbi__write_context s = { 0 }; - if (stbi__start_write_file(&s,filename)) { - int r = stbi_write_bmp_core(&s, x, y, comp, data); - stbi__end_write_file(&s); - return r; - } else - return 0; +STBIWDEF int stbi_write_bmp(const char* filename, int x, int y, int comp, const void* data) { + stbi__write_context s = {0}; + if (stbi__start_write_file(&s, filename)) { + int r = stbi_write_bmp_core(&s, x, y, comp, data); + stbi__end_write_file(&s); + return r; + } else + return 0; } -#endif //!STBI_WRITE_NO_STDIO - -static int stbi_write_tga_core(stbi__write_context *s, int x, int y, int comp, void *data) -{ - int has_alpha = (comp == 2 || comp == 4); - int colorbytes = has_alpha ? comp-1 : comp; - int format = colorbytes < 2 ? 3 : 2; // 3 color channels (RGB/RGBA) = 2, 1 color channel (Y/YA) = 3 - - if (y < 0 || x < 0) - return 0; - - if (!stbi_write_tga_with_rle) { - return stbiw__outfile(s, -1, -1, x, y, comp, 0, (void *) data, has_alpha, 0, - "111 221 2222 11", 0, 0, format, 0, 0, 0, 0, 0, x, y, (colorbytes + has_alpha) * 8, has_alpha * 8); - } else { - int i,j,k; - int jend, jdir; - - stbiw__writef(s, "111 221 2222 11", 0,0,format+8, 0,0,0, 0,0,x,y, (colorbytes + has_alpha) * 8, has_alpha * 8); - - if (stbi__flip_vertically_on_write) { - j = 0; - jend = y; - jdir = 1; - } else { - j = y-1; - jend = -1; - jdir = -1; - } - for (; j != jend; j += jdir) { - unsigned char *row = (unsigned char *) data + j * x * comp; - int len; - - for (i = 0; i < x; i += len) { - unsigned char *begin = row + i * comp; - int diff = 1; - len = 1; - - if (i < x - 1) { - ++len; - diff = memcmp(begin, row + (i + 1) * comp, comp); - if (diff) { - const unsigned char *prev = begin; - for (k = i + 2; k < x && len < 128; ++k) { - if (memcmp(prev, row + k * comp, comp)) { - prev += comp; - ++len; - } else { - --len; - break; - } - } - } else { - for (k = i + 2; k < x && len < 128; ++k) { - if (!memcmp(begin, row + k * comp, comp)) { - ++len; - } else { - break; - } - } - } +#endif //! STBI_WRITE_NO_STDIO + +static int stbi_write_tga_core(stbi__write_context* s, int x, int y, int comp, void* data) { + int has_alpha = (comp == 2 || comp == 4); + int colorbytes = has_alpha ? comp - 1 : comp; + int format = + colorbytes < 2 ? 3 : 2; // 3 color channels (RGB/RGBA) = 2, 1 color channel (Y/YA) = 3 + + if (y < 0 || x < 0) + return 0; + + if (!stbi_write_tga_with_rle) { + return stbiw__outfile(s, + -1, + -1, + x, + y, + comp, + 0, + (void*)data, + has_alpha, + 0, + "111 221 2222 11", + 0, + 0, + format, + 0, + 0, + 0, + 0, + 0, + x, + y, + (colorbytes + has_alpha) * 8, + has_alpha * 8); + } else { + int i, j, k; + int jend, jdir; + + stbiw__writef(s, + "111 221 2222 11", + 0, + 0, + format + 8, + 0, + 0, + 0, + 0, + 0, + x, + y, + (colorbytes + has_alpha) * 8, + has_alpha * 8); + + if (stbi__flip_vertically_on_write) { + j = 0; + jend = y; + jdir = 1; + } else { + j = y - 1; + jend = -1; + jdir = -1; + } + for (; j != jend; j += jdir) { + unsigned char* row = (unsigned char*)data + j * x * comp; + int len; + + for (i = 0; i < x; i += len) { + unsigned char* begin = row + i * comp; + int diff = 1; + len = 1; + + if (i < x - 1) { + ++len; + diff = memcmp(begin, row + (i + 1) * comp, comp); + if (diff) { + const unsigned char* prev = begin; + for (k = i + 2; k < x && len < 128; ++k) { + if (memcmp(prev, row + k * comp, comp)) { + prev += comp; + ++len; + } else { + --len; + break; + } } - - if (diff) { - unsigned char header = STBIW_UCHAR(len - 1); - stbiw__write1(s, header); - for (k = 0; k < len; ++k) { - stbiw__write_pixel(s, -1, comp, has_alpha, 0, begin + k * comp); - } - } else { - unsigned char header = STBIW_UCHAR(len - 129); - stbiw__write1(s, header); - stbiw__write_pixel(s, -1, comp, has_alpha, 0, begin); + } else { + for (k = i + 2; k < x && len < 128; ++k) { + if (!memcmp(begin, row + k * comp, comp)) { + ++len; + } else { + break; + } } - } + } + } + + if (diff) { + unsigned char header = STBIW_UCHAR(len - 1); + stbiw__write1(s, header); + for (k = 0; k < len; ++k) { + stbiw__write_pixel(s, -1, comp, has_alpha, 0, begin + k * comp); + } + } else { + unsigned char header = STBIW_UCHAR(len - 129); + stbiw__write1(s, header); + stbiw__write_pixel(s, -1, comp, has_alpha, 0, begin); + } } - stbiw__write_flush(s); - } - return 1; + } + stbiw__write_flush(s); + } + return 1; } -STBIWDEF int stbi_write_tga_to_func(stbi_write_func *func, void *context, int x, int y, int comp, const void *data) -{ - stbi__write_context s = { 0 }; - stbi__start_write_callbacks(&s, func, context); - return stbi_write_tga_core(&s, x, y, comp, (void *) data); +STBIWDEF int stbi_write_tga_to_func( + stbi_write_func* func, void* context, int x, int y, int comp, const void* data) { + stbi__write_context s = {0}; + stbi__start_write_callbacks(&s, func, context); + return stbi_write_tga_core(&s, x, y, comp, (void*)data); } #ifndef STBI_WRITE_NO_STDIO -STBIWDEF int stbi_write_tga(char const *filename, int x, int y, int comp, const void *data) -{ - stbi__write_context s = { 0 }; - if (stbi__start_write_file(&s,filename)) { - int r = stbi_write_tga_core(&s, x, y, comp, (void *) data); - stbi__end_write_file(&s); - return r; - } else - return 0; +STBIWDEF int stbi_write_tga(const char* filename, int x, int y, int comp, const void* data) { + stbi__write_context s = {0}; + if (stbi__start_write_file(&s, filename)) { + int r = stbi_write_tga_core(&s, x, y, comp, (void*)data); + stbi__end_write_file(&s); + return r; + } else + return 0; } #endif @@ -620,978 +719,1234 @@ STBIWDEF int stbi_write_tga(char const *filename, int x, int y, int comp, const // Radiance RGBE HDR writer // by Baldur Karlsson -#define stbiw__max(a, b) ((a) > (b) ? (a) : (b)) +#define stbiw__max(a, b) ((a) > (b) ? (a) : (b)) -static void stbiw__linear_to_rgbe(unsigned char *rgbe, float *linear) -{ - int exponent; - float maxcomp = stbiw__max(linear[0], stbiw__max(linear[1], linear[2])); +static void stbiw__linear_to_rgbe(unsigned char* rgbe, float* linear) { + int exponent; + float maxcomp = stbiw__max(linear[0], stbiw__max(linear[1], linear[2])); - if (maxcomp < 1e-32f) { - rgbe[0] = rgbe[1] = rgbe[2] = rgbe[3] = 0; - } else { - float normalize = (float) frexp(maxcomp, &exponent) * 256.0f/maxcomp; + if (maxcomp < 1e-32f) { + rgbe[0] = rgbe[1] = rgbe[2] = rgbe[3] = 0; + } else { + float normalize = (float)frexp(maxcomp, &exponent) * 256.0f / maxcomp; - rgbe[0] = (unsigned char)(linear[0] * normalize); - rgbe[1] = (unsigned char)(linear[1] * normalize); - rgbe[2] = (unsigned char)(linear[2] * normalize); - rgbe[3] = (unsigned char)(exponent + 128); - } + rgbe[0] = (unsigned char)(linear[0] * normalize); + rgbe[1] = (unsigned char)(linear[1] * normalize); + rgbe[2] = (unsigned char)(linear[2] * normalize); + rgbe[3] = (unsigned char)(exponent + 128); + } } -static void stbiw__write_run_data(stbi__write_context *s, int length, unsigned char databyte) -{ - unsigned char lengthbyte = STBIW_UCHAR(length+128); - STBIW_ASSERT(length+128 <= 255); - s->func(s->context, &lengthbyte, 1); - s->func(s->context, &databyte, 1); +static void stbiw__write_run_data(stbi__write_context* s, int length, unsigned char databyte) { + unsigned char lengthbyte = STBIW_UCHAR(length + 128); + STBIW_ASSERT(length + 128 <= 255); + s->func(s->context, &lengthbyte, 1); + s->func(s->context, &databyte, 1); } -static void stbiw__write_dump_data(stbi__write_context *s, int length, unsigned char *data) -{ - unsigned char lengthbyte = STBIW_UCHAR(length); - STBIW_ASSERT(length <= 128); // inconsistent with spec but consistent with official code - s->func(s->context, &lengthbyte, 1); - s->func(s->context, data, length); +static void stbiw__write_dump_data(stbi__write_context* s, int length, unsigned char* data) { + unsigned char lengthbyte = STBIW_UCHAR(length); + STBIW_ASSERT(length <= 128); // inconsistent with spec but consistent with official code + s->func(s->context, &lengthbyte, 1); + s->func(s->context, data, length); } -static void stbiw__write_hdr_scanline(stbi__write_context *s, int width, int ncomp, unsigned char *scratch, float *scanline) -{ - unsigned char scanlineheader[4] = { 2, 2, 0, 0 }; - unsigned char rgbe[4]; - float linear[3]; - int x; - - scanlineheader[2] = (width&0xff00)>>8; - scanlineheader[3] = (width&0x00ff); - - /* skip RLE for images too small or large */ - if (width < 8 || width >= 32768) { - for (x=0; x < width; x++) { - switch (ncomp) { - case 4: /* fallthrough */ - case 3: linear[2] = scanline[x*ncomp + 2]; - linear[1] = scanline[x*ncomp + 1]; - linear[0] = scanline[x*ncomp + 0]; - break; - default: - linear[0] = linear[1] = linear[2] = scanline[x*ncomp + 0]; - break; - } - stbiw__linear_to_rgbe(rgbe, linear); - s->func(s->context, rgbe, 4); +static void stbiw__write_hdr_scanline( + stbi__write_context* s, int width, int ncomp, unsigned char* scratch, float* scanline) { + unsigned char scanlineheader[4] = {2, 2, 0, 0}; + unsigned char rgbe[4]; + float linear[3]; + int x; + + scanlineheader[2] = (width & 0xff00) >> 8; + scanlineheader[3] = (width & 0x00ff); + + /* skip RLE for images too small or large */ + if (width < 8 || width >= 32768) { + for (x = 0; x < width; x++) { + switch (ncomp) { + case 4: /* fallthrough */ + case 3: + linear[2] = scanline[x * ncomp + 2]; + linear[1] = scanline[x * ncomp + 1]; + linear[0] = scanline[x * ncomp + 0]; + break; + default: + linear[0] = linear[1] = linear[2] = scanline[x * ncomp + 0]; + break; } - } else { - int c,r; - /* encode into scratch buffer */ - for (x=0; x < width; x++) { - switch(ncomp) { - case 4: /* fallthrough */ - case 3: linear[2] = scanline[x*ncomp + 2]; - linear[1] = scanline[x*ncomp + 1]; - linear[0] = scanline[x*ncomp + 0]; - break; - default: - linear[0] = linear[1] = linear[2] = scanline[x*ncomp + 0]; - break; - } - stbiw__linear_to_rgbe(rgbe, linear); - scratch[x + width*0] = rgbe[0]; - scratch[x + width*1] = rgbe[1]; - scratch[x + width*2] = rgbe[2]; - scratch[x + width*3] = rgbe[3]; + stbiw__linear_to_rgbe(rgbe, linear); + s->func(s->context, rgbe, 4); + } + } else { + int c, r; + /* encode into scratch buffer */ + for (x = 0; x < width; x++) { + switch (ncomp) { + case 4: /* fallthrough */ + case 3: + linear[2] = scanline[x * ncomp + 2]; + linear[1] = scanline[x * ncomp + 1]; + linear[0] = scanline[x * ncomp + 0]; + break; + default: + linear[0] = linear[1] = linear[2] = scanline[x * ncomp + 0]; + break; } - - s->func(s->context, scanlineheader, 4); - - /* RLE each component separately */ - for (c=0; c < 4; c++) { - unsigned char *comp = &scratch[width*c]; - - x = 0; - while (x < width) { - // find first run - r = x; - while (r+2 < width) { - if (comp[r] == comp[r+1] && comp[r] == comp[r+2]) - break; - ++r; - } - if (r+2 >= width) - r = width; - // dump up to first run - while (x < r) { - int len = r-x; - if (len > 128) len = 128; - stbiw__write_dump_data(s, len, &comp[x]); - x += len; - } - // if there's a run, output it - if (r+2 < width) { // same test as what we break out of in search loop, so only true if we break'd - // find next byte after run - while (r < width && comp[r] == comp[x]) - ++r; - // output run up to r - while (x < r) { - int len = r-x; - if (len > 127) len = 127; - stbiw__write_run_data(s, len, comp[x]); - x += len; - } - } - } + stbiw__linear_to_rgbe(rgbe, linear); + scratch[x + width * 0] = rgbe[0]; + scratch[x + width * 1] = rgbe[1]; + scratch[x + width * 2] = rgbe[2]; + scratch[x + width * 3] = rgbe[3]; + } + + s->func(s->context, scanlineheader, 4); + + /* RLE each component separately */ + for (c = 0; c < 4; c++) { + unsigned char* comp = &scratch[width * c]; + + x = 0; + while (x < width) { + // find first run + r = x; + while (r + 2 < width) { + if (comp[r] == comp[r + 1] && comp[r] == comp[r + 2]) + break; + ++r; + } + if (r + 2 >= width) + r = width; + // dump up to first run + while (x < r) { + int len = r - x; + if (len > 128) + len = 128; + stbiw__write_dump_data(s, len, &comp[x]); + x += len; + } + // if there's a run, output it + if (r + 2 < width) { // same test as what we break out of in search + // loop, so only true if we break'd + // find next byte after run + while (r < width && comp[r] == comp[x]) + ++r; + // output run up to r + while (x < r) { + int len = r - x; + if (len > 127) + len = 127; + stbiw__write_run_data(s, len, comp[x]); + x += len; + } + } } - } + } + } } -static int stbi_write_hdr_core(stbi__write_context *s, int x, int y, int comp, float *data) -{ - if (y <= 0 || x <= 0 || data == NULL) - return 0; - else { - // Each component is stored separately. Allocate scratch space for full output scanline. - unsigned char *scratch = (unsigned char *) STBIW_MALLOC(x*4); - int i, len; - char buffer[128]; - char header[] = "#?RADIANCE\n# Written by stb_image_write.h\nFORMAT=32-bit_rle_rgbe\n"; - s->func(s->context, header, sizeof(header)-1); +static int stbi_write_hdr_core(stbi__write_context* s, int x, int y, int comp, float* data) { + if (y <= 0 || x <= 0 || data == NULL) + return 0; + else { + // Each component is stored separately. Allocate scratch space for full + // output scanline. + unsigned char* scratch = (unsigned char*)STBIW_MALLOC(x * 4); + int i, len; + char buffer[128]; + char header[] = "#?RADIANCE\n# Written by stb_image_write.h\nFORMAT=32-bit_rle_rgbe\n"; + s->func(s->context, header, sizeof(header) - 1); #ifdef __STDC_WANT_SECURE_LIB__ - len = sprintf_s(buffer, sizeof(buffer), "EXPOSURE= 1.0000000000000\n\n-Y %d +X %d\n", y, x); + len = sprintf_s( + buffer, sizeof(buffer), "EXPOSURE= 1.0000000000000\n\n-Y %d +X %d\n", y, x); #else - len = sprintf(buffer, "EXPOSURE= 1.0000000000000\n\n-Y %d +X %d\n", y, x); + len = sprintf(buffer, "EXPOSURE= 1.0000000000000\n\n-Y %d +X %d\n", y, x); #endif - s->func(s->context, buffer, len); - - for(i=0; i < y; i++) - stbiw__write_hdr_scanline(s, x, comp, scratch, data + comp*x*(stbi__flip_vertically_on_write ? y-1-i : i)); - STBIW_FREE(scratch); - return 1; - } + s->func(s->context, buffer, len); + + for (i = 0; i < y; i++) + stbiw__write_hdr_scanline( + s, x, comp, scratch, data + comp * x * (stbi__flip_vertically_on_write ? y - 1 - i : i)); + STBIW_FREE(scratch); + return 1; + } } -STBIWDEF int stbi_write_hdr_to_func(stbi_write_func *func, void *context, int x, int y, int comp, const float *data) -{ - stbi__write_context s = { 0 }; - stbi__start_write_callbacks(&s, func, context); - return stbi_write_hdr_core(&s, x, y, comp, (float *) data); +STBIWDEF int stbi_write_hdr_to_func( + stbi_write_func* func, void* context, int x, int y, int comp, const float* data) { + stbi__write_context s = {0}; + stbi__start_write_callbacks(&s, func, context); + return stbi_write_hdr_core(&s, x, y, comp, (float*)data); } #ifndef STBI_WRITE_NO_STDIO -STBIWDEF int stbi_write_hdr(char const *filename, int x, int y, int comp, const float *data) -{ - stbi__write_context s = { 0 }; - if (stbi__start_write_file(&s,filename)) { - int r = stbi_write_hdr_core(&s, x, y, comp, (float *) data); - stbi__end_write_file(&s); - return r; - } else - return 0; +STBIWDEF int stbi_write_hdr(const char* filename, int x, int y, int comp, const float* data) { + stbi__write_context s = {0}; + if (stbi__start_write_file(&s, filename)) { + int r = stbi_write_hdr_core(&s, x, y, comp, (float*)data); + stbi__end_write_file(&s); + return r; + } else + return 0; } #endif // STBI_WRITE_NO_STDIO - ////////////////////////////////////////////////////////////////////////////// // // PNG writer // #ifndef STBIW_ZLIB_COMPRESS -// stretchy buffer; stbiw__sbpush() == vector<>::push_back() -- stbiw__sbcount() == vector<>::size() -#define stbiw__sbraw(a) ((int *) (void *) (a) - 2) -#define stbiw__sbm(a) stbiw__sbraw(a)[0] -#define stbiw__sbn(a) stbiw__sbraw(a)[1] - -#define stbiw__sbneedgrow(a,n) ((a)==0 || stbiw__sbn(a)+n >= stbiw__sbm(a)) -#define stbiw__sbmaybegrow(a,n) (stbiw__sbneedgrow(a,(n)) ? stbiw__sbgrow(a,n) : 0) -#define stbiw__sbgrow(a,n) stbiw__sbgrowf((void **) &(a), (n), sizeof(*(a))) - -#define stbiw__sbpush(a, v) (stbiw__sbmaybegrow(a,1), (a)[stbiw__sbn(a)++] = (v)) -#define stbiw__sbcount(a) ((a) ? stbiw__sbn(a) : 0) -#define stbiw__sbfree(a) ((a) ? STBIW_FREE(stbiw__sbraw(a)),0 : 0) - -static void *stbiw__sbgrowf(void **arr, int increment, int itemsize) -{ - int m = *arr ? 2*stbiw__sbm(*arr)+increment : increment+1; - void *p = STBIW_REALLOC_SIZED(*arr ? stbiw__sbraw(*arr) : 0, *arr ? (stbiw__sbm(*arr)*itemsize + sizeof(int)*2) : 0, itemsize * m + sizeof(int)*2); - STBIW_ASSERT(p); - if (p) { - if (!*arr) ((int *) p)[1] = 0; - *arr = (void *) ((int *) p + 2); - stbiw__sbm(*arr) = m; - } - return *arr; +// stretchy buffer; stbiw__sbpush() == vector<>::push_back() -- stbiw__sbcount() +// == vector<>::size() +#define stbiw__sbraw(a) ((int*)(void*)(a) - 2) +#define stbiw__sbm(a) stbiw__sbraw(a)[0] +#define stbiw__sbn(a) stbiw__sbraw(a)[1] + +#define stbiw__sbneedgrow(a, n) ((a) == 0 || stbiw__sbn(a) + n >= stbiw__sbm(a)) +#define stbiw__sbmaybegrow(a, n) (stbiw__sbneedgrow(a, (n)) ? stbiw__sbgrow(a, n) : 0) +#define stbiw__sbgrow(a, n) stbiw__sbgrowf((void**)&(a), (n), sizeof(*(a))) + +#define stbiw__sbpush(a, v) (stbiw__sbmaybegrow(a, 1), (a)[stbiw__sbn(a)++] = (v)) +#define stbiw__sbcount(a) ((a) ? stbiw__sbn(a) : 0) +#define stbiw__sbfree(a) ((a) ? STBIW_FREE(stbiw__sbraw(a)), 0 : 0) + +static void* stbiw__sbgrowf(void** arr, int increment, int itemsize) { + int m = *arr ? 2 * stbiw__sbm(*arr) + increment : increment + 1; + void* p = STBIW_REALLOC_SIZED(*arr ? stbiw__sbraw(*arr) : 0, + *arr ? (stbiw__sbm(*arr) * itemsize + sizeof(int) * 2) : 0, + itemsize * m + sizeof(int) * 2); + STBIW_ASSERT(p); + if (p) { + if (!*arr) + ((int*)p)[1] = 0; + *arr = (void*)((int*)p + 2); + stbiw__sbm(*arr) = m; + } + return *arr; } -static unsigned char *stbiw__zlib_flushf(unsigned char *data, unsigned int *bitbuffer, int *bitcount) -{ - while (*bitcount >= 8) { - stbiw__sbpush(data, STBIW_UCHAR(*bitbuffer)); - *bitbuffer >>= 8; - *bitcount -= 8; - } - return data; +static unsigned char* + stbiw__zlib_flushf(unsigned char* data, unsigned int* bitbuffer, int* bitcount) { + while (*bitcount >= 8) { + stbiw__sbpush(data, STBIW_UCHAR(*bitbuffer)); + *bitbuffer >>= 8; + *bitcount -= 8; + } + return data; } -static int stbiw__zlib_bitrev(int code, int codebits) -{ - int res=0; - while (codebits--) { - res = (res << 1) | (code & 1); - code >>= 1; - } - return res; +static int stbiw__zlib_bitrev(int code, int codebits) { + int res = 0; + while (codebits--) { + res = (res << 1) | (code & 1); + code >>= 1; + } + return res; } -static unsigned int stbiw__zlib_countm(unsigned char *a, unsigned char *b, int limit) -{ - int i; - for (i=0; i < limit && i < 258; ++i) - if (a[i] != b[i]) break; - return i; +static unsigned int stbiw__zlib_countm(unsigned char* a, unsigned char* b, int limit) { + int i; + for (i = 0; i < limit && i < 258; ++i) + if (a[i] != b[i]) + break; + return i; } -static unsigned int stbiw__zhash(unsigned char *data) -{ - stbiw_uint32 hash = data[0] + (data[1] << 8) + (data[2] << 16); - hash ^= hash << 3; - hash += hash >> 5; - hash ^= hash << 4; - hash += hash >> 17; - hash ^= hash << 25; - hash += hash >> 6; - return hash; +static unsigned int stbiw__zhash(unsigned char* data) { + stbiw_uint32 hash = data[0] + (data[1] << 8) + (data[2] << 16); + hash ^= hash << 3; + hash += hash >> 5; + hash ^= hash << 4; + hash += hash >> 17; + hash ^= hash << 25; + hash += hash >> 6; + return hash; } #define stbiw__zlib_flush() (out = stbiw__zlib_flushf(out, &bitbuf, &bitcount)) -#define stbiw__zlib_add(code,codebits) \ - (bitbuf |= (code) << bitcount, bitcount += (codebits), stbiw__zlib_flush()) -#define stbiw__zlib_huffa(b,c) stbiw__zlib_add(stbiw__zlib_bitrev(b,c),c) +#define stbiw__zlib_add(code, codebits) \ + (bitbuf |= (code) << bitcount, bitcount += (codebits), stbiw__zlib_flush()) +#define stbiw__zlib_huffa(b, c) stbiw__zlib_add(stbiw__zlib_bitrev(b, c), c) // default huffman tables -#define stbiw__zlib_huff1(n) stbiw__zlib_huffa(0x30 + (n), 8) -#define stbiw__zlib_huff2(n) stbiw__zlib_huffa(0x190 + (n)-144, 9) -#define stbiw__zlib_huff3(n) stbiw__zlib_huffa(0 + (n)-256,7) -#define stbiw__zlib_huff4(n) stbiw__zlib_huffa(0xc0 + (n)-280,8) -#define stbiw__zlib_huff(n) ((n) <= 143 ? stbiw__zlib_huff1(n) : (n) <= 255 ? stbiw__zlib_huff2(n) : (n) <= 279 ? stbiw__zlib_huff3(n) : stbiw__zlib_huff4(n)) +#define stbiw__zlib_huff1(n) stbiw__zlib_huffa(0x30 + (n), 8) +#define stbiw__zlib_huff2(n) stbiw__zlib_huffa(0x190 + (n) - 144, 9) +#define stbiw__zlib_huff3(n) stbiw__zlib_huffa(0 + (n) - 256, 7) +#define stbiw__zlib_huff4(n) stbiw__zlib_huffa(0xc0 + (n) - 280, 8) +#define stbiw__zlib_huff(n) \ + ((n) <= 143 ? stbiw__zlib_huff1(n) \ + : (n) <= 255 ? stbiw__zlib_huff2(n) \ + : (n) <= 279 ? stbiw__zlib_huff3(n) \ + : stbiw__zlib_huff4(n)) #define stbiw__zlib_huffb(n) ((n) <= 143 ? stbiw__zlib_huff1(n) : stbiw__zlib_huff2(n)) -#define stbiw__ZHASH 16384 +#define stbiw__ZHASH 16384 #endif // STBIW_ZLIB_COMPRESS -STBIWDEF unsigned char * stbi_zlib_compress(unsigned char *data, int data_len, int *out_len, int quality) -{ +STBIWDEF unsigned char* + stbi_zlib_compress(unsigned char* data, int data_len, int* out_len, int quality) { #ifdef STBIW_ZLIB_COMPRESS - // user provided a zlib compress implementation, use that - return STBIW_ZLIB_COMPRESS(data, data_len, out_len, quality); -#else // use builtin - static unsigned short lengthc[] = { 3,4,5,6,7,8,9,10,11,13,15,17,19,23,27,31,35,43,51,59,67,83,99,115,131,163,195,227,258, 259 }; - static unsigned char lengtheb[]= { 0,0,0,0,0,0,0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3, 4, 4, 4, 4, 5, 5, 5, 5, 0 }; - static unsigned short distc[] = { 1,2,3,4,5,7,9,13,17,25,33,49,65,97,129,193,257,385,513,769,1025,1537,2049,3073,4097,6145,8193,12289,16385,24577, 32768 }; - static unsigned char disteb[] = { 0,0,0,0,1,1,2,2,3,3,4,4,5,5,6,6,7,7,8,8,9,9,10,10,11,11,12,12,13,13 }; - unsigned int bitbuf=0; - int i,j, bitcount=0; - unsigned char *out = NULL; - unsigned char ***hash_table = (unsigned char***) STBIW_MALLOC(stbiw__ZHASH * sizeof(unsigned char**)); - if (hash_table == NULL) - return NULL; - if (quality < 5) quality = 5; - - stbiw__sbpush(out, 0x78); // DEFLATE 32K window - stbiw__sbpush(out, 0x5e); // FLEVEL = 1 - stbiw__zlib_add(1,1); // BFINAL = 1 - stbiw__zlib_add(1,2); // BTYPE = 1 -- fixed huffman - - for (i=0; i < stbiw__ZHASH; ++i) - hash_table[i] = NULL; - - i=0; - while (i < data_len-3) { - // hash next 3 bytes of data to be compressed - int h = stbiw__zhash(data+i)&(stbiw__ZHASH-1), best=3; - unsigned char *bestloc = 0; - unsigned char **hlist = hash_table[h]; - int n = stbiw__sbcount(hlist); - for (j=0; j < n; ++j) { - if (hlist[j]-data > i-32768) { // if entry lies within window - int d = stbiw__zlib_countm(hlist[j], data+i, data_len-i); - if (d >= best) { best=d; bestloc=hlist[j]; } - } + // user provided a zlib compress implementation, use that + return STBIW_ZLIB_COMPRESS(data, data_len, out_len, quality); +#else // use builtin + static unsigned short lengthc[] = {3, 4, 5, 6, 7, 8, 9, 10, 11, 13, + 15, 17, 19, 23, 27, 31, 35, 43, 51, 59, + 67, 83, 99, 115, 131, 163, 195, 227, 258, 259}; + static unsigned char lengtheb[] = {0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, + 2, 3, 3, 3, 3, 4, 4, 4, 4, 5, 5, 5, 5, 0}; + static unsigned short distc[] = { + 1, 2, 3, 4, 5, 7, 9, 13, 17, 25, 33, 49, 65, 97, 129, 193, + 257, 385, 513, 769, 1025, 1537, 2049, 3073, 4097, 6145, 8193, 12289, 16385, 24577, 32768}; + static unsigned char disteb[] = {0, 0, 0, 0, 1, 1, 2, 2, 3, 3, 4, 4, 5, 5, 6, + 6, 7, 7, 8, 8, 9, 9, 10, 10, 11, 11, 12, 12, 13, 13}; + unsigned int bitbuf = 0; + int i, j, bitcount = 0; + unsigned char* out = NULL; + unsigned char*** hash_table = + (unsigned char***)STBIW_MALLOC(stbiw__ZHASH * sizeof(unsigned char**)); + if (hash_table == NULL) + return NULL; + if (quality < 5) + quality = 5; + + stbiw__sbpush(out, 0x78); // DEFLATE 32K window + stbiw__sbpush(out, 0x5e); // FLEVEL = 1 + stbiw__zlib_add(1, 1); // BFINAL = 1 + stbiw__zlib_add(1, 2); // BTYPE = 1 -- fixed huffman + + for (i = 0; i < stbiw__ZHASH; ++i) + hash_table[i] = NULL; + + i = 0; + while (i < data_len - 3) { + // hash next 3 bytes of data to be compressed + int h = stbiw__zhash(data + i) & (stbiw__ZHASH - 1), best = 3; + unsigned char* bestloc = 0; + unsigned char** hlist = hash_table[h]; + int n = stbiw__sbcount(hlist); + for (j = 0; j < n; ++j) { + if (hlist[j] - data > i - 32768) { // if entry lies within window + int d = stbiw__zlib_countm(hlist[j], data + i, data_len - i); + if (d >= best) { + best = d; + bestloc = hlist[j]; + } } - // when hash table entry is too long, delete half the entries - if (hash_table[h] && stbiw__sbn(hash_table[h]) == 2*quality) { - STBIW_MEMMOVE(hash_table[h], hash_table[h]+quality, sizeof(hash_table[h][0])*quality); - stbiw__sbn(hash_table[h]) = quality; - } - stbiw__sbpush(hash_table[h],data+i); - - if (bestloc) { - // "lazy matching" - check match at *next* byte, and if it's better, do cur byte as literal - h = stbiw__zhash(data+i+1)&(stbiw__ZHASH-1); - hlist = hash_table[h]; - n = stbiw__sbcount(hlist); - for (j=0; j < n; ++j) { - if (hlist[j]-data > i-32767) { - int e = stbiw__zlib_countm(hlist[j], data+i+1, data_len-i-1); - if (e > best) { // if next match is better, bail on current match - bestloc = NULL; - break; - } - } - } - } - - if (bestloc) { - int d = (int) (data+i - bestloc); // distance back - STBIW_ASSERT(d <= 32767 && best <= 258); - for (j=0; best > lengthc[j+1]-1; ++j); - stbiw__zlib_huff(j+257); - if (lengtheb[j]) stbiw__zlib_add(best - lengthc[j], lengtheb[j]); - for (j=0; d > distc[j+1]-1; ++j); - stbiw__zlib_add(stbiw__zlib_bitrev(j,5),5); - if (disteb[j]) stbiw__zlib_add(d - distc[j], disteb[j]); - i += best; - } else { - stbiw__zlib_huffb(data[i]); - ++i; + } + // when hash table entry is too long, delete half the entries + if (hash_table[h] && stbiw__sbn(hash_table[h]) == 2 * quality) { + STBIW_MEMMOVE(hash_table[h], hash_table[h] + quality, sizeof(hash_table[h][0]) * quality); + stbiw__sbn(hash_table[h]) = quality; + } + stbiw__sbpush(hash_table[h], data + i); + + if (bestloc) { + // "lazy matching" - check match at *next* byte, and if it's better, do + // cur byte as literal + h = stbiw__zhash(data + i + 1) & (stbiw__ZHASH - 1); + hlist = hash_table[h]; + n = stbiw__sbcount(hlist); + for (j = 0; j < n; ++j) { + if (hlist[j] - data > i - 32767) { + int e = stbiw__zlib_countm(hlist[j], data + i + 1, data_len - i - 1); + if (e > best) { // if next match is better, bail on current match + bestloc = NULL; + break; + } + } } - } - // write out final bytes - for (;i < data_len; ++i) + } + + if (bestloc) { + int d = (int)(data + i - bestloc); // distance back + STBIW_ASSERT(d <= 32767 && best <= 258); + for (j = 0; best > lengthc[j + 1] - 1; ++j) + ; + stbiw__zlib_huff(j + 257); + if (lengtheb[j]) + stbiw__zlib_add(best - lengthc[j], lengtheb[j]); + for (j = 0; d > distc[j + 1] - 1; ++j) + ; + stbiw__zlib_add(stbiw__zlib_bitrev(j, 5), 5); + if (disteb[j]) + stbiw__zlib_add(d - distc[j], disteb[j]); + i += best; + } else { stbiw__zlib_huffb(data[i]); - stbiw__zlib_huff(256); // end of block - // pad with 0 bits to byte boundary - while (bitcount) - stbiw__zlib_add(0,1); - - for (i=0; i < stbiw__ZHASH; ++i) - (void) stbiw__sbfree(hash_table[i]); - STBIW_FREE(hash_table); - - { - // compute adler32 on input - unsigned int s1=1, s2=0; - int blocklen = (int) (data_len % 5552); - j=0; - while (j < data_len) { - for (i=0; i < blocklen; ++i) { s1 += data[j+i]; s2 += s1; } - s1 %= 65521; s2 %= 65521; - j += blocklen; - blocklen = 5552; + ++i; + } + } + // write out final bytes + for (; i < data_len; ++i) + stbiw__zlib_huffb(data[i]); + stbiw__zlib_huff(256); // end of block + // pad with 0 bits to byte boundary + while (bitcount) + stbiw__zlib_add(0, 1); + + for (i = 0; i < stbiw__ZHASH; ++i) + (void)stbiw__sbfree(hash_table[i]); + STBIW_FREE(hash_table); + + { + // compute adler32 on input + unsigned int s1 = 1, s2 = 0; + int blocklen = (int)(data_len % 5552); + j = 0; + while (j < data_len) { + for (i = 0; i < blocklen; ++i) { + s1 += data[j + i]; + s2 += s1; } - stbiw__sbpush(out, STBIW_UCHAR(s2 >> 8)); - stbiw__sbpush(out, STBIW_UCHAR(s2)); - stbiw__sbpush(out, STBIW_UCHAR(s1 >> 8)); - stbiw__sbpush(out, STBIW_UCHAR(s1)); - } - *out_len = stbiw__sbn(out); - // make returned pointer freeable - STBIW_MEMMOVE(stbiw__sbraw(out), out, *out_len); - return (unsigned char *) stbiw__sbraw(out); + s1 %= 65521; + s2 %= 65521; + j += blocklen; + blocklen = 5552; + } + stbiw__sbpush(out, STBIW_UCHAR(s2 >> 8)); + stbiw__sbpush(out, STBIW_UCHAR(s2)); + stbiw__sbpush(out, STBIW_UCHAR(s1 >> 8)); + stbiw__sbpush(out, STBIW_UCHAR(s1)); + } + *out_len = stbiw__sbn(out); + // make returned pointer freeable + STBIW_MEMMOVE(stbiw__sbraw(out), out, *out_len); + return (unsigned char*)stbiw__sbraw(out); #endif // STBIW_ZLIB_COMPRESS } -static unsigned int stbiw__crc32(unsigned char *buffer, int len) -{ +static unsigned int stbiw__crc32(unsigned char* buffer, int len) { #ifdef STBIW_CRC32 - return STBIW_CRC32(buffer, len); + return STBIW_CRC32(buffer, len); #else - static unsigned int crc_table[256] = - { - 0x00000000, 0x77073096, 0xEE0E612C, 0x990951BA, 0x076DC419, 0x706AF48F, 0xE963A535, 0x9E6495A3, - 0x0eDB8832, 0x79DCB8A4, 0xE0D5E91E, 0x97D2D988, 0x09B64C2B, 0x7EB17CBD, 0xE7B82D07, 0x90BF1D91, - 0x1DB71064, 0x6AB020F2, 0xF3B97148, 0x84BE41DE, 0x1ADAD47D, 0x6DDDE4EB, 0xF4D4B551, 0x83D385C7, - 0x136C9856, 0x646BA8C0, 0xFD62F97A, 0x8A65C9EC, 0x14015C4F, 0x63066CD9, 0xFA0F3D63, 0x8D080DF5, - 0x3B6E20C8, 0x4C69105E, 0xD56041E4, 0xA2677172, 0x3C03E4D1, 0x4B04D447, 0xD20D85FD, 0xA50AB56B, - 0x35B5A8FA, 0x42B2986C, 0xDBBBC9D6, 0xACBCF940, 0x32D86CE3, 0x45DF5C75, 0xDCD60DCF, 0xABD13D59, - 0x26D930AC, 0x51DE003A, 0xC8D75180, 0xBFD06116, 0x21B4F4B5, 0x56B3C423, 0xCFBA9599, 0xB8BDA50F, - 0x2802B89E, 0x5F058808, 0xC60CD9B2, 0xB10BE924, 0x2F6F7C87, 0x58684C11, 0xC1611DAB, 0xB6662D3D, - 0x76DC4190, 0x01DB7106, 0x98D220BC, 0xEFD5102A, 0x71B18589, 0x06B6B51F, 0x9FBFE4A5, 0xE8B8D433, - 0x7807C9A2, 0x0F00F934, 0x9609A88E, 0xE10E9818, 0x7F6A0DBB, 0x086D3D2D, 0x91646C97, 0xE6635C01, - 0x6B6B51F4, 0x1C6C6162, 0x856530D8, 0xF262004E, 0x6C0695ED, 0x1B01A57B, 0x8208F4C1, 0xF50FC457, - 0x65B0D9C6, 0x12B7E950, 0x8BBEB8EA, 0xFCB9887C, 0x62DD1DDF, 0x15DA2D49, 0x8CD37CF3, 0xFBD44C65, - 0x4DB26158, 0x3AB551CE, 0xA3BC0074, 0xD4BB30E2, 0x4ADFA541, 0x3DD895D7, 0xA4D1C46D, 0xD3D6F4FB, - 0x4369E96A, 0x346ED9FC, 0xAD678846, 0xDA60B8D0, 0x44042D73, 0x33031DE5, 0xAA0A4C5F, 0xDD0D7CC9, - 0x5005713C, 0x270241AA, 0xBE0B1010, 0xC90C2086, 0x5768B525, 0x206F85B3, 0xB966D409, 0xCE61E49F, - 0x5EDEF90E, 0x29D9C998, 0xB0D09822, 0xC7D7A8B4, 0x59B33D17, 0x2EB40D81, 0xB7BD5C3B, 0xC0BA6CAD, - 0xEDB88320, 0x9ABFB3B6, 0x03B6E20C, 0x74B1D29A, 0xEAD54739, 0x9DD277AF, 0x04DB2615, 0x73DC1683, - 0xE3630B12, 0x94643B84, 0x0D6D6A3E, 0x7A6A5AA8, 0xE40ECF0B, 0x9309FF9D, 0x0A00AE27, 0x7D079EB1, - 0xF00F9344, 0x8708A3D2, 0x1E01F268, 0x6906C2FE, 0xF762575D, 0x806567CB, 0x196C3671, 0x6E6B06E7, - 0xFED41B76, 0x89D32BE0, 0x10DA7A5A, 0x67DD4ACC, 0xF9B9DF6F, 0x8EBEEFF9, 0x17B7BE43, 0x60B08ED5, - 0xD6D6A3E8, 0xA1D1937E, 0x38D8C2C4, 0x4FDFF252, 0xD1BB67F1, 0xA6BC5767, 0x3FB506DD, 0x48B2364B, - 0xD80D2BDA, 0xAF0A1B4C, 0x36034AF6, 0x41047A60, 0xDF60EFC3, 0xA867DF55, 0x316E8EEF, 0x4669BE79, - 0xCB61B38C, 0xBC66831A, 0x256FD2A0, 0x5268E236, 0xCC0C7795, 0xBB0B4703, 0x220216B9, 0x5505262F, - 0xC5BA3BBE, 0xB2BD0B28, 0x2BB45A92, 0x5CB36A04, 0xC2D7FFA7, 0xB5D0CF31, 0x2CD99E8B, 0x5BDEAE1D, - 0x9B64C2B0, 0xEC63F226, 0x756AA39C, 0x026D930A, 0x9C0906A9, 0xEB0E363F, 0x72076785, 0x05005713, - 0x95BF4A82, 0xE2B87A14, 0x7BB12BAE, 0x0CB61B38, 0x92D28E9B, 0xE5D5BE0D, 0x7CDCEFB7, 0x0BDBDF21, - 0x86D3D2D4, 0xF1D4E242, 0x68DDB3F8, 0x1FDA836E, 0x81BE16CD, 0xF6B9265B, 0x6FB077E1, 0x18B74777, - 0x88085AE6, 0xFF0F6A70, 0x66063BCA, 0x11010B5C, 0x8F659EFF, 0xF862AE69, 0x616BFFD3, 0x166CCF45, - 0xA00AE278, 0xD70DD2EE, 0x4E048354, 0x3903B3C2, 0xA7672661, 0xD06016F7, 0x4969474D, 0x3E6E77DB, - 0xAED16A4A, 0xD9D65ADC, 0x40DF0B66, 0x37D83BF0, 0xA9BCAE53, 0xDEBB9EC5, 0x47B2CF7F, 0x30B5FFE9, - 0xBDBDF21C, 0xCABAC28A, 0x53B39330, 0x24B4A3A6, 0xBAD03605, 0xCDD70693, 0x54DE5729, 0x23D967BF, - 0xB3667A2E, 0xC4614AB8, 0x5D681B02, 0x2A6F2B94, 0xB40BBE37, 0xC30C8EA1, 0x5A05DF1B, 0x2D02EF8D - }; - - unsigned int crc = ~0u; - int i; - for (i=0; i < len; ++i) - crc = (crc >> 8) ^ crc_table[buffer[i] ^ (crc & 0xff)]; - return ~crc; + static unsigned int crc_table[256] = { + 0x00000000, 0x77073096, 0xEE0E612C, 0x990951BA, 0x076DC419, 0x706AF48F, 0xE963A535, + 0x9E6495A3, 0x0eDB8832, 0x79DCB8A4, 0xE0D5E91E, 0x97D2D988, 0x09B64C2B, 0x7EB17CBD, + 0xE7B82D07, 0x90BF1D91, 0x1DB71064, 0x6AB020F2, 0xF3B97148, 0x84BE41DE, 0x1ADAD47D, + 0x6DDDE4EB, 0xF4D4B551, 0x83D385C7, 0x136C9856, 0x646BA8C0, 0xFD62F97A, 0x8A65C9EC, + 0x14015C4F, 0x63066CD9, 0xFA0F3D63, 0x8D080DF5, 0x3B6E20C8, 0x4C69105E, 0xD56041E4, + 0xA2677172, 0x3C03E4D1, 0x4B04D447, 0xD20D85FD, 0xA50AB56B, 0x35B5A8FA, 0x42B2986C, + 0xDBBBC9D6, 0xACBCF940, 0x32D86CE3, 0x45DF5C75, 0xDCD60DCF, 0xABD13D59, 0x26D930AC, + 0x51DE003A, 0xC8D75180, 0xBFD06116, 0x21B4F4B5, 0x56B3C423, 0xCFBA9599, 0xB8BDA50F, + 0x2802B89E, 0x5F058808, 0xC60CD9B2, 0xB10BE924, 0x2F6F7C87, 0x58684C11, 0xC1611DAB, + 0xB6662D3D, 0x76DC4190, 0x01DB7106, 0x98D220BC, 0xEFD5102A, 0x71B18589, 0x06B6B51F, + 0x9FBFE4A5, 0xE8B8D433, 0x7807C9A2, 0x0F00F934, 0x9609A88E, 0xE10E9818, 0x7F6A0DBB, + 0x086D3D2D, 0x91646C97, 0xE6635C01, 0x6B6B51F4, 0x1C6C6162, 0x856530D8, 0xF262004E, + 0x6C0695ED, 0x1B01A57B, 0x8208F4C1, 0xF50FC457, 0x65B0D9C6, 0x12B7E950, 0x8BBEB8EA, + 0xFCB9887C, 0x62DD1DDF, 0x15DA2D49, 0x8CD37CF3, 0xFBD44C65, 0x4DB26158, 0x3AB551CE, + 0xA3BC0074, 0xD4BB30E2, 0x4ADFA541, 0x3DD895D7, 0xA4D1C46D, 0xD3D6F4FB, 0x4369E96A, + 0x346ED9FC, 0xAD678846, 0xDA60B8D0, 0x44042D73, 0x33031DE5, 0xAA0A4C5F, 0xDD0D7CC9, + 0x5005713C, 0x270241AA, 0xBE0B1010, 0xC90C2086, 0x5768B525, 0x206F85B3, 0xB966D409, + 0xCE61E49F, 0x5EDEF90E, 0x29D9C998, 0xB0D09822, 0xC7D7A8B4, 0x59B33D17, 0x2EB40D81, + 0xB7BD5C3B, 0xC0BA6CAD, 0xEDB88320, 0x9ABFB3B6, 0x03B6E20C, 0x74B1D29A, 0xEAD54739, + 0x9DD277AF, 0x04DB2615, 0x73DC1683, 0xE3630B12, 0x94643B84, 0x0D6D6A3E, 0x7A6A5AA8, + 0xE40ECF0B, 0x9309FF9D, 0x0A00AE27, 0x7D079EB1, 0xF00F9344, 0x8708A3D2, 0x1E01F268, + 0x6906C2FE, 0xF762575D, 0x806567CB, 0x196C3671, 0x6E6B06E7, 0xFED41B76, 0x89D32BE0, + 0x10DA7A5A, 0x67DD4ACC, 0xF9B9DF6F, 0x8EBEEFF9, 0x17B7BE43, 0x60B08ED5, 0xD6D6A3E8, + 0xA1D1937E, 0x38D8C2C4, 0x4FDFF252, 0xD1BB67F1, 0xA6BC5767, 0x3FB506DD, 0x48B2364B, + 0xD80D2BDA, 0xAF0A1B4C, 0x36034AF6, 0x41047A60, 0xDF60EFC3, 0xA867DF55, 0x316E8EEF, + 0x4669BE79, 0xCB61B38C, 0xBC66831A, 0x256FD2A0, 0x5268E236, 0xCC0C7795, 0xBB0B4703, + 0x220216B9, 0x5505262F, 0xC5BA3BBE, 0xB2BD0B28, 0x2BB45A92, 0x5CB36A04, 0xC2D7FFA7, + 0xB5D0CF31, 0x2CD99E8B, 0x5BDEAE1D, 0x9B64C2B0, 0xEC63F226, 0x756AA39C, 0x026D930A, + 0x9C0906A9, 0xEB0E363F, 0x72076785, 0x05005713, 0x95BF4A82, 0xE2B87A14, 0x7BB12BAE, + 0x0CB61B38, 0x92D28E9B, 0xE5D5BE0D, 0x7CDCEFB7, 0x0BDBDF21, 0x86D3D2D4, 0xF1D4E242, + 0x68DDB3F8, 0x1FDA836E, 0x81BE16CD, 0xF6B9265B, 0x6FB077E1, 0x18B74777, 0x88085AE6, + 0xFF0F6A70, 0x66063BCA, 0x11010B5C, 0x8F659EFF, 0xF862AE69, 0x616BFFD3, 0x166CCF45, + 0xA00AE278, 0xD70DD2EE, 0x4E048354, 0x3903B3C2, 0xA7672661, 0xD06016F7, 0x4969474D, + 0x3E6E77DB, 0xAED16A4A, 0xD9D65ADC, 0x40DF0B66, 0x37D83BF0, 0xA9BCAE53, 0xDEBB9EC5, + 0x47B2CF7F, 0x30B5FFE9, 0xBDBDF21C, 0xCABAC28A, 0x53B39330, 0x24B4A3A6, 0xBAD03605, + 0xCDD70693, 0x54DE5729, 0x23D967BF, 0xB3667A2E, 0xC4614AB8, 0x5D681B02, 0x2A6F2B94, + 0xB40BBE37, 0xC30C8EA1, 0x5A05DF1B, 0x2D02EF8D}; + + unsigned int crc = ~0u; + int i; + for (i = 0; i < len; ++i) + crc = (crc >> 8) ^ crc_table[buffer[i] ^ (crc & 0xff)]; + return ~crc; #endif } -#define stbiw__wpng4(o,a,b,c,d) ((o)[0]=STBIW_UCHAR(a),(o)[1]=STBIW_UCHAR(b),(o)[2]=STBIW_UCHAR(c),(o)[3]=STBIW_UCHAR(d),(o)+=4) -#define stbiw__wp32(data,v) stbiw__wpng4(data, (v)>>24,(v)>>16,(v)>>8,(v)); -#define stbiw__wptag(data,s) stbiw__wpng4(data, s[0],s[1],s[2],s[3]) - -static void stbiw__wpcrc(unsigned char **data, int len) -{ - unsigned int crc = stbiw__crc32(*data - len - 4, len+4); - stbiw__wp32(*data, crc); +#define stbiw__wpng4(o, a, b, c, d) \ + ((o)[0] = STBIW_UCHAR(a), \ + (o)[1] = STBIW_UCHAR(b), \ + (o)[2] = STBIW_UCHAR(c), \ + (o)[3] = STBIW_UCHAR(d), \ + (o) += 4) +#define stbiw__wp32(data, v) stbiw__wpng4(data, (v) >> 24, (v) >> 16, (v) >> 8, (v)); +#define stbiw__wptag(data, s) stbiw__wpng4(data, s[0], s[1], s[2], s[3]) + +static void stbiw__wpcrc(unsigned char** data, int len) { + unsigned int crc = stbiw__crc32(*data - len - 4, len + 4); + stbiw__wp32(*data, crc); } -static unsigned char stbiw__paeth(int a, int b, int c) -{ - int p = a + b - c, pa = abs(p-a), pb = abs(p-b), pc = abs(p-c); - if (pa <= pb && pa <= pc) return STBIW_UCHAR(a); - if (pb <= pc) return STBIW_UCHAR(b); - return STBIW_UCHAR(c); +static unsigned char stbiw__paeth(int a, int b, int c) { + int p = a + b - c, pa = abs(p - a), pb = abs(p - b), pc = abs(p - c); + if (pa <= pb && pa <= pc) + return STBIW_UCHAR(a); + if (pb <= pc) + return STBIW_UCHAR(b); + return STBIW_UCHAR(c); } // @OPTIMIZE: provide an option that always forces left-predict or paeth predict -static void stbiw__encode_png_line(unsigned char *pixels, int stride_bytes, int width, int height, int y, int n, int filter_type, signed char *line_buffer) -{ - static int mapping[] = { 0,1,2,3,4 }; - static int firstmap[] = { 0,1,0,5,6 }; - int *mymap = (y != 0) ? mapping : firstmap; - int i; - int type = mymap[filter_type]; - unsigned char *z = pixels + stride_bytes * (stbi__flip_vertically_on_write ? height-1-y : y); - int signed_stride = stbi__flip_vertically_on_write ? -stride_bytes : stride_bytes; - - if (type==0) { - memcpy(line_buffer, z, width*n); - return; - } - - // first loop isn't optimized since it's just one pixel - for (i = 0; i < n; ++i) { - switch (type) { - case 1: line_buffer[i] = z[i]; break; - case 2: line_buffer[i] = z[i] - z[i-signed_stride]; break; - case 3: line_buffer[i] = z[i] - (z[i-signed_stride]>>1); break; - case 4: line_buffer[i] = (signed char) (z[i] - stbiw__paeth(0,z[i-signed_stride],0)); break; - case 5: line_buffer[i] = z[i]; break; - case 6: line_buffer[i] = z[i]; break; - } - } - switch (type) { - case 1: for (i=n; i < width*n; ++i) line_buffer[i] = z[i] - z[i-n]; break; - case 2: for (i=n; i < width*n; ++i) line_buffer[i] = z[i] - z[i-signed_stride]; break; - case 3: for (i=n; i < width*n; ++i) line_buffer[i] = z[i] - ((z[i-n] + z[i-signed_stride])>>1); break; - case 4: for (i=n; i < width*n; ++i) line_buffer[i] = z[i] - stbiw__paeth(z[i-n], z[i-signed_stride], z[i-signed_stride-n]); break; - case 5: for (i=n; i < width*n; ++i) line_buffer[i] = z[i] - (z[i-n]>>1); break; - case 6: for (i=n; i < width*n; ++i) line_buffer[i] = z[i] - stbiw__paeth(z[i-n], 0,0); break; - } +static void stbiw__encode_png_line(unsigned char* pixels, + int stride_bytes, + int width, + int height, + int y, + int n, + int filter_type, + signed char* line_buffer) { + static int mapping[] = {0, 1, 2, 3, 4}; + static int firstmap[] = {0, 1, 0, 5, 6}; + int* mymap = (y != 0) ? mapping : firstmap; + int i; + int type = mymap[filter_type]; + unsigned char* z = pixels + stride_bytes * (stbi__flip_vertically_on_write ? height - 1 - y : y); + int signed_stride = stbi__flip_vertically_on_write ? -stride_bytes : stride_bytes; + + if (type == 0) { + memcpy(line_buffer, z, width * n); + return; + } + + // first loop isn't optimized since it's just one pixel + for (i = 0; i < n; ++i) { + switch (type) { + case 1: + line_buffer[i] = z[i]; + break; + case 2: + line_buffer[i] = z[i] - z[i - signed_stride]; + break; + case 3: + line_buffer[i] = z[i] - (z[i - signed_stride] >> 1); + break; + case 4: + line_buffer[i] = (signed char)(z[i] - stbiw__paeth(0, z[i - signed_stride], 0)); + break; + case 5: + line_buffer[i] = z[i]; + break; + case 6: + line_buffer[i] = z[i]; + break; + } + } + switch (type) { + case 1: + for (i = n; i < width * n; ++i) + line_buffer[i] = z[i] - z[i - n]; + break; + case 2: + for (i = n; i < width * n; ++i) + line_buffer[i] = z[i] - z[i - signed_stride]; + break; + case 3: + for (i = n; i < width * n; ++i) + line_buffer[i] = z[i] - ((z[i - n] + z[i - signed_stride]) >> 1); + break; + case 4: + for (i = n; i < width * n; ++i) + line_buffer[i] = + z[i] - stbiw__paeth(z[i - n], z[i - signed_stride], z[i - signed_stride - n]); + break; + case 5: + for (i = n; i < width * n; ++i) + line_buffer[i] = z[i] - (z[i - n] >> 1); + break; + case 6: + for (i = n; i < width * n; ++i) + line_buffer[i] = z[i] - stbiw__paeth(z[i - n], 0, 0); + break; + } } -STBIWDEF unsigned char *stbi_write_png_to_mem(const unsigned char *pixels, int stride_bytes, int x, int y, int n, int *out_len) -{ - int force_filter = stbi_write_force_png_filter; - int ctype[5] = { -1, 0, 4, 2, 6 }; - unsigned char sig[8] = { 137,80,78,71,13,10,26,10 }; - unsigned char *out,*o, *filt, *zlib; - signed char *line_buffer; - int j,zlen; - - if (stride_bytes == 0) - stride_bytes = x * n; - - if (force_filter >= 5) { - force_filter = -1; - } - - filt = (unsigned char *) STBIW_MALLOC((x*n+1) * y); if (!filt) return 0; - line_buffer = (signed char *) STBIW_MALLOC(x * n); if (!line_buffer) { STBIW_FREE(filt); return 0; } - for (j=0; j < y; ++j) { - int filter_type; - if (force_filter > -1) { - filter_type = force_filter; - stbiw__encode_png_line((unsigned char*)(pixels), stride_bytes, x, y, j, n, force_filter, line_buffer); - } else { // Estimate the best filter by running through all of them: - int best_filter = 0, best_filter_val = 0x7fffffff, est, i; - for (filter_type = 0; filter_type < 5; filter_type++) { - stbiw__encode_png_line((unsigned char*)(pixels), stride_bytes, x, y, j, n, filter_type, line_buffer); - - // Estimate the entropy of the line using this filter; the less, the better. - est = 0; - for (i = 0; i < x*n; ++i) { - est += abs((signed char) line_buffer[i]); - } - if (est < best_filter_val) { - best_filter_val = est; - best_filter = filter_type; - } - } - if (filter_type != best_filter) { // If the last iteration already got us the best filter, don't redo it - stbiw__encode_png_line((unsigned char*)(pixels), stride_bytes, x, y, j, n, best_filter, line_buffer); - filter_type = best_filter; - } +STBIWDEF unsigned char* stbi_write_png_to_mem( + const unsigned char* pixels, int stride_bytes, int x, int y, int n, int* out_len) { + int force_filter = stbi_write_force_png_filter; + int ctype[5] = {-1, 0, 4, 2, 6}; + unsigned char sig[8] = {137, 80, 78, 71, 13, 10, 26, 10}; + unsigned char *out, *o, *filt, *zlib; + signed char* line_buffer; + int j, zlen; + + if (stride_bytes == 0) + stride_bytes = x * n; + + if (force_filter >= 5) { + force_filter = -1; + } + + filt = (unsigned char*)STBIW_MALLOC((x * n + 1) * y); + if (!filt) + return 0; + line_buffer = (signed char*)STBIW_MALLOC(x * n); + if (!line_buffer) { + STBIW_FREE(filt); + return 0; + } + for (j = 0; j < y; ++j) { + int filter_type; + if (force_filter > -1) { + filter_type = force_filter; + stbiw__encode_png_line( + (unsigned char*)(pixels), stride_bytes, x, y, j, n, force_filter, line_buffer); + } else { // Estimate the best filter by running through all of them: + int best_filter = 0, best_filter_val = 0x7fffffff, est, i; + for (filter_type = 0; filter_type < 5; filter_type++) { + stbiw__encode_png_line( + (unsigned char*)(pixels), stride_bytes, x, y, j, n, filter_type, line_buffer); + + // Estimate the entropy of the line using this filter; the less, the + // better. + est = 0; + for (i = 0; i < x * n; ++i) { + est += abs((signed char)line_buffer[i]); + } + if (est < best_filter_val) { + best_filter_val = est; + best_filter = filter_type; + } + } + if (filter_type != best_filter) { // If the last iteration already got us + // the best filter, don't redo it + stbiw__encode_png_line( + (unsigned char*)(pixels), stride_bytes, x, y, j, n, best_filter, line_buffer); + filter_type = best_filter; } - // when we get here, filter_type contains the filter type, and line_buffer contains the data - filt[j*(x*n+1)] = (unsigned char) filter_type; - STBIW_MEMMOVE(filt+j*(x*n+1)+1, line_buffer, x*n); - } - STBIW_FREE(line_buffer); - zlib = stbi_zlib_compress(filt, y*( x*n+1), &zlen, stbi_write_png_compression_level); - STBIW_FREE(filt); - if (!zlib) return 0; - - // each tag requires 12 bytes of overhead - out = (unsigned char *) STBIW_MALLOC(8 + 12+13 + 12+zlen + 12); - if (!out) return 0; - *out_len = 8 + 12+13 + 12+zlen + 12; - - o=out; - STBIW_MEMMOVE(o,sig,8); o+= 8; - stbiw__wp32(o, 13); // header length - stbiw__wptag(o, "IHDR"); - stbiw__wp32(o, x); - stbiw__wp32(o, y); - *o++ = 8; - *o++ = STBIW_UCHAR(ctype[n]); - *o++ = 0; - *o++ = 0; - *o++ = 0; - stbiw__wpcrc(&o,13); - - stbiw__wp32(o, zlen); - stbiw__wptag(o, "IDAT"); - STBIW_MEMMOVE(o, zlib, zlen); - o += zlen; - STBIW_FREE(zlib); - stbiw__wpcrc(&o, zlen); - - stbiw__wp32(o,0); - stbiw__wptag(o, "IEND"); - stbiw__wpcrc(&o,0); - - STBIW_ASSERT(o == out + *out_len); - - return out; + } + // when we get here, filter_type contains the filter type, and line_buffer + // contains the data + filt[j * (x * n + 1)] = (unsigned char)filter_type; + STBIW_MEMMOVE(filt + j * (x * n + 1) + 1, line_buffer, x * n); + } + STBIW_FREE(line_buffer); + zlib = stbi_zlib_compress(filt, y * (x * n + 1), &zlen, stbi_write_png_compression_level); + STBIW_FREE(filt); + if (!zlib) + return 0; + + // each tag requires 12 bytes of overhead + out = (unsigned char*)STBIW_MALLOC(8 + 12 + 13 + 12 + zlen + 12); + if (!out) + return 0; + *out_len = 8 + 12 + 13 + 12 + zlen + 12; + + o = out; + STBIW_MEMMOVE(o, sig, 8); + o += 8; + stbiw__wp32(o, 13); // header length + stbiw__wptag(o, "IHDR"); + stbiw__wp32(o, x); + stbiw__wp32(o, y); + *o++ = 8; + *o++ = STBIW_UCHAR(ctype[n]); + *o++ = 0; + *o++ = 0; + *o++ = 0; + stbiw__wpcrc(&o, 13); + + stbiw__wp32(o, zlen); + stbiw__wptag(o, "IDAT"); + STBIW_MEMMOVE(o, zlib, zlen); + o += zlen; + STBIW_FREE(zlib); + stbiw__wpcrc(&o, zlen); + + stbiw__wp32(o, 0); + stbiw__wptag(o, "IEND"); + stbiw__wpcrc(&o, 0); + + STBIW_ASSERT(o == out + *out_len); + + return out; } #ifndef STBI_WRITE_NO_STDIO -STBIWDEF int stbi_write_png(char const *filename, int x, int y, int comp, const void *data, int stride_bytes) -{ - FILE *f; - int len; - unsigned char *png = stbi_write_png_to_mem((const unsigned char *) data, stride_bytes, x, y, comp, &len); - if (png == NULL) return 0; - - f = stbiw__fopen(filename, "wb"); - if (!f) { STBIW_FREE(png); return 0; } - fwrite(png, 1, len, f); - fclose(f); - STBIW_FREE(png); - return 1; +STBIWDEF int stbi_write_png( + const char* filename, int x, int y, int comp, const void* data, int stride_bytes) { + FILE* f; + int len; + unsigned char* png = + stbi_write_png_to_mem((const unsigned char*)data, stride_bytes, x, y, comp, &len); + if (png == NULL) + return 0; + + f = stbiw__fopen(filename, "wb"); + if (!f) { + STBIW_FREE(png); + return 0; + } + fwrite(png, 1, len, f); + fclose(f); + STBIW_FREE(png); + return 1; } #endif -STBIWDEF int stbi_write_png_to_func(stbi_write_func *func, void *context, int x, int y, int comp, const void *data, int stride_bytes) -{ - int len; - unsigned char *png = stbi_write_png_to_mem((const unsigned char *) data, stride_bytes, x, y, comp, &len); - if (png == NULL) return 0; - func(context, png, len); - STBIW_FREE(png); - return 1; +STBIWDEF int stbi_write_png_to_func(stbi_write_func* func, + void* context, + int x, + int y, + int comp, + const void* data, + int stride_bytes) { + int len; + unsigned char* png = + stbi_write_png_to_mem((const unsigned char*)data, stride_bytes, x, y, comp, &len); + if (png == NULL) + return 0; + func(context, png, len); + STBIW_FREE(png); + return 1; } - /* *************************************************************************** * * JPEG writer * * This is based on Jon Olick's jo_jpeg.cpp: - * public domain Simple, Minimalistic JPEG writer - http://www.jonolick.com/code.html + * public domain Simple, Minimalistic JPEG writer - + * http://www.jonolick.com/code.html */ -static const unsigned char stbiw__jpg_ZigZag[] = { 0,1,5,6,14,15,27,28,2,4,7,13,16,26,29,42,3,8,12,17,25,30,41,43,9,11,18, - 24,31,40,44,53,10,19,23,32,39,45,52,54,20,22,33,38,46,51,55,60,21,34,37,47,50,56,59,61,35,36,48,49,57,58,62,63 }; - -static void stbiw__jpg_writeBits(stbi__write_context *s, int *bitBufP, int *bitCntP, const unsigned short *bs) { - int bitBuf = *bitBufP, bitCnt = *bitCntP; - bitCnt += bs[1]; - bitBuf |= bs[0] << (24 - bitCnt); - while(bitCnt >= 8) { - unsigned char c = (bitBuf >> 16) & 255; - stbiw__putc(s, c); - if(c == 255) { - stbiw__putc(s, 0); - } - bitBuf <<= 8; - bitCnt -= 8; - } - *bitBufP = bitBuf; - *bitCntP = bitCnt; +static const unsigned char stbiw__jpg_ZigZag[] = { + 0, 1, 5, 6, 14, 15, 27, 28, 2, 4, 7, 13, 16, 26, 29, 42, 3, 8, 12, 17, 25, 30, + 41, 43, 9, 11, 18, 24, 31, 40, 44, 53, 10, 19, 23, 32, 39, 45, 52, 54, 20, 22, 33, 38, + 46, 51, 55, 60, 21, 34, 37, 47, 50, 56, 59, 61, 35, 36, 48, 49, 57, 58, 62, 63}; + +static void stbiw__jpg_writeBits(stbi__write_context* s, + int* bitBufP, + int* bitCntP, + const unsigned short* bs) { + int bitBuf = *bitBufP, bitCnt = *bitCntP; + bitCnt += bs[1]; + bitBuf |= bs[0] << (24 - bitCnt); + while (bitCnt >= 8) { + unsigned char c = (bitBuf >> 16) & 255; + stbiw__putc(s, c); + if (c == 255) { + stbiw__putc(s, 0); + } + bitBuf <<= 8; + bitCnt -= 8; + } + *bitBufP = bitBuf; + *bitCntP = bitCnt; } -static void stbiw__jpg_DCT(float *d0p, float *d1p, float *d2p, float *d3p, float *d4p, float *d5p, float *d6p, float *d7p) { - float d0 = *d0p, d1 = *d1p, d2 = *d2p, d3 = *d3p, d4 = *d4p, d5 = *d5p, d6 = *d6p, d7 = *d7p; - float z1, z2, z3, z4, z5, z11, z13; - - float tmp0 = d0 + d7; - float tmp7 = d0 - d7; - float tmp1 = d1 + d6; - float tmp6 = d1 - d6; - float tmp2 = d2 + d5; - float tmp5 = d2 - d5; - float tmp3 = d3 + d4; - float tmp4 = d3 - d4; - - // Even part - float tmp10 = tmp0 + tmp3; // phase 2 - float tmp13 = tmp0 - tmp3; - float tmp11 = tmp1 + tmp2; - float tmp12 = tmp1 - tmp2; - - d0 = tmp10 + tmp11; // phase 3 - d4 = tmp10 - tmp11; - - z1 = (tmp12 + tmp13) * 0.707106781f; // c4 - d2 = tmp13 + z1; // phase 5 - d6 = tmp13 - z1; - - // Odd part - tmp10 = tmp4 + tmp5; // phase 2 - tmp11 = tmp5 + tmp6; - tmp12 = tmp6 + tmp7; - - // The rotator is modified from fig 4-8 to avoid extra negations. - z5 = (tmp10 - tmp12) * 0.382683433f; // c6 - z2 = tmp10 * 0.541196100f + z5; // c2-c6 - z4 = tmp12 * 1.306562965f + z5; // c2+c6 - z3 = tmp11 * 0.707106781f; // c4 - - z11 = tmp7 + z3; // phase 5 - z13 = tmp7 - z3; - - *d5p = z13 + z2; // phase 6 - *d3p = z13 - z2; - *d1p = z11 + z4; - *d7p = z11 - z4; - - *d0p = d0; *d2p = d2; *d4p = d4; *d6p = d6; +static void stbiw__jpg_DCT(float* d0p, + float* d1p, + float* d2p, + float* d3p, + float* d4p, + float* d5p, + float* d6p, + float* d7p) { + float d0 = *d0p, d1 = *d1p, d2 = *d2p, d3 = *d3p, d4 = *d4p, d5 = *d5p, d6 = *d6p, d7 = *d7p; + float z1, z2, z3, z4, z5, z11, z13; + + float tmp0 = d0 + d7; + float tmp7 = d0 - d7; + float tmp1 = d1 + d6; + float tmp6 = d1 - d6; + float tmp2 = d2 + d5; + float tmp5 = d2 - d5; + float tmp3 = d3 + d4; + float tmp4 = d3 - d4; + + // Even part + float tmp10 = tmp0 + tmp3; // phase 2 + float tmp13 = tmp0 - tmp3; + float tmp11 = tmp1 + tmp2; + float tmp12 = tmp1 - tmp2; + + d0 = tmp10 + tmp11; // phase 3 + d4 = tmp10 - tmp11; + + z1 = (tmp12 + tmp13) * 0.707106781f; // c4 + d2 = tmp13 + z1; // phase 5 + d6 = tmp13 - z1; + + // Odd part + tmp10 = tmp4 + tmp5; // phase 2 + tmp11 = tmp5 + tmp6; + tmp12 = tmp6 + tmp7; + + // The rotator is modified from fig 4-8 to avoid extra negations. + z5 = (tmp10 - tmp12) * 0.382683433f; // c6 + z2 = tmp10 * 0.541196100f + z5; // c2-c6 + z4 = tmp12 * 1.306562965f + z5; // c2+c6 + z3 = tmp11 * 0.707106781f; // c4 + + z11 = tmp7 + z3; // phase 5 + z13 = tmp7 - z3; + + *d5p = z13 + z2; // phase 6 + *d3p = z13 - z2; + *d1p = z11 + z4; + *d7p = z11 - z4; + + *d0p = d0; + *d2p = d2; + *d4p = d4; + *d6p = d6; } static void stbiw__jpg_calcBits(int val, unsigned short bits[2]) { - int tmp1 = val < 0 ? -val : val; - val = val < 0 ? val-1 : val; - bits[1] = 1; - while(tmp1 >>= 1) { - ++bits[1]; - } - bits[0] = val & ((1<>= 1) { + ++bits[1]; + } + bits[0] = val & ((1 << bits[1]) - 1); } -static int stbiw__jpg_processDU(stbi__write_context *s, int *bitBuf, int *bitCnt, float *CDU, int du_stride, float *fdtbl, int DC, const unsigned short HTDC[256][2], const unsigned short HTAC[256][2]) { - const unsigned short EOB[2] = { HTAC[0x00][0], HTAC[0x00][1] }; - const unsigned short M16zeroes[2] = { HTAC[0xF0][0], HTAC[0xF0][1] }; - int dataOff, i, j, n, diff, end0pos, x, y; - int DU[64]; - - // DCT rows - for(dataOff=0, n=du_stride*8; dataOff0)&&(DU[end0pos]==0); --end0pos) { - } - // end0pos = first element in reverse order !=0 - if(end0pos == 0) { - stbiw__jpg_writeBits(s, bitBuf, bitCnt, EOB); - return DU[0]; - } - for(i = 1; i <= end0pos; ++i) { - int startpos = i; - int nrzeroes; - unsigned short bits[2]; - for (; DU[i]==0 && i<=end0pos; ++i) { - } - nrzeroes = i-startpos; - if ( nrzeroes >= 16 ) { - int lng = nrzeroes>>4; - int nrmarker; - for (nrmarker=1; nrmarker <= lng; ++nrmarker) - stbiw__jpg_writeBits(s, bitBuf, bitCnt, M16zeroes); - nrzeroes &= 15; - } - stbiw__jpg_calcBits(DU[i], bits); - stbiw__jpg_writeBits(s, bitBuf, bitCnt, HTAC[(nrzeroes<<4)+bits[1]]); - stbiw__jpg_writeBits(s, bitBuf, bitCnt, bits); - } - if(end0pos != 63) { - stbiw__jpg_writeBits(s, bitBuf, bitCnt, EOB); - } - return DU[0]; +static int stbiw__jpg_processDU(stbi__write_context* s, + int* bitBuf, + int* bitCnt, + float* CDU, + int du_stride, + float* fdtbl, + int DC, + const unsigned short HTDC[256][2], + const unsigned short HTAC[256][2]) { + const unsigned short EOB[2] = {HTAC[0x00][0], HTAC[0x00][1]}; + const unsigned short M16zeroes[2] = {HTAC[0xF0][0], HTAC[0xF0][1]}; + int dataOff, i, j, n, diff, end0pos, x, y; + int DU[64]; + + // DCT rows + for (dataOff = 0, n = du_stride * 8; dataOff < n; dataOff += du_stride) { + stbiw__jpg_DCT(&CDU[dataOff], + &CDU[dataOff + 1], + &CDU[dataOff + 2], + &CDU[dataOff + 3], + &CDU[dataOff + 4], + &CDU[dataOff + 5], + &CDU[dataOff + 6], + &CDU[dataOff + 7]); + } + // DCT columns + for (dataOff = 0; dataOff < 8; ++dataOff) { + stbiw__jpg_DCT(&CDU[dataOff], + &CDU[dataOff + du_stride], + &CDU[dataOff + du_stride * 2], + &CDU[dataOff + du_stride * 3], + &CDU[dataOff + du_stride * 4], + &CDU[dataOff + du_stride * 5], + &CDU[dataOff + du_stride * 6], + &CDU[dataOff + du_stride * 7]); + } + // Quantize/descale/zigzag the coefficients + for (y = 0, j = 0; y < 8; ++y) { + for (x = 0; x < 8; ++x, ++j) { + float v; + i = y * du_stride + x; + v = CDU[i] * fdtbl[j]; + // DU[stbiw__jpg_ZigZag[j]] = (int)(v < 0 ? ceilf(v - 0.5f) : floorf(v + + // 0.5f)); ceilf() and floorf() are C99, not C89, but I /think/ they're + // not needed here anyway? + DU[stbiw__jpg_ZigZag[j]] = (int)(v < 0 ? v - 0.5f : v + 0.5f); + } + } + + // Encode DC + diff = DU[0] - DC; + if (diff == 0) { + stbiw__jpg_writeBits(s, bitBuf, bitCnt, HTDC[0]); + } else { + unsigned short bits[2]; + stbiw__jpg_calcBits(diff, bits); + stbiw__jpg_writeBits(s, bitBuf, bitCnt, HTDC[bits[1]]); + stbiw__jpg_writeBits(s, bitBuf, bitCnt, bits); + } + // Encode ACs + end0pos = 63; + for (; (end0pos > 0) && (DU[end0pos] == 0); --end0pos) { + } + // end0pos = first element in reverse order !=0 + if (end0pos == 0) { + stbiw__jpg_writeBits(s, bitBuf, bitCnt, EOB); + return DU[0]; + } + for (i = 1; i <= end0pos; ++i) { + int startpos = i; + int nrzeroes; + unsigned short bits[2]; + for (; DU[i] == 0 && i <= end0pos; ++i) { + } + nrzeroes = i - startpos; + if (nrzeroes >= 16) { + int lng = nrzeroes >> 4; + int nrmarker; + for (nrmarker = 1; nrmarker <= lng; ++nrmarker) + stbiw__jpg_writeBits(s, bitBuf, bitCnt, M16zeroes); + nrzeroes &= 15; + } + stbiw__jpg_calcBits(DU[i], bits); + stbiw__jpg_writeBits(s, bitBuf, bitCnt, HTAC[(nrzeroes << 4) + bits[1]]); + stbiw__jpg_writeBits(s, bitBuf, bitCnt, bits); + } + if (end0pos != 63) { + stbiw__jpg_writeBits(s, bitBuf, bitCnt, EOB); + } + return DU[0]; } -static int stbi_write_jpg_core(stbi__write_context *s, int width, int height, int comp, const void* data, int quality) { - // Constants that don't pollute global namespace - static const unsigned char std_dc_luminance_nrcodes[] = {0,0,1,5,1,1,1,1,1,1,0,0,0,0,0,0,0}; - static const unsigned char std_dc_luminance_values[] = {0,1,2,3,4,5,6,7,8,9,10,11}; - static const unsigned char std_ac_luminance_nrcodes[] = {0,0,2,1,3,3,2,4,3,5,5,4,4,0,0,1,0x7d}; - static const unsigned char std_ac_luminance_values[] = { - 0x01,0x02,0x03,0x00,0x04,0x11,0x05,0x12,0x21,0x31,0x41,0x06,0x13,0x51,0x61,0x07,0x22,0x71,0x14,0x32,0x81,0x91,0xa1,0x08, - 0x23,0x42,0xb1,0xc1,0x15,0x52,0xd1,0xf0,0x24,0x33,0x62,0x72,0x82,0x09,0x0a,0x16,0x17,0x18,0x19,0x1a,0x25,0x26,0x27,0x28, - 0x29,0x2a,0x34,0x35,0x36,0x37,0x38,0x39,0x3a,0x43,0x44,0x45,0x46,0x47,0x48,0x49,0x4a,0x53,0x54,0x55,0x56,0x57,0x58,0x59, - 0x5a,0x63,0x64,0x65,0x66,0x67,0x68,0x69,0x6a,0x73,0x74,0x75,0x76,0x77,0x78,0x79,0x7a,0x83,0x84,0x85,0x86,0x87,0x88,0x89, - 0x8a,0x92,0x93,0x94,0x95,0x96,0x97,0x98,0x99,0x9a,0xa2,0xa3,0xa4,0xa5,0xa6,0xa7,0xa8,0xa9,0xaa,0xb2,0xb3,0xb4,0xb5,0xb6, - 0xb7,0xb8,0xb9,0xba,0xc2,0xc3,0xc4,0xc5,0xc6,0xc7,0xc8,0xc9,0xca,0xd2,0xd3,0xd4,0xd5,0xd6,0xd7,0xd8,0xd9,0xda,0xe1,0xe2, - 0xe3,0xe4,0xe5,0xe6,0xe7,0xe8,0xe9,0xea,0xf1,0xf2,0xf3,0xf4,0xf5,0xf6,0xf7,0xf8,0xf9,0xfa - }; - static const unsigned char std_dc_chrominance_nrcodes[] = {0,0,3,1,1,1,1,1,1,1,1,1,0,0,0,0,0}; - static const unsigned char std_dc_chrominance_values[] = {0,1,2,3,4,5,6,7,8,9,10,11}; - static const unsigned char std_ac_chrominance_nrcodes[] = {0,0,2,1,2,4,4,3,4,7,5,4,4,0,1,2,0x77}; - static const unsigned char std_ac_chrominance_values[] = { - 0x00,0x01,0x02,0x03,0x11,0x04,0x05,0x21,0x31,0x06,0x12,0x41,0x51,0x07,0x61,0x71,0x13,0x22,0x32,0x81,0x08,0x14,0x42,0x91, - 0xa1,0xb1,0xc1,0x09,0x23,0x33,0x52,0xf0,0x15,0x62,0x72,0xd1,0x0a,0x16,0x24,0x34,0xe1,0x25,0xf1,0x17,0x18,0x19,0x1a,0x26, - 0x27,0x28,0x29,0x2a,0x35,0x36,0x37,0x38,0x39,0x3a,0x43,0x44,0x45,0x46,0x47,0x48,0x49,0x4a,0x53,0x54,0x55,0x56,0x57,0x58, - 0x59,0x5a,0x63,0x64,0x65,0x66,0x67,0x68,0x69,0x6a,0x73,0x74,0x75,0x76,0x77,0x78,0x79,0x7a,0x82,0x83,0x84,0x85,0x86,0x87, - 0x88,0x89,0x8a,0x92,0x93,0x94,0x95,0x96,0x97,0x98,0x99,0x9a,0xa2,0xa3,0xa4,0xa5,0xa6,0xa7,0xa8,0xa9,0xaa,0xb2,0xb3,0xb4, - 0xb5,0xb6,0xb7,0xb8,0xb9,0xba,0xc2,0xc3,0xc4,0xc5,0xc6,0xc7,0xc8,0xc9,0xca,0xd2,0xd3,0xd4,0xd5,0xd6,0xd7,0xd8,0xd9,0xda, - 0xe2,0xe3,0xe4,0xe5,0xe6,0xe7,0xe8,0xe9,0xea,0xf2,0xf3,0xf4,0xf5,0xf6,0xf7,0xf8,0xf9,0xfa - }; - // Huffman tables - static const unsigned short YDC_HT[256][2] = { {0,2},{2,3},{3,3},{4,3},{5,3},{6,3},{14,4},{30,5},{62,6},{126,7},{254,8},{510,9}}; - static const unsigned short UVDC_HT[256][2] = { {0,2},{1,2},{2,2},{6,3},{14,4},{30,5},{62,6},{126,7},{254,8},{510,9},{1022,10},{2046,11}}; - static const unsigned short YAC_HT[256][2] = { - {10,4},{0,2},{1,2},{4,3},{11,4},{26,5},{120,7},{248,8},{1014,10},{65410,16},{65411,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {12,4},{27,5},{121,7},{502,9},{2038,11},{65412,16},{65413,16},{65414,16},{65415,16},{65416,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {28,5},{249,8},{1015,10},{4084,12},{65417,16},{65418,16},{65419,16},{65420,16},{65421,16},{65422,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {58,6},{503,9},{4085,12},{65423,16},{65424,16},{65425,16},{65426,16},{65427,16},{65428,16},{65429,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {59,6},{1016,10},{65430,16},{65431,16},{65432,16},{65433,16},{65434,16},{65435,16},{65436,16},{65437,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {122,7},{2039,11},{65438,16},{65439,16},{65440,16},{65441,16},{65442,16},{65443,16},{65444,16},{65445,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {123,7},{4086,12},{65446,16},{65447,16},{65448,16},{65449,16},{65450,16},{65451,16},{65452,16},{65453,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {250,8},{4087,12},{65454,16},{65455,16},{65456,16},{65457,16},{65458,16},{65459,16},{65460,16},{65461,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {504,9},{32704,15},{65462,16},{65463,16},{65464,16},{65465,16},{65466,16},{65467,16},{65468,16},{65469,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {505,9},{65470,16},{65471,16},{65472,16},{65473,16},{65474,16},{65475,16},{65476,16},{65477,16},{65478,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {506,9},{65479,16},{65480,16},{65481,16},{65482,16},{65483,16},{65484,16},{65485,16},{65486,16},{65487,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {1017,10},{65488,16},{65489,16},{65490,16},{65491,16},{65492,16},{65493,16},{65494,16},{65495,16},{65496,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {1018,10},{65497,16},{65498,16},{65499,16},{65500,16},{65501,16},{65502,16},{65503,16},{65504,16},{65505,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {2040,11},{65506,16},{65507,16},{65508,16},{65509,16},{65510,16},{65511,16},{65512,16},{65513,16},{65514,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {65515,16},{65516,16},{65517,16},{65518,16},{65519,16},{65520,16},{65521,16},{65522,16},{65523,16},{65524,16},{0,0},{0,0},{0,0},{0,0},{0,0}, - {2041,11},{65525,16},{65526,16},{65527,16},{65528,16},{65529,16},{65530,16},{65531,16},{65532,16},{65533,16},{65534,16},{0,0},{0,0},{0,0},{0,0},{0,0} - }; - static const unsigned short UVAC_HT[256][2] = { - {0,2},{1,2},{4,3},{10,4},{24,5},{25,5},{56,6},{120,7},{500,9},{1014,10},{4084,12},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {11,4},{57,6},{246,8},{501,9},{2038,11},{4085,12},{65416,16},{65417,16},{65418,16},{65419,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {26,5},{247,8},{1015,10},{4086,12},{32706,15},{65420,16},{65421,16},{65422,16},{65423,16},{65424,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {27,5},{248,8},{1016,10},{4087,12},{65425,16},{65426,16},{65427,16},{65428,16},{65429,16},{65430,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {58,6},{502,9},{65431,16},{65432,16},{65433,16},{65434,16},{65435,16},{65436,16},{65437,16},{65438,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {59,6},{1017,10},{65439,16},{65440,16},{65441,16},{65442,16},{65443,16},{65444,16},{65445,16},{65446,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {121,7},{2039,11},{65447,16},{65448,16},{65449,16},{65450,16},{65451,16},{65452,16},{65453,16},{65454,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {122,7},{2040,11},{65455,16},{65456,16},{65457,16},{65458,16},{65459,16},{65460,16},{65461,16},{65462,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {249,8},{65463,16},{65464,16},{65465,16},{65466,16},{65467,16},{65468,16},{65469,16},{65470,16},{65471,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {503,9},{65472,16},{65473,16},{65474,16},{65475,16},{65476,16},{65477,16},{65478,16},{65479,16},{65480,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {504,9},{65481,16},{65482,16},{65483,16},{65484,16},{65485,16},{65486,16},{65487,16},{65488,16},{65489,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {505,9},{65490,16},{65491,16},{65492,16},{65493,16},{65494,16},{65495,16},{65496,16},{65497,16},{65498,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {506,9},{65499,16},{65500,16},{65501,16},{65502,16},{65503,16},{65504,16},{65505,16},{65506,16},{65507,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {2041,11},{65508,16},{65509,16},{65510,16},{65511,16},{65512,16},{65513,16},{65514,16},{65515,16},{65516,16},{0,0},{0,0},{0,0},{0,0},{0,0},{0,0}, - {16352,14},{65517,16},{65518,16},{65519,16},{65520,16},{65521,16},{65522,16},{65523,16},{65524,16},{65525,16},{0,0},{0,0},{0,0},{0,0},{0,0}, - {1018,10},{32707,15},{65526,16},{65527,16},{65528,16},{65529,16},{65530,16},{65531,16},{65532,16},{65533,16},{65534,16},{0,0},{0,0},{0,0},{0,0},{0,0} - }; - static const int YQT[] = {16,11,10,16,24,40,51,61,12,12,14,19,26,58,60,55,14,13,16,24,40,57,69,56,14,17,22,29,51,87,80,62,18,22, - 37,56,68,109,103,77,24,35,55,64,81,104,113,92,49,64,78,87,103,121,120,101,72,92,95,98,112,100,103,99}; - static const int UVQT[] = {17,18,24,47,99,99,99,99,18,21,26,66,99,99,99,99,24,26,56,99,99,99,99,99,47,66,99,99,99,99,99,99, - 99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99}; - static const float aasf[] = { 1.0f * 2.828427125f, 1.387039845f * 2.828427125f, 1.306562965f * 2.828427125f, 1.175875602f * 2.828427125f, - 1.0f * 2.828427125f, 0.785694958f * 2.828427125f, 0.541196100f * 2.828427125f, 0.275899379f * 2.828427125f }; - - int row, col, i, k, subsample; - float fdtbl_Y[64], fdtbl_UV[64]; - unsigned char YTable[64], UVTable[64]; - - if(!data || !width || !height || comp > 4 || comp < 1) { - return 0; - } - - quality = quality ? quality : 90; - subsample = quality <= 90 ? 1 : 0; - quality = quality < 1 ? 1 : quality > 100 ? 100 : quality; - quality = quality < 50 ? 5000 / quality : 200 - quality * 2; - - for(i = 0; i < 64; ++i) { - int uvti, yti = (YQT[i]*quality+50)/100; - YTable[stbiw__jpg_ZigZag[i]] = (unsigned char) (yti < 1 ? 1 : yti > 255 ? 255 : yti); - uvti = (UVQT[i]*quality+50)/100; - UVTable[stbiw__jpg_ZigZag[i]] = (unsigned char) (uvti < 1 ? 1 : uvti > 255 ? 255 : uvti); - } - - for(row = 0, k = 0; row < 8; ++row) { - for(col = 0; col < 8; ++col, ++k) { - fdtbl_Y[k] = 1 / (YTable [stbiw__jpg_ZigZag[k]] * aasf[row] * aasf[col]); - fdtbl_UV[k] = 1 / (UVTable[stbiw__jpg_ZigZag[k]] * aasf[row] * aasf[col]); - } - } - - // Write Headers - { - static const unsigned char head0[] = { 0xFF,0xD8,0xFF,0xE0,0,0x10,'J','F','I','F',0,1,1,0,0,1,0,1,0,0,0xFF,0xDB,0,0x84,0 }; - static const unsigned char head2[] = { 0xFF,0xDA,0,0xC,3,1,0,2,0x11,3,0x11,0,0x3F,0 }; - const unsigned char head1[] = { 0xFF,0xC0,0,0x11,8,(unsigned char)(height>>8),STBIW_UCHAR(height),(unsigned char)(width>>8),STBIW_UCHAR(width), - 3,1,(unsigned char)(subsample?0x22:0x11),0,2,0x11,1,3,0x11,1,0xFF,0xC4,0x01,0xA2,0 }; - s->func(s->context, (void*)head0, sizeof(head0)); - s->func(s->context, (void*)YTable, sizeof(YTable)); - stbiw__putc(s, 1); - s->func(s->context, UVTable, sizeof(UVTable)); - s->func(s->context, (void*)head1, sizeof(head1)); - s->func(s->context, (void*)(std_dc_luminance_nrcodes+1), sizeof(std_dc_luminance_nrcodes)-1); - s->func(s->context, (void*)std_dc_luminance_values, sizeof(std_dc_luminance_values)); - stbiw__putc(s, 0x10); // HTYACinfo - s->func(s->context, (void*)(std_ac_luminance_nrcodes+1), sizeof(std_ac_luminance_nrcodes)-1); - s->func(s->context, (void*)std_ac_luminance_values, sizeof(std_ac_luminance_values)); - stbiw__putc(s, 1); // HTUDCinfo - s->func(s->context, (void*)(std_dc_chrominance_nrcodes+1), sizeof(std_dc_chrominance_nrcodes)-1); - s->func(s->context, (void*)std_dc_chrominance_values, sizeof(std_dc_chrominance_values)); - stbiw__putc(s, 0x11); // HTUACinfo - s->func(s->context, (void*)(std_ac_chrominance_nrcodes+1), sizeof(std_ac_chrominance_nrcodes)-1); - s->func(s->context, (void*)std_ac_chrominance_values, sizeof(std_ac_chrominance_values)); - s->func(s->context, (void*)head2, sizeof(head2)); - } - - // Encode 8x8 macroblocks - { - static const unsigned short fillBits[] = {0x7F, 7}; - int DCY=0, DCU=0, DCV=0; - int bitBuf=0, bitCnt=0; - // comp == 2 is grey+alpha (alpha is ignored) - int ofsG = comp > 2 ? 1 : 0, ofsB = comp > 2 ? 2 : 0; - const unsigned char *dataR = (const unsigned char *)data; - const unsigned char *dataG = dataR + ofsG; - const unsigned char *dataB = dataR + ofsB; - int x, y, pos; - if(subsample) { - for(y = 0; y < height; y += 16) { - for(x = 0; x < width; x += 16) { - float Y[256], U[256], V[256]; - for(row = y, pos = 0; row < y+16; ++row) { - // row >= height => use last input row - int clamped_row = (row < height) ? row : height - 1; - int base_p = (stbi__flip_vertically_on_write ? (height-1-clamped_row) : clamped_row)*width*comp; - for(col = x; col < x+16; ++col, ++pos) { - // if col >= width => use pixel from last input column - int p = base_p + ((col < width) ? col : (width-1))*comp; - float r = dataR[p], g = dataG[p], b = dataB[p]; - Y[pos]= +0.29900f*r + 0.58700f*g + 0.11400f*b - 128; - U[pos]= -0.16874f*r - 0.33126f*g + 0.50000f*b; - V[pos]= +0.50000f*r - 0.41869f*g - 0.08131f*b; - } - } - DCY = stbiw__jpg_processDU(s, &bitBuf, &bitCnt, Y+0, 16, fdtbl_Y, DCY, YDC_HT, YAC_HT); - DCY = stbiw__jpg_processDU(s, &bitBuf, &bitCnt, Y+8, 16, fdtbl_Y, DCY, YDC_HT, YAC_HT); - DCY = stbiw__jpg_processDU(s, &bitBuf, &bitCnt, Y+128, 16, fdtbl_Y, DCY, YDC_HT, YAC_HT); - DCY = stbiw__jpg_processDU(s, &bitBuf, &bitCnt, Y+136, 16, fdtbl_Y, DCY, YDC_HT, YAC_HT); - - // subsample U,V - { - float subU[64], subV[64]; - int yy, xx; - for(yy = 0, pos = 0; yy < 8; ++yy) { - for(xx = 0; xx < 8; ++xx, ++pos) { - int j = yy*32+xx*2; - subU[pos] = (U[j+0] + U[j+1] + U[j+16] + U[j+17]) * 0.25f; - subV[pos] = (V[j+0] + V[j+1] + V[j+16] + V[j+17]) * 0.25f; - } - } - DCU = stbiw__jpg_processDU(s, &bitBuf, &bitCnt, subU, 8, fdtbl_UV, DCU, UVDC_HT, UVAC_HT); - DCV = stbiw__jpg_processDU(s, &bitBuf, &bitCnt, subV, 8, fdtbl_UV, DCV, UVDC_HT, UVAC_HT); - } +static int stbi_write_jpg_core( + stbi__write_context* s, int width, int height, int comp, const void* data, int quality) { + // Constants that don't pollute global namespace + static const unsigned char std_dc_luminance_nrcodes[] = { + 0, 0, 1, 5, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0}; + static const unsigned char std_dc_luminance_values[] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11}; + static const unsigned char std_ac_luminance_nrcodes[] = { + 0, 0, 2, 1, 3, 3, 2, 4, 3, 5, 5, 4, 4, 0, 0, 1, 0x7d}; + static const unsigned char std_ac_luminance_values[] = { + 0x01, 0x02, 0x03, 0x00, 0x04, 0x11, 0x05, 0x12, 0x21, 0x31, 0x41, 0x06, 0x13, 0x51, 0x61, + 0x07, 0x22, 0x71, 0x14, 0x32, 0x81, 0x91, 0xa1, 0x08, 0x23, 0x42, 0xb1, 0xc1, 0x15, 0x52, + 0xd1, 0xf0, 0x24, 0x33, 0x62, 0x72, 0x82, 0x09, 0x0a, 0x16, 0x17, 0x18, 0x19, 0x1a, 0x25, + 0x26, 0x27, 0x28, 0x29, 0x2a, 0x34, 0x35, 0x36, 0x37, 0x38, 0x39, 0x3a, 0x43, 0x44, 0x45, + 0x46, 0x47, 0x48, 0x49, 0x4a, 0x53, 0x54, 0x55, 0x56, 0x57, 0x58, 0x59, 0x5a, 0x63, 0x64, + 0x65, 0x66, 0x67, 0x68, 0x69, 0x6a, 0x73, 0x74, 0x75, 0x76, 0x77, 0x78, 0x79, 0x7a, 0x83, + 0x84, 0x85, 0x86, 0x87, 0x88, 0x89, 0x8a, 0x92, 0x93, 0x94, 0x95, 0x96, 0x97, 0x98, 0x99, + 0x9a, 0xa2, 0xa3, 0xa4, 0xa5, 0xa6, 0xa7, 0xa8, 0xa9, 0xaa, 0xb2, 0xb3, 0xb4, 0xb5, 0xb6, + 0xb7, 0xb8, 0xb9, 0xba, 0xc2, 0xc3, 0xc4, 0xc5, 0xc6, 0xc7, 0xc8, 0xc9, 0xca, 0xd2, 0xd3, + 0xd4, 0xd5, 0xd6, 0xd7, 0xd8, 0xd9, 0xda, 0xe1, 0xe2, 0xe3, 0xe4, 0xe5, 0xe6, 0xe7, 0xe8, + 0xe9, 0xea, 0xf1, 0xf2, 0xf3, 0xf4, 0xf5, 0xf6, 0xf7, 0xf8, 0xf9, 0xfa}; + static const unsigned char std_dc_chrominance_nrcodes[] = { + 0, 0, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0}; + static const unsigned char std_dc_chrominance_values[] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11}; + static const unsigned char std_ac_chrominance_nrcodes[] = { + 0, 0, 2, 1, 2, 4, 4, 3, 4, 7, 5, 4, 4, 0, 1, 2, 0x77}; + static const unsigned char std_ac_chrominance_values[] = { + 0x00, 0x01, 0x02, 0x03, 0x11, 0x04, 0x05, 0x21, 0x31, 0x06, 0x12, 0x41, 0x51, 0x07, 0x61, + 0x71, 0x13, 0x22, 0x32, 0x81, 0x08, 0x14, 0x42, 0x91, 0xa1, 0xb1, 0xc1, 0x09, 0x23, 0x33, + 0x52, 0xf0, 0x15, 0x62, 0x72, 0xd1, 0x0a, 0x16, 0x24, 0x34, 0xe1, 0x25, 0xf1, 0x17, 0x18, + 0x19, 0x1a, 0x26, 0x27, 0x28, 0x29, 0x2a, 0x35, 0x36, 0x37, 0x38, 0x39, 0x3a, 0x43, 0x44, + 0x45, 0x46, 0x47, 0x48, 0x49, 0x4a, 0x53, 0x54, 0x55, 0x56, 0x57, 0x58, 0x59, 0x5a, 0x63, + 0x64, 0x65, 0x66, 0x67, 0x68, 0x69, 0x6a, 0x73, 0x74, 0x75, 0x76, 0x77, 0x78, 0x79, 0x7a, + 0x82, 0x83, 0x84, 0x85, 0x86, 0x87, 0x88, 0x89, 0x8a, 0x92, 0x93, 0x94, 0x95, 0x96, 0x97, + 0x98, 0x99, 0x9a, 0xa2, 0xa3, 0xa4, 0xa5, 0xa6, 0xa7, 0xa8, 0xa9, 0xaa, 0xb2, 0xb3, 0xb4, + 0xb5, 0xb6, 0xb7, 0xb8, 0xb9, 0xba, 0xc2, 0xc3, 0xc4, 0xc5, 0xc6, 0xc7, 0xc8, 0xc9, 0xca, + 0xd2, 0xd3, 0xd4, 0xd5, 0xd6, 0xd7, 0xd8, 0xd9, 0xda, 0xe2, 0xe3, 0xe4, 0xe5, 0xe6, 0xe7, + 0xe8, 0xe9, 0xea, 0xf2, 0xf3, 0xf4, 0xf5, 0xf6, 0xf7, 0xf8, 0xf9, 0xfa}; + // Huffman tables + static const unsigned short YDC_HT[256][2] = {{0, 2}, + {2, 3}, + {3, 3}, + {4, 3}, + {5, 3}, + {6, 3}, + {14, 4}, + {30, 5}, + {62, 6}, + {126, 7}, + {254, 8}, + {510, 9}}; + static const unsigned short UVDC_HT[256][2] = {{0, 2}, + {1, 2}, + {2, 2}, + {6, 3}, + {14, 4}, + {30, 5}, + {62, 6}, + {126, 7}, + {254, 8}, + {510, 9}, + {1022, 10}, + {2046, 11}}; + static const unsigned short YAC_HT[256][2] = { + {10, 4}, {0, 2}, {1, 2}, {4, 3}, {11, 4}, {26, 5}, {120, 7}, + {248, 8}, {1014, 10}, {65410, 16}, {65411, 16}, {0, 0}, {0, 0}, {0, 0}, + {0, 0}, {0, 0}, {0, 0}, {12, 4}, {27, 5}, {121, 7}, {502, 9}, + {2038, 11}, {65412, 16}, {65413, 16}, {65414, 16}, {65415, 16}, {65416, 16}, {0, 0}, + {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {28, 5}, {249, 8}, + {1015, 10}, {4084, 12}, {65417, 16}, {65418, 16}, {65419, 16}, {65420, 16}, {65421, 16}, + {65422, 16}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, + {58, 6}, {503, 9}, {4085, 12}, {65423, 16}, {65424, 16}, {65425, 16}, {65426, 16}, + {65427, 16}, {65428, 16}, {65429, 16}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, + {0, 0}, {0, 0}, {59, 6}, {1016, 10}, {65430, 16}, {65431, 16}, {65432, 16}, + {65433, 16}, {65434, 16}, {65435, 16}, {65436, 16}, {65437, 16}, {0, 0}, {0, 0}, + {0, 0}, {0, 0}, {0, 0}, {0, 0}, {122, 7}, {2039, 11}, {65438, 16}, + {65439, 16}, {65440, 16}, {65441, 16}, {65442, 16}, {65443, 16}, {65444, 16}, {65445, 16}, + {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {123, 7}, + {4086, 12}, {65446, 16}, {65447, 16}, {65448, 16}, {65449, 16}, {65450, 16}, {65451, 16}, + {65452, 16}, {65453, 16}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, + {0, 0}, {250, 8}, {4087, 12}, {65454, 16}, {65455, 16}, {65456, 16}, {65457, 16}, + {65458, 16}, {65459, 16}, {65460, 16}, {65461, 16}, {0, 0}, {0, 0}, {0, 0}, + {0, 0}, {0, 0}, {0, 0}, {504, 9}, {32704, 15}, {65462, 16}, {65463, 16}, + {65464, 16}, {65465, 16}, {65466, 16}, {65467, 16}, {65468, 16}, {65469, 16}, {0, 0}, + {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {505, 9}, {65470, 16}, + {65471, 16}, {65472, 16}, {65473, 16}, {65474, 16}, {65475, 16}, {65476, 16}, {65477, 16}, + {65478, 16}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, + {506, 9}, {65479, 16}, {65480, 16}, {65481, 16}, {65482, 16}, {65483, 16}, {65484, 16}, + {65485, 16}, {65486, 16}, {65487, 16}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, + {0, 0}, {0, 0}, {1017, 10}, {65488, 16}, {65489, 16}, {65490, 16}, {65491, 16}, + {65492, 16}, {65493, 16}, {65494, 16}, {65495, 16}, {65496, 16}, {0, 0}, {0, 0}, + {0, 0}, {0, 0}, {0, 0}, {0, 0}, {1018, 10}, {65497, 16}, {65498, 16}, + {65499, 16}, {65500, 16}, {65501, 16}, {65502, 16}, {65503, 16}, {65504, 16}, {65505, 16}, + {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {2040, 11}, + {65506, 16}, {65507, 16}, {65508, 16}, {65509, 16}, {65510, 16}, {65511, 16}, {65512, 16}, + {65513, 16}, {65514, 16}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, + {0, 0}, {65515, 16}, {65516, 16}, {65517, 16}, {65518, 16}, {65519, 16}, {65520, 16}, + {65521, 16}, {65522, 16}, {65523, 16}, {65524, 16}, {0, 0}, {0, 0}, {0, 0}, + {0, 0}, {0, 0}, {2041, 11}, {65525, 16}, {65526, 16}, {65527, 16}, {65528, 16}, + {65529, 16}, {65530, 16}, {65531, 16}, {65532, 16}, {65533, 16}, {65534, 16}, {0, 0}, + {0, 0}, {0, 0}, {0, 0}, {0, 0}}; + static const unsigned short UVAC_HT[256][2] = { + {0, 2}, {1, 2}, {4, 3}, {10, 4}, {24, 5}, {25, 5}, {56, 6}, + {120, 7}, {500, 9}, {1014, 10}, {4084, 12}, {0, 0}, {0, 0}, {0, 0}, + {0, 0}, {0, 0}, {0, 0}, {11, 4}, {57, 6}, {246, 8}, {501, 9}, + {2038, 11}, {4085, 12}, {65416, 16}, {65417, 16}, {65418, 16}, {65419, 16}, {0, 0}, + {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {26, 5}, {247, 8}, + {1015, 10}, {4086, 12}, {32706, 15}, {65420, 16}, {65421, 16}, {65422, 16}, {65423, 16}, + {65424, 16}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, + {27, 5}, {248, 8}, {1016, 10}, {4087, 12}, {65425, 16}, {65426, 16}, {65427, 16}, + {65428, 16}, {65429, 16}, {65430, 16}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, + {0, 0}, {0, 0}, {58, 6}, {502, 9}, {65431, 16}, {65432, 16}, {65433, 16}, + {65434, 16}, {65435, 16}, {65436, 16}, {65437, 16}, {65438, 16}, {0, 0}, {0, 0}, + {0, 0}, {0, 0}, {0, 0}, {0, 0}, {59, 6}, {1017, 10}, {65439, 16}, + {65440, 16}, {65441, 16}, {65442, 16}, {65443, 16}, {65444, 16}, {65445, 16}, {65446, 16}, + {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {121, 7}, + {2039, 11}, {65447, 16}, {65448, 16}, {65449, 16}, {65450, 16}, {65451, 16}, {65452, 16}, + {65453, 16}, {65454, 16}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, + {0, 0}, {122, 7}, {2040, 11}, {65455, 16}, {65456, 16}, {65457, 16}, {65458, 16}, + {65459, 16}, {65460, 16}, {65461, 16}, {65462, 16}, {0, 0}, {0, 0}, {0, 0}, + {0, 0}, {0, 0}, {0, 0}, {249, 8}, {65463, 16}, {65464, 16}, {65465, 16}, + {65466, 16}, {65467, 16}, {65468, 16}, {65469, 16}, {65470, 16}, {65471, 16}, {0, 0}, + {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {503, 9}, {65472, 16}, + {65473, 16}, {65474, 16}, {65475, 16}, {65476, 16}, {65477, 16}, {65478, 16}, {65479, 16}, + {65480, 16}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, + {504, 9}, {65481, 16}, {65482, 16}, {65483, 16}, {65484, 16}, {65485, 16}, {65486, 16}, + {65487, 16}, {65488, 16}, {65489, 16}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, + {0, 0}, {0, 0}, {505, 9}, {65490, 16}, {65491, 16}, {65492, 16}, {65493, 16}, + {65494, 16}, {65495, 16}, {65496, 16}, {65497, 16}, {65498, 16}, {0, 0}, {0, 0}, + {0, 0}, {0, 0}, {0, 0}, {0, 0}, {506, 9}, {65499, 16}, {65500, 16}, + {65501, 16}, {65502, 16}, {65503, 16}, {65504, 16}, {65505, 16}, {65506, 16}, {65507, 16}, + {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {2041, 11}, + {65508, 16}, {65509, 16}, {65510, 16}, {65511, 16}, {65512, 16}, {65513, 16}, {65514, 16}, + {65515, 16}, {65516, 16}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, {0, 0}, + {0, 0}, {16352, 14}, {65517, 16}, {65518, 16}, {65519, 16}, {65520, 16}, {65521, 16}, + {65522, 16}, {65523, 16}, {65524, 16}, {65525, 16}, {0, 0}, {0, 0}, {0, 0}, + {0, 0}, {0, 0}, {1018, 10}, {32707, 15}, {65526, 16}, {65527, 16}, {65528, 16}, + {65529, 16}, {65530, 16}, {65531, 16}, {65532, 16}, {65533, 16}, {65534, 16}, {0, 0}, + {0, 0}, {0, 0}, {0, 0}, {0, 0}}; + static const int YQT[] = {16, 11, 10, 16, 24, 40, 51, 61, 12, 12, 14, 19, 26, 58, 60, 55, + 14, 13, 16, 24, 40, 57, 69, 56, 14, 17, 22, 29, 51, 87, 80, 62, + 18, 22, 37, 56, 68, 109, 103, 77, 24, 35, 55, 64, 81, 104, 113, 92, + 49, 64, 78, 87, 103, 121, 120, 101, 72, 92, 95, 98, 112, 100, 103, 99}; + static const int UVQT[] = {17, 18, 24, 47, 99, 99, 99, 99, 18, 21, 26, 66, 99, 99, 99, 99, + 24, 26, 56, 99, 99, 99, 99, 99, 47, 66, 99, 99, 99, 99, 99, 99, + 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, + 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99}; + static const float aasf[] = {1.0f * 2.828427125f, + 1.387039845f * 2.828427125f, + 1.306562965f * 2.828427125f, + 1.175875602f * 2.828427125f, + 1.0f * 2.828427125f, + 0.785694958f * 2.828427125f, + 0.541196100f * 2.828427125f, + 0.275899379f * 2.828427125f}; + + int row, col, i, k, subsample; + float fdtbl_Y[64], fdtbl_UV[64]; + unsigned char YTable[64], UVTable[64]; + + if (!data || !width || !height || comp > 4 || comp < 1) { + return 0; + } + + quality = quality ? quality : 90; + subsample = quality <= 90 ? 1 : 0; + quality = quality < 1 ? 1 : quality > 100 ? 100 : quality; + quality = quality < 50 ? 5000 / quality : 200 - quality * 2; + + for (i = 0; i < 64; ++i) { + int uvti, yti = (YQT[i] * quality + 50) / 100; + YTable[stbiw__jpg_ZigZag[i]] = (unsigned char)(yti < 1 ? 1 : yti > 255 ? 255 : yti); + uvti = (UVQT[i] * quality + 50) / 100; + UVTable[stbiw__jpg_ZigZag[i]] = (unsigned char)(uvti < 1 ? 1 : uvti > 255 ? 255 : uvti); + } + + for (row = 0, k = 0; row < 8; ++row) { + for (col = 0; col < 8; ++col, ++k) { + fdtbl_Y[k] = 1 / (YTable[stbiw__jpg_ZigZag[k]] * aasf[row] * aasf[col]); + fdtbl_UV[k] = 1 / (UVTable[stbiw__jpg_ZigZag[k]] * aasf[row] * aasf[col]); + } + } + + // Write Headers + { + static const unsigned char head0[] = {0xFF, 0xD8, 0xFF, 0xE0, 0, 0x10, 'J', 'F', 'I', + 'F', 0, 1, 1, 0, 0, 1, 0, 1, + 0, 0, 0xFF, 0xDB, 0, 0x84, 0}; + static const unsigned char head2[] = { + 0xFF, 0xDA, 0, 0xC, 3, 1, 0, 2, 0x11, 3, 0x11, 0, 0x3F, 0}; + const unsigned char head1[] = {0xFF, + 0xC0, + 0, + 0x11, + 8, + (unsigned char)(height >> 8), + STBIW_UCHAR(height), + (unsigned char)(width >> 8), + STBIW_UCHAR(width), + 3, + 1, + (unsigned char)(subsample ? 0x22 : 0x11), + 0, + 2, + 0x11, + 1, + 3, + 0x11, + 1, + 0xFF, + 0xC4, + 0x01, + 0xA2, + 0}; + s->func(s->context, (void*)head0, sizeof(head0)); + s->func(s->context, (void*)YTable, sizeof(YTable)); + stbiw__putc(s, 1); + s->func(s->context, UVTable, sizeof(UVTable)); + s->func(s->context, (void*)head1, sizeof(head1)); + s->func( + s->context, (void*)(std_dc_luminance_nrcodes + 1), sizeof(std_dc_luminance_nrcodes) - 1); + s->func(s->context, (void*)std_dc_luminance_values, sizeof(std_dc_luminance_values)); + stbiw__putc(s, 0x10); // HTYACinfo + s->func( + s->context, (void*)(std_ac_luminance_nrcodes + 1), sizeof(std_ac_luminance_nrcodes) - 1); + s->func(s->context, (void*)std_ac_luminance_values, sizeof(std_ac_luminance_values)); + stbiw__putc(s, 1); // HTUDCinfo + s->func(s->context, + (void*)(std_dc_chrominance_nrcodes + 1), + sizeof(std_dc_chrominance_nrcodes) - 1); + s->func(s->context, (void*)std_dc_chrominance_values, sizeof(std_dc_chrominance_values)); + stbiw__putc(s, 0x11); // HTUACinfo + s->func(s->context, + (void*)(std_ac_chrominance_nrcodes + 1), + sizeof(std_ac_chrominance_nrcodes) - 1); + s->func(s->context, (void*)std_ac_chrominance_values, sizeof(std_ac_chrominance_values)); + s->func(s->context, (void*)head2, sizeof(head2)); + } + + // Encode 8x8 macroblocks + { + static const unsigned short fillBits[] = {0x7F, 7}; + int DCY = 0, DCU = 0, DCV = 0; + int bitBuf = 0, bitCnt = 0; + // comp == 2 is grey+alpha (alpha is ignored) + int ofsG = comp > 2 ? 1 : 0, ofsB = comp > 2 ? 2 : 0; + const unsigned char* dataR = (const unsigned char*)data; + const unsigned char* dataG = dataR + ofsG; + const unsigned char* dataB = dataR + ofsB; + int x, y, pos; + if (subsample) { + for (y = 0; y < height; y += 16) { + for (x = 0; x < width; x += 16) { + float Y[256], U[256], V[256]; + for (row = y, pos = 0; row < y + 16; ++row) { + // row >= height => use last input row + int clamped_row = (row < height) ? row : height - 1; + int base_p = + (stbi__flip_vertically_on_write ? (height - 1 - clamped_row) : clamped_row) * + width * comp; + for (col = x; col < x + 16; ++col, ++pos) { + // if col >= width => use pixel from last input column + int p = base_p + ((col < width) ? col : (width - 1)) * comp; + float r = dataR[p], g = dataG[p], b = dataB[p]; + Y[pos] = +0.29900f * r + 0.58700f * g + 0.11400f * b - 128; + U[pos] = -0.16874f * r - 0.33126f * g + 0.50000f * b; + V[pos] = +0.50000f * r - 0.41869f * g - 0.08131f * b; } - } - } else { - for(y = 0; y < height; y += 8) { - for(x = 0; x < width; x += 8) { - float Y[64], U[64], V[64]; - for(row = y, pos = 0; row < y+8; ++row) { - // row >= height => use last input row - int clamped_row = (row < height) ? row : height - 1; - int base_p = (stbi__flip_vertically_on_write ? (height-1-clamped_row) : clamped_row)*width*comp; - for(col = x; col < x+8; ++col, ++pos) { - // if col >= width => use pixel from last input column - int p = base_p + ((col < width) ? col : (width-1))*comp; - float r = dataR[p], g = dataG[p], b = dataB[p]; - Y[pos]= +0.29900f*r + 0.58700f*g + 0.11400f*b - 128; - U[pos]= -0.16874f*r - 0.33126f*g + 0.50000f*b; - V[pos]= +0.50000f*r - 0.41869f*g - 0.08131f*b; - } - } - - DCY = stbiw__jpg_processDU(s, &bitBuf, &bitCnt, Y, 8, fdtbl_Y, DCY, YDC_HT, YAC_HT); - DCU = stbiw__jpg_processDU(s, &bitBuf, &bitCnt, U, 8, fdtbl_UV, DCU, UVDC_HT, UVAC_HT); - DCV = stbiw__jpg_processDU(s, &bitBuf, &bitCnt, V, 8, fdtbl_UV, DCV, UVDC_HT, UVAC_HT); + } + DCY = stbiw__jpg_processDU(s, &bitBuf, &bitCnt, Y + 0, 16, fdtbl_Y, DCY, YDC_HT, YAC_HT); + DCY = stbiw__jpg_processDU(s, &bitBuf, &bitCnt, Y + 8, 16, fdtbl_Y, DCY, YDC_HT, YAC_HT); + DCY = + stbiw__jpg_processDU(s, &bitBuf, &bitCnt, Y + 128, 16, fdtbl_Y, DCY, YDC_HT, YAC_HT); + DCY = + stbiw__jpg_processDU(s, &bitBuf, &bitCnt, Y + 136, 16, fdtbl_Y, DCY, YDC_HT, YAC_HT); + + // subsample U,V + { + float subU[64], subV[64]; + int yy, xx; + for (yy = 0, pos = 0; yy < 8; ++yy) { + for (xx = 0; xx < 8; ++xx, ++pos) { + int j = yy * 32 + xx * 2; + subU[pos] = (U[j + 0] + U[j + 1] + U[j + 16] + U[j + 17]) * 0.25f; + subV[pos] = (V[j + 0] + V[j + 1] + V[j + 16] + V[j + 17]) * 0.25f; + } } - } + DCU = + stbiw__jpg_processDU(s, &bitBuf, &bitCnt, subU, 8, fdtbl_UV, DCU, UVDC_HT, UVAC_HT); + DCV = + stbiw__jpg_processDU(s, &bitBuf, &bitCnt, subV, 8, fdtbl_UV, DCV, UVDC_HT, UVAC_HT); + } + } } + } else { + for (y = 0; y < height; y += 8) { + for (x = 0; x < width; x += 8) { + float Y[64], U[64], V[64]; + for (row = y, pos = 0; row < y + 8; ++row) { + // row >= height => use last input row + int clamped_row = (row < height) ? row : height - 1; + int base_p = + (stbi__flip_vertically_on_write ? (height - 1 - clamped_row) : clamped_row) * + width * comp; + for (col = x; col < x + 8; ++col, ++pos) { + // if col >= width => use pixel from last input column + int p = base_p + ((col < width) ? col : (width - 1)) * comp; + float r = dataR[p], g = dataG[p], b = dataB[p]; + Y[pos] = +0.29900f * r + 0.58700f * g + 0.11400f * b - 128; + U[pos] = -0.16874f * r - 0.33126f * g + 0.50000f * b; + V[pos] = +0.50000f * r - 0.41869f * g - 0.08131f * b; + } + } - // Do the bit alignment of the EOI marker - stbiw__jpg_writeBits(s, &bitBuf, &bitCnt, fillBits); - } + DCY = stbiw__jpg_processDU(s, &bitBuf, &bitCnt, Y, 8, fdtbl_Y, DCY, YDC_HT, YAC_HT); + DCU = stbiw__jpg_processDU(s, &bitBuf, &bitCnt, U, 8, fdtbl_UV, DCU, UVDC_HT, UVAC_HT); + DCV = stbiw__jpg_processDU(s, &bitBuf, &bitCnt, V, 8, fdtbl_UV, DCV, UVDC_HT, UVAC_HT); + } + } + } - // EOI - stbiw__putc(s, 0xFF); - stbiw__putc(s, 0xD9); + // Do the bit alignment of the EOI marker + stbiw__jpg_writeBits(s, &bitBuf, &bitCnt, fillBits); + } - return 1; -} + // EOI + stbiw__putc(s, 0xFF); + stbiw__putc(s, 0xD9); -STBIWDEF int stbi_write_jpg_to_func(stbi_write_func *func, void *context, int x, int y, int comp, const void *data, int quality) -{ - stbi__write_context s = { 0 }; - stbi__start_write_callbacks(&s, func, context); - return stbi_write_jpg_core(&s, x, y, comp, (void *) data, quality); + return 1; } +STBIWDEF int stbi_write_jpg_to_func( + stbi_write_func* func, void* context, int x, int y, int comp, const void* data, int quality) { + stbi__write_context s = {0}; + stbi__start_write_callbacks(&s, func, context); + return stbi_write_jpg_core(&s, x, y, comp, (void*)data, quality); +} #ifndef STBI_WRITE_NO_STDIO -STBIWDEF int stbi_write_jpg(char const *filename, int x, int y, int comp, const void *data, int quality) -{ - stbi__write_context s = { 0 }; - if (stbi__start_write_file(&s,filename)) { - int r = stbi_write_jpg_core(&s, x, y, comp, data, quality); - stbi__end_write_file(&s); - return r; - } else - return 0; +STBIWDEF int + stbi_write_jpg(const char* filename, int x, int y, int comp, const void* data, int quality) { + stbi__write_context s = {0}; + if (stbi__start_write_file(&s, filename)) { + int r = stbi_write_jpg_core(&s, x, y, comp, data, quality); + stbi__end_write_file(&s); + return r; + } else + return 0; } #endif @@ -1608,24 +1963,13 @@ STBIWDEF int stbi_write_jpg(char const *filename, int x, int y, int comp, const 1.09 (2018-02-11) fix typo in zlib quality API, improve STB_I_W_STATIC in C++ 1.08 (2018-01-29) - add stbi__flip_vertically_on_write, external zlib, zlib quality, choose PNG filter - 1.07 (2017-07-24) - doc fix - 1.06 (2017-07-23) - writing JPEG (using Jon Olick's code) - 1.05 ??? - 1.04 (2017-03-03) - monochrome BMP expansion - 1.03 ??? - 1.02 (2016-04-02) - avoid allocating large structures on the stack - 1.01 (2016-01-16) - STBIW_REALLOC_SIZED: support allocators with no realloc support - avoid race-condition in crc initialization - minor compile issues - 1.00 (2015-09-14) - installable file IO function - 0.99 (2015-09-13) + add stbi__flip_vertically_on_write, external zlib, zlib quality, + choose PNG filter 1.07 (2017-07-24) doc fix 1.06 (2017-07-23) writing JPEG + (using Jon Olick's code) 1.05 ??? 1.04 (2017-03-03) monochrome BMP + expansion 1.03 ??? 1.02 (2016-04-02) avoid allocating large structures on + the stack 1.01 (2016-01-16) STBIW_REALLOC_SIZED: support allocators with no + realloc support avoid race-condition in crc initialization minor compile + issues 1.00 (2015-09-14) installable file IO function 0.99 (2015-09-13) warning fixes; TGA rle support 0.98 (2015-04-08) added STBIW_MALLOC, STBIW_ASSERT etc @@ -1635,7 +1979,7 @@ STBIWDEF int stbi_write_jpg(char const *filename, int x, int y, int comp, const add HDR output fix monochrome BMP 0.95 (2014-08-17) - add monochrome TGA output + add monochrome TGA output 0.94 (2014-05-31) rename private functions to avoid conflicts with stb_image.h 0.93 (2014-05-27) diff --git a/seissolxdmf/README.md b/seissolxdmf/README.md index cabd7c5..bf34e10 100644 --- a/seissolxdmf/README.md +++ b/seissolxdmf/README.md @@ -1,7 +1,9 @@ seissolxdmf =============== + A python reader for SeisSol xdmf output (posix or hdf5) and hdf5 meshes. Below is an simple example, illustrating the use of the module: + ```python import seissolxdmf fn = 'test-fault.xdmf' diff --git a/seissolxdmf/seissolxdmf/seissolxdmf.py b/seissolxdmf/seissolxdmf/seissolxdmf.py index f180f83..f89f5f0 100644 --- a/seissolxdmf/seissolxdmf/seissolxdmf.py +++ b/seissolxdmf/seissolxdmf/seissolxdmf.py @@ -55,7 +55,7 @@ def GetDtype(self, data_prec, isInt): def ReadSimpleBinaryFile(self, absolute_path, MemDimension, data_prec, isInt, idt=-1): """Read block of data in binary format (posix) - idt!=-1 loads only one time step + idt!=-1 loads only one time step this function is a special case of ReadSimpleBinaryFileChunk but kept for performance reasons """ oneDtMem = True if idt != -1 else False diff --git a/seissolxdmfwriter/README.md b/seissolxdmfwriter/README.md index 4c6903d..7dd387a 100644 --- a/seissolxdmfwriter/README.md +++ b/seissolxdmfwriter/README.md @@ -1,5 +1,6 @@ seissolxdmfwriter ================= + A python module to write seissol outputs readable by paraview. ```python @@ -8,7 +9,9 @@ import seissolxdmf as sx import numpy as np fn = "test-fault.xdmf" + # Read data from input file using seissolxdmf + sx = sx.seissolxdmf(fn) geom = sx.ReadGeometry() connect = sx.ReadConnect() @@ -31,8 +34,8 @@ sxw.write( backend="hdf5", ) -# Finally, the module can be use to write data directly from seissolxdmf, limiting -# the memory requirements +# Finally, the module can be use to write data directly from seissolxdmf, +# limiting the memory requirements sxw.write_from_seissol_output( 'test-fault-sx', @@ -46,11 +49,22 @@ sxw.write_from_seissol_output( ``` -The module also encapsulates `seissol_output_extractor`, which can be used to extract and process data from SeisSol output files, allowing selection of variables, time steps, spatial ranges, and output format. +The module also encapsulates `seissol_output_extractor`, +which can be used to extract and process data from SeisSol +output files, allowing selection of variables, time steps, +spatial ranges, and output format. Here is an example of use: ```bash -# extracts PSR, Vr and partition, at 2nd and 4th time steps and at simulation time 0.5, from test-fault.xdmf and write into test_new-fault.xdmf -# use seissol_output_extractor --h for additionnal info about the arguments -seissol_output_extractor test-fault.xdmf --time "i2,i4,0.5" --variable PSR Vr partition --add2prefix "_new" +# extracts PSR, Vr and partition, at 2nd and 4th time +# steps and at simulation time 0.5, from test-fault.xdmf +# and write into test_new-fault.xdmf +# +# use seissol_output_extractor --h for additionnal +# info about the arguments + +seissol_output_extractor test-fault.xdmf \ + --time "i2,i4,0.5" \ + --variable PSR Vr partition \ + --add2prefix "_new" ``` diff --git a/tag_faultside/SConstruct b/tag_faultside/SConstruct index 1eb9680..4c3b6c1 100644 --- a/tag_faultside/SConstruct +++ b/tag_faultside/SConstruct @@ -1,4 +1,3 @@ -#!/usr/bin/python ## # @file # This file is part of SeisSol. @@ -76,4 +75,4 @@ Export('env') SConscript('src/SConscript', variant_dir='#/build/'+buildDir, src_dir='#/', duplicate=0) Import('env') -env.Program('#/build/tag-faultside', env.sourceFiles) \ No newline at end of file +env.Program('#/build/tag-faultside', env.sourceFiles) diff --git a/tag_faultside/src/KDTree.cpp b/tag_faultside/src/KDTree.cpp index ae90dfc..50bd11f 100644 --- a/tag_faultside/src/KDTree.cpp +++ b/tag_faultside/src/KDTree.cpp @@ -2,8 +2,10 @@ * @file * This file is part of SeisSol. * - * @author Carsten Uphoff (c.uphoff AT tum.de, http://www5.in.tum.de/wiki/index.php/Carsten_Uphoff,_M.Sc.) - * @author Sebastian Rettenberger (sebastian.rettenberger AT tum.de, http://www5.in.tum.de/wiki/index.php/Sebastian_Rettenberger) + * @author Carsten Uphoff (c.uphoff AT tum.de, + * http://www5.in.tum.de/wiki/index.php/Carsten_Uphoff,_M.Sc.) + * @author Sebastian Rettenberger (sebastian.rettenberger AT tum.de, + * http://www5.in.tum.de/wiki/index.php/Sebastian_Rettenberger) * * @section LICENSE * Copyright (c) 2016, SeisSol Group @@ -41,108 +43,102 @@ #include "KDTree.h" #include -#include #include +#include -KDTree::KDTree(std::unordered_set const& points, int maxLeafSize, bool split[3]) - : maxLeafN(maxLeafSize) -{ - memcpy(this->split, split, sizeof(bool)*3); +KDTree::KDTree(const std::unordered_set& points, int maxLeafSize, bool split[3]) + : maxLeafN(maxLeafSize) { + memcpy(this->split, split, sizeof(bool) * 3); - int n = points.size(); + int n = points.size(); - // Copy data and change row-major to column-major storage. - data = new Point[n]; - Point* next = data; - for (std::unordered_set::const_iterator it = points.begin(); - it != points.end(); it++) { - memcpy(next->coords, it->coords, sizeof(double)*3); - next++; - } + // Copy data and change row-major to column-major storage. + data = new Point[n]; + Point* next = data; + for (std::unordered_set::const_iterator it = points.begin(); it != points.end(); it++) { + memcpy(next->coords, it->coords, sizeof(double) * 3); + next++; + } - idx = new int[n]; - for (int i = 0; i < n; ++i) { - idx[i] = i; - } + idx = new int[n]; + for (int i = 0; i < n; ++i) { + idx[i] = i; + } - int maxHeight = 1 + ceil(log2(n / static_cast(maxLeafN))); - int maxNodes = (1 << maxHeight) - 1; - nodes = new Node[maxNodes]; - nodes[0].start = 0; - nodes[0].n = n; + int maxHeight = 1 + ceil(log2(n / static_cast(maxLeafN))); + int maxNodes = (1 << maxHeight) - 1; + nodes = new Node[maxNodes]; + nodes[0].start = 0; + nodes[0].n = n; - buildTree(0, 0); + buildTree(0, 0); } -KDTree::~KDTree() -{ - delete[] data; - delete[] idx; - delete[] nodes; +KDTree::~KDTree() { + delete[] data; + delete[] idx; + delete[] nodes; } -void KDTree::swap(int i, int j) -{ - if (i != j) { - std::swap(idx[i], idx[j]); +void KDTree::swap(int i, int j) { + if (i != j) { + std::swap(idx[i], idx[j]); std::swap(data[i], data[j]); - } + } } -int KDTree::partition(int left, int right, int pivotIdx, int splitdim) -{ - double pivot = data[pivotIdx].coords[splitdim]; - int st = left; - for (int i = left; i < right; ++i) { - if (data[i].coords[splitdim] < pivot) { - swap(st, i); - ++st; - } - } - swap(st, right); - return st; +int KDTree::partition(int left, int right, int pivotIdx, int splitdim) { + double pivot = data[pivotIdx].coords[splitdim]; + int st = left; + for (int i = left; i < right; ++i) { + if (data[i].coords[splitdim] < pivot) { + swap(st, i); + ++st; + } + } + swap(st, right); + return st; } -void KDTree::buildTree(int k, int splitdim) -{ - Node& node = nodes[k]; - node.splitdim = splitdim; - if (node.n > maxLeafN) { - int half = (node.n % 2 != 0) ? (node.n + 1)/2 : node.n/2; - int l = node.start; - int r = l + node.n - 1; - int median_idx = l + half; - if (l != r) { - int pivotIdx; - while (true) { - pivotIdx = r; - pivotIdx = partition(l, r, pivotIdx, splitdim); - if (median_idx == pivotIdx) { - break; - } else if (median_idx < pivotIdx) { - r = pivotIdx - 1; - } else { - l = pivotIdx + 1; - } - } - } - node.pivot = data[median_idx].coords[splitdim]; +void KDTree::buildTree(int k, int splitdim) { + Node& node = nodes[k]; + node.splitdim = splitdim; + if (node.n > maxLeafN) { + int half = (node.n % 2 != 0) ? (node.n + 1) / 2 : node.n / 2; + int l = node.start; + int r = l + node.n - 1; + int median_idx = l + half; + if (l != r) { + int pivotIdx; + while (true) { + pivotIdx = r; + pivotIdx = partition(l, r, pivotIdx, splitdim); + if (median_idx == pivotIdx) { + break; + } else if (median_idx < pivotIdx) { + r = pivotIdx - 1; + } else { + l = pivotIdx + 1; + } + } + } + node.pivot = data[median_idx].coords[splitdim]; - Node& left = nodes[leftChild(k)]; - left.start = node.start; - left.n = half; + Node& left = nodes[leftChild(k)]; + left.start = node.start; + left.n = half; - Node& right = nodes[rightChild(k)]; - right.start = median_idx; - right.n = node.n - half; + Node& right = nodes[rightChild(k)]; + right.start = median_idx; + right.n = node.n - half; - int nextSplitdim = splitdim; - do { - nextSplitdim = (nextSplitdim + 1) % 3; - } while (!split[nextSplitdim]); - buildTree(leftChild(k), nextSplitdim); - buildTree(rightChild(k), nextSplitdim); - } else { - node.isLeaf = true; - } + int nextSplitdim = splitdim; + do { + nextSplitdim = (nextSplitdim + 1) % 3; + } while (!split[nextSplitdim]); + buildTree(leftChild(k), nextSplitdim); + buildTree(rightChild(k), nextSplitdim); + } else { + node.isLeaf = true; + } } diff --git a/tag_faultside/src/KDTree.h b/tag_faultside/src/KDTree.h index 528f8f9..75c7c18 100644 --- a/tag_faultside/src/KDTree.h +++ b/tag_faultside/src/KDTree.h @@ -2,8 +2,10 @@ * @file * This file is part of SeisSol. * - * @author Carsten Uphoff (c.uphoff AT tum.de, http://www5.in.tum.de/wiki/index.php/Carsten_Uphoff,_M.Sc.) - * @author Sebastian Rettenberger (sebastian.rettenberger AT tum.de, http://www5.in.tum.de/wiki/index.php/Sebastian_Rettenberger) + * @author Carsten Uphoff (c.uphoff AT tum.de, + * http://www5.in.tum.de/wiki/index.php/Carsten_Uphoff,_M.Sc.) + * @author Sebastian Rettenberger (sebastian.rettenberger AT tum.de, + * http://www5.in.tum.de/wiki/index.php/Sebastian_Rettenberger) * * @section LICENSE * Copyright (c) 2016, SeisSol Group @@ -53,107 +55,97 @@ union Point { double z; }; - bool operator==(const Point& other) const { - return x == other.x & y == other.y & z == other.z; - } + bool operator==(const Point& other) const { return x == other.x & y == other.y & z == other.z; } }; // a little helper that should IMHO be standardized -template -std::size_t make_hash(const T& v) -{ - return std::hash()(v); +template +std::size_t make_hash(const T& v) { + return std::hash()(v); } // adapted from boost::hash_combine -inline void hash_combine(std::size_t& h, const std::size_t& v) -{ - h ^= v + 0x9e3779b9 + (h << 6) + (h >> 2); +inline void hash_combine(std::size_t& h, const std::size_t& v) { + h ^= v + 0x9e3779b9 + (h << 6) + (h >> 2); } -namespace std -{ +namespace std { // support for Point -template<> -struct hash -{ - size_t operator()(const Point& v) const - { - size_t x = make_hash(v.x); - size_t y = make_hash(v.y); - size_t z = make_hash(v.z); - hash_combine(x, y); - hash_combine(x, z); - return x; - } +template <> +struct hash { + size_t operator()(const Point& v) const { + size_t x = make_hash(v.x); + size_t y = make_hash(v.y); + size_t z = make_hash(v.z); + hash_combine(x, y); + hash_combine(x, z); + return x; + } }; -} +} // namespace std class KDTree { -public: - /** - * @param split The dimensions that should be used for spliting. - * At least one of the 3 values has to be true. - */ - KDTree(std::unordered_set const& points, int maxLeafSize, - bool split[3]); - ~KDTree(); - - template - void search(Support const& support, Action& action) const - { - searchTree(0, support, action); - } - - inline Point const* points() const { return data; } - inline int index(int r) const { return idx[r]; } - inline int numPoints() const { return nodes[0].n;} - -private: - int leftChild(int k) const { return 2*k + 1; } - int rightChild(int k) const { return 2*k + 2; } - - struct Node { + public: + /** + * @param split The dimensions that should be used for spliting. + * At least one of the 3 values has to be true. + */ + KDTree(const std::unordered_set& points, int maxLeafSize, bool split[3]); + ~KDTree(); + + template + void search(const Support& support, Action& action) const { + searchTree(0, support, action); + } + + inline const Point* points() const { return data; } + inline int index(int r) const { return idx[r]; } + inline int numPoints() const { return nodes[0].n; } + + private: + int leftChild(int k) const { return 2 * k + 1; } + int rightChild(int k) const { return 2 * k + 2; } + + struct Node { Node() : isLeaf(false) {} - double pivot; - int start; - int n; - int splitdim; - bool isLeaf; - }; - Node* nodes; - - void swap(int i, int j); - int partition(int left, int right, int pivotIdx, int splitdim); - void buildTree(int k, int splitdim); - - template - void searchTree(int k, Support const& support, Action& action) const; - - Point* data; - int* idx; - bool split[3]; - int maxLeafN; + double pivot; + int start; + int n; + int splitdim; + bool isLeaf; + }; + Node* nodes; + + void swap(int i, int j); + int partition(int left, int right, int pivotIdx, int splitdim); + void buildTree(int k, int splitdim); + + template + void searchTree(int k, const Support& support, Action& action) const; + + Point* data; + int* idx; + bool split[3]; + int maxLeafN; }; -template -void KDTree::searchTree(int k, Support const& support, Action& action) const -{ - Node& node = nodes[k]; - if (node.isLeaf) { - for (int i = node.start; i < node.start + node.n; ++i) { +template +void KDTree::searchTree(int k, const Support& support, Action& action) const { + Node& node = nodes[k]; + if (node.isLeaf) { + for (int i = node.start; i < node.start + node.n; ++i) { action(data[i]); - } - } else { - if (support(node.splitdim, 0) <= node.pivot) { - searchTree(leftChild(k), support, action); - } - if (support(node.splitdim, 1) >= node.pivot) { - searchTree(rightChild(k), support, action); - } - } + } + } else { + if (support(node.splitdim, 0) <= node.pivot) { + searchTree(leftChild(k), support, action); + } + if (support(node.splitdim, 1) >= node.pivot) { + searchTree(rightChild(k), support, action); + } + } } #endif diff --git a/tag_faultside/src/SConscript b/tag_faultside/src/SConscript index b3ddfad..7593cf2 100644 --- a/tag_faultside/src/SConscript +++ b/tag_faultside/src/SConscript @@ -1,4 +1,3 @@ -#!/usr/bin/python ## # @file # This file is part of SeisSol. @@ -44,4 +43,4 @@ sourceFiles = ['main.cpp', 'KDTree.cpp'] for f in sourceFiles: env.sourceFiles.append(env.Object(f)) -Export('env') \ No newline at end of file +Export('env') diff --git a/tag_faultside/src/main.cpp b/tag_faultside/src/main.cpp index 844fa43..d286ca9 100644 --- a/tag_faultside/src/main.cpp +++ b/tag_faultside/src/main.cpp @@ -2,7 +2,8 @@ * @file * This file is part of SeisSol. * - * @author Sebastian Rettenberger (sebastian.rettenberger AT tum.de, http://www5.in.tum.de/wiki/index.php/Sebastian_Rettenberger) + * @author Sebastian Rettenberger (sebastian.rettenberger AT tum.de, + * http://www5.in.tum.de/wiki/index.php/Sebastian_Rettenberger) * * @section LICENSE * Copyright (c) 2016, SeisSol Group @@ -37,34 +38,29 @@ * @section DESCRIPTION */ +#include "KDTree.h" +#include "utils/args.h" +#include "utils/logger.h" + #include #include +#include #include +#include #include #include -#include - -#include - -#include "utils/args.h" -#include "utils/logger.h" - -#include "KDTree.h" - const static int FACE2NODES[4][3] = {{0, 2, 1}, {0, 1, 3}, {0, 3, 2}, {1, 2, 3}}; -template -static void checkH5Err(T status) -{ - if (status < 0) - logError() << "An HDF5 error occurred"; +template +static void checkH5Err(T status) { + if (status < 0) + logError() << "An HDF5 error occurred"; } -static void checkNcError(int error) -{ - if (error != NC_NOERR) - logError() << "An netCDF error occurred:" << nc_strerror(error); +static void checkNcError(int error) { + if (error != NC_NOERR) + logError() << "An netCDF error occurred:" << nc_strerror(error); } struct Support { @@ -82,272 +78,279 @@ struct Support { limits[2][1] = -std::numeric_limits::infinity(); } - double operator()(int splitdim, int side) const { - return limits[splitdim][side]; - } + double operator()(int splitdim, int side) const { return limits[splitdim][side]; } }; struct Action { - std::vector points; + std::vector points; - double operator()(Point& p) { - points.push_back(p); - } + double operator()(Point& p) { points.push_back(p); } }; -int main(int argc, char* argv[]) -{ - utils::Args args; -// args.addOption("level", 'l', "gzip compressen level [0-9]", utils::Args::Required, false); - args.addAdditionalOption("mesh", "mesh file"); - args.addAdditionalOption("output", "output file"); +int main(int argc, char* argv[]) { + utils::Args args; + // args.addOption("level", 'l', "gzip compressen level [0-9]", + // utils::Args::Required, false); + args.addAdditionalOption("mesh", "mesh file"); + args.addAdditionalOption("output", "output file"); + + switch (args.parse(argc, argv)) { + case utils::Args::Help: + return 0; + case utils::Args::Error: + return 1; + } + + const unsigned int direction = 0; + + std::string mesh = args.getAdditionalArgument("mesh"); + std::string output = args.getAdditionalArgument("output"); + + // Read the fault vertices + logInfo() << "Build fault K-D-tree..."; + int ncFile; + checkNcError(nc_open(mesh.c_str(), NC_NOWRITE, &ncFile)); + + int ncDimPart; + checkNcError(nc_inq_dimid(ncFile, "partitions", &ncDimPart)); + size_t partitions; + checkNcError(nc_inq_dimlen(ncFile, ncDimPart, &partitions)); + + int ncDimElem; + checkNcError(nc_inq_dimid(ncFile, "elements", &ncDimElem)); + size_t maxElements; + checkNcError(nc_inq_dimlen(ncFile, ncDimElem, &maxElements)); - switch (args.parse(argc, argv)) { - case utils::Args::Help: - return 0; - case utils::Args::Error: - return 1; - } + int ncDimVrtx; + checkNcError(nc_inq_dimid(ncFile, "vertices", &ncDimVrtx)); + size_t maxVertices; + checkNcError(nc_inq_dimlen(ncFile, ncDimVrtx, &maxVertices)); - const unsigned int direction = 0; + int ncElemSize; + checkNcError(nc_inq_varid(ncFile, "element_size", &ncElemSize)); - std::string mesh = args.getAdditionalArgument("mesh"); - std::string output = args.getAdditionalArgument("output"); + int ncElemVertices; + checkNcError(nc_inq_varid(ncFile, "element_vertices", &ncElemVertices)); - // Read the fault vertices - logInfo() << "Build fault K-D-tree..."; - int ncFile; - checkNcError(nc_open(mesh.c_str(), NC_NOWRITE, &ncFile)); + int ncElemBoundaries; + checkNcError(nc_inq_varid(ncFile, "element_boundaries", &ncElemBoundaries)); - int ncDimPart; - checkNcError(nc_inq_dimid(ncFile, "partitions", &ncDimPart)); - size_t partitions; - checkNcError(nc_inq_dimlen(ncFile, ncDimPart, &partitions)); + int ncVrtxSize; + checkNcError(nc_inq_varid(ncFile, "vertex_size", &ncVrtxSize)); - int ncDimElem; - checkNcError(nc_inq_dimid(ncFile, "elements", &ncDimElem)); - size_t maxElements; - checkNcError(nc_inq_dimlen(ncFile, ncDimElem, &maxElements)); + int ncVrtxCoords; + checkNcError(nc_inq_varid(ncFile, "vertex_coordinates", &ncVrtxCoords)); - int ncDimVrtx; - checkNcError(nc_inq_dimid(ncFile, "vertices", &ncDimVrtx)); - size_t maxVertices; - checkNcError(nc_inq_dimlen(ncFile, ncDimVrtx, &maxVertices)); + int* elementSize = new int[partitions]; + checkNcError(nc_get_var_int(ncFile, ncElemSize, elementSize)); - int ncElemSize; - checkNcError(nc_inq_varid(ncFile, "element_size", &ncElemSize)); + int* vertexSize = new int[partitions]; + checkNcError(nc_get_var_int(ncFile, ncVrtxSize, vertexSize)); - int ncElemVertices; - checkNcError(nc_inq_varid(ncFile, "element_vertices", &ncElemVertices)); + int* elementVertices = new int[maxElements * 4]; + int* elementBoundaries = new int[maxElements * 4]; - int ncElemBoundaries; - checkNcError(nc_inq_varid(ncFile, "element_boundaries", &ncElemBoundaries)); + double* vertexCoordinates = new double[maxVertices * 3]; - int ncVrtxSize; - checkNcError(nc_inq_varid(ncFile, "vertex_size", &ncVrtxSize)); + std::unordered_set faultPoints; - int ncVrtxCoords; - checkNcError(nc_inq_varid(ncFile, "vertex_coordinates", &ncVrtxCoords)); + unsigned long totalElements = 0; - int* elementSize = new int[partitions]; - checkNcError(nc_get_var_int(ncFile, ncElemSize, elementSize)); + double minfault = std::numeric_limits::infinity(); + double maxfault = -std::numeric_limits::infinity(); - int* vertexSize = new int[partitions]; - checkNcError(nc_get_var_int(ncFile, ncVrtxSize, vertexSize)); + double minfault2 = std::numeric_limits::infinity(); + double maxfault2 = -std::numeric_limits::infinity(); - int* elementVertices = new int[maxElements*4]; - int* elementBoundaries = new int[maxElements*4]; + double outsideborder1 = 0; + double outsideborder2 = 0; - double* vertexCoordinates = new double[maxVertices*3]; + for (size_t p = 0; p < partitions; p++) { + totalElements += elementSize[p]; - std::unordered_set faultPoints; + size_t offset[3] = {p, 0, 0}; + size_t size[3] = {1, static_cast(elementSize[p]), 4}; - unsigned long totalElements = 0; + checkNcError(nc_get_vara_int(ncFile, ncElemVertices, offset, size, elementVertices)); + checkNcError(nc_get_vara_int(ncFile, ncElemBoundaries, offset, size, elementBoundaries)); - double minfault = std::numeric_limits::infinity(); - double maxfault = -std::numeric_limits::infinity(); + size[1] = vertexSize[p]; + size[2] = 3; + // logInfo() << "size" << size[0] << size[1] << size[2]; + checkNcError(nc_get_vara_double(ncFile, ncVrtxCoords, offset, size, vertexCoordinates)); - double minfault2 = std::numeric_limits::infinity(); - double maxfault2 = -std::numeric_limits::infinity(); + for (unsigned int i = 0; i < elementSize[p]; i++) { + for (unsigned int j = 0; j < 4; j++) { + if (elementBoundaries[i * 4 + j] == 3) { + // Is a fault boundary + for (unsigned int k = 0; k < 3; k++) { + Point point; + memcpy(point.coords, + &vertexCoordinates[elementVertices[i * 4 + FACE2NODES[j][k]] * 3], + sizeof(double) * 3); + faultPoints.insert(point); + + minfault = std::min(minfault, point.coords[direction]); + maxfault = std::max(maxfault, point.coords[direction]); + + if (point.coords[1 - direction] < minfault2) { + minfault2 = point.coords[1 - direction]; + outsideborder1 = point.coords[direction]; + } + if (point.coords[1 - direction] > maxfault2) { + maxfault2 = point.coords[1 - direction]; + outsideborder2 = point.coords[direction]; + } + } + } + } + } + } + + logInfo() << "Cells <" << minfault2 << "and >" << maxfault2 << "have the be handled manually"; + logInfo() << "A good choice might be" << outsideborder1 << "reps." << outsideborder2; + + bool split[3] = {true, true, true}; + split[direction] = false; + KDTree kdtree(faultPoints, 4, split); + + faultPoints.clear(); // Free memory + + unsigned int* isLeft = new unsigned int[totalElements]; + std::fill(isLeft, isLeft + totalElements, 1); + + unsigned long globElement = 0; + for (size_t p = 0; p < partitions; p++) { + logInfo() << "Processing partition" << utils::nospace << p << "..."; + + size_t offset[3] = {p, 0, 0}; + size_t size[3] = {1, static_cast(elementSize[p]), 4}; + + checkNcError(nc_get_vara_int(ncFile, ncElemVertices, offset, size, elementVertices)); + + size[1] = vertexSize[p]; + size[2] = 3; + checkNcError(nc_get_vara_double(ncFile, ncVrtxCoords, offset, size, vertexCoordinates)); + + for (unsigned int i = 0; i < elementSize[p]; i++) { + double avg = 0; + + Support sup; + for (unsigned int j = 0; j < 4; j++) { + sup.limits[0][0] = + std::min(sup.limits[0][0], vertexCoordinates[elementVertices[i * 4 + j] * 3]); + sup.limits[0][1] = + std::max(sup.limits[0][1], vertexCoordinates[elementVertices[i * 4 + j] * 3]); + sup.limits[1][0] = + std::min(sup.limits[1][0], vertexCoordinates[elementVertices[i * 4 + j] * 3 + 1]); + sup.limits[1][1] = + std::max(sup.limits[1][1], vertexCoordinates[elementVertices[i * 4 + j] * 3 + 1]); + sup.limits[2][0] = + std::min(sup.limits[1][0], vertexCoordinates[elementVertices[i * 4 + j] * 3 + 2]); + sup.limits[2][1] = + std::max(sup.limits[1][1], vertexCoordinates[elementVertices[i * 4 + j] * 3 + 2]); + + avg += vertexCoordinates[elementVertices[i * 4 + j] * 3 + direction]; + } + + avg /= 4; + + if (avg < minfault) { + // Do nothing + } else if (avg > maxfault) { + isLeft[globElement] = 0; + } else { + Action act; + kdtree.search(sup, act); + + bool hasPoint = false; + + for (std::vector::const_iterator it = act.points.begin(); it != act.points.end(); + it++) { + if (it->coords[1 - direction] >= sup.limits[1 - direction][0] && + it->coords[1 - direction] <= sup.limits[1 - direction][1] && + it->z >= sup.limits[2][0] && it->z <= sup.limits[2][1]) { + hasPoint = true; + + if (it->coords[direction] < avg) { + isLeft[globElement] = 0; + break; + } + } + } + + if (!hasPoint) { + if (avg > minfault) { + isLeft[globElement] = 0; + } + } + } + + globElement++; + } + } + + checkNcError(nc_close(ncFile)); + + delete[] elementVertices; + delete[] elementBoundaries; + delete[] vertexCoordinates; + delete[] elementSize; + delete[] vertexSize; + + hid_t h5File = H5Fopen(output.c_str(), H5F_ACC_RDWR, H5P_DEFAULT); + checkH5Err(h5File); + + hsize_t chunkSize = 256 * 1024 * 1024 / sizeof(unsigned int); + + hid_t h5var = H5Dopen(h5File, "/is_left", H5P_DEFAULT); + hid_t h5space; + if (h5var >= 0) { + logWarning() << "Overwriting old dataset from HDF5 file"; + + h5space = H5Dget_space(h5var); + checkH5Err(h5space); + int ndims = H5Sget_simple_extent_ndims(h5space); + checkH5Err(ndims); + if (ndims != 1) + logError() << "Old dataset has wrong dimension"; + + hsize_t extent; + checkH5Err(H5Sget_simple_extent_dims(h5space, &extent, 0L)); + if (extent != totalElements) + logError() << "Old dataset has wrong size"; + } else { + // Create new dataset + hsize_t dim = totalElements; + h5space = H5Screate_simple(1, &dim, 0L); + checkH5Err(h5space); + hid_t h5pcreate = H5Pcreate(H5P_DATASET_CREATE); + checkH5Err(h5pcreate); + hsize_t chunkDim = std::min(chunkSize, dim); + checkH5Err(H5Pset_chunk(h5pcreate, 1, &chunkDim)); + // checkH5Err(H5Pset_szip(h5pcreate, H5_SZIP_NN_OPTION_MASK, 4)); + checkH5Err(H5Pset_deflate(h5pcreate, 5)); + h5var = + H5Dcreate(h5File, "/is_left", H5T_STD_U32LE, h5space, H5P_DEFAULT, h5pcreate, H5P_DEFAULT); + checkH5Err(h5var); + checkH5Err(H5Pclose(h5pcreate)); + } - double outsideborder1 = 0; - double outsideborder2 = 0; + hsize_t dim = totalElements; + hid_t h5memspace = H5Screate_simple(1, &dim, 0L); + checkH5Err(h5memspace); - for (size_t p = 0; p < partitions; p++) { - totalElements += elementSize[p]; + checkH5Err(H5Sselect_all(h5memspace)); - size_t offset[3] = {p, 0, 0}; - size_t size[3] = {1, static_cast(elementSize[p]), 4}; + checkH5Err(H5Dwrite(h5var, H5T_NATIVE_UINT32, h5memspace, h5space, H5P_DEFAULT, isLeft)); - checkNcError(nc_get_vara_int(ncFile, ncElemVertices, offset, size, elementVertices)); - checkNcError(nc_get_vara_int(ncFile, ncElemBoundaries, offset, size, elementBoundaries)); - - size[1] = vertexSize[p]; size[2] = 3; -// logInfo() << "size" << size[0] << size[1] << size[2]; - checkNcError(nc_get_vara_double(ncFile, ncVrtxCoords, offset, size, vertexCoordinates)); - - for (unsigned int i = 0; i < elementSize[p]; i++) { - for (unsigned int j = 0; j < 4; j++) { - if (elementBoundaries[i*4+j] == 3) { - // Is a fault boundary - for (unsigned int k = 0; k < 3; k++) { - Point point; - memcpy(point.coords, &vertexCoordinates[elementVertices[i*4 + FACE2NODES[j][k]]*3], sizeof(double)*3); - faultPoints.insert(point); + checkH5Err(H5Sclose(h5memspace)); + checkH5Err(H5Dclose(h5var)); + checkH5Err(H5Sclose(h5space)); - minfault = std::min(minfault, point.coords[direction]); - maxfault = std::max(maxfault, point.coords[direction]); + delete[] isLeft; - if (point.coords[1-direction] < minfault2) { - minfault2 = point.coords[1-direction]; - outsideborder1 = point.coords[direction]; - } - if (point.coords[1-direction] > maxfault2) { - maxfault2 = point.coords[1-direction]; - outsideborder2 = point.coords[direction]; - } - } - } - } - } - } + checkH5Err(H5Fclose(h5File)); - logInfo() << "Cells <" << minfault2 << "and >" << maxfault2 << "have the be handled manually"; - logInfo() << "A good choice might be" << outsideborder1 << "reps." << outsideborder2; - - bool split[3] = {true, true, true}; - split[direction] = false; - KDTree kdtree(faultPoints, 4, split); - - faultPoints.clear(); // Free memory - - unsigned int* isLeft = new unsigned int[totalElements]; - std::fill(isLeft, isLeft+totalElements, 1); - - unsigned long globElement = 0; - for (size_t p = 0; p < partitions; p++) { - logInfo() << "Processing partition" << utils::nospace << p << "..."; - - size_t offset[3] = {p, 0, 0}; - size_t size[3] = {1, static_cast(elementSize[p]), 4}; - - checkNcError(nc_get_vara_int(ncFile, ncElemVertices, offset, size, elementVertices)); - - size[1] = vertexSize[p]; size[2] = 3; - checkNcError(nc_get_vara_double(ncFile, ncVrtxCoords, offset, size, vertexCoordinates)); - - for (unsigned int i = 0; i < elementSize[p]; i++) { - double avg = 0; - - Support sup; - for (unsigned int j = 0; j < 4; j++) { - sup.limits[0][0] = std::min(sup.limits[0][0], vertexCoordinates[elementVertices[i*4+j]*3]); - sup.limits[0][1] = std::max(sup.limits[0][1], vertexCoordinates[elementVertices[i*4+j]*3]); - sup.limits[1][0] = std::min(sup.limits[1][0], vertexCoordinates[elementVertices[i*4+j]*3+1]); - sup.limits[1][1] = std::max(sup.limits[1][1], vertexCoordinates[elementVertices[i*4+j]*3+1]); - sup.limits[2][0] = std::min(sup.limits[1][0], vertexCoordinates[elementVertices[i*4+j]*3+2]); - sup.limits[2][1] = std::max(sup.limits[1][1], vertexCoordinates[elementVertices[i*4+j]*3+2]); - - avg += vertexCoordinates[elementVertices[i*4+j]*3 + direction]; - } - - avg /= 4; - - if (avg < minfault) { - // Do nothing - } else if (avg > maxfault) { - isLeft[globElement] = 0; - } else { - Action act; - kdtree.search(sup, act); - - bool hasPoint = false; - - for (std::vector::const_iterator it = act.points.begin(); - it != act.points.end(); it++) { - if (it->coords[1-direction] >= sup.limits[1-direction][0] && it->coords[1-direction] <= sup.limits[1-direction][1] - && it->z >= sup.limits[2][0] && it->z <= sup.limits[2][1]) { - hasPoint = true; - - if (it->coords[direction] < avg) { - isLeft[globElement] = 0; - break; - } - } - } - - if (!hasPoint) { - if (avg > minfault) { - isLeft[globElement] = 0; - } - } - } - - globElement++; - } - } - - checkNcError(nc_close(ncFile)); - - delete [] elementVertices; - delete [] elementBoundaries; - delete [] vertexCoordinates; - delete [] elementSize; - delete [] vertexSize; - - hid_t h5File = H5Fopen(output.c_str(), H5F_ACC_RDWR, H5P_DEFAULT); - checkH5Err(h5File); - - hsize_t chunkSize = 256*1024*1024 / sizeof(unsigned int); - - hid_t h5var = H5Dopen(h5File, "/is_left", H5P_DEFAULT); - hid_t h5space; - if (h5var >= 0) { - logWarning() << "Overwriting old dataset from HDF5 file"; - - h5space = H5Dget_space(h5var); - checkH5Err(h5space); - int ndims = H5Sget_simple_extent_ndims(h5space); - checkH5Err(ndims); - if (ndims != 1) - logError() << "Old dataset has wrong dimension"; - - hsize_t extent; - checkH5Err(H5Sget_simple_extent_dims(h5space, &extent, 0L)); - if (extent != totalElements) - logError() << "Old dataset has wrong size"; - } else { - // Create new dataset - hsize_t dim = totalElements; - h5space = H5Screate_simple(1, &dim, 0L); - checkH5Err(h5space); - hid_t h5pcreate = H5Pcreate(H5P_DATASET_CREATE); - checkH5Err(h5pcreate); - hsize_t chunkDim = std::min(chunkSize, dim); - checkH5Err(H5Pset_chunk(h5pcreate, 1, &chunkDim)); -// checkH5Err(H5Pset_szip(h5pcreate, H5_SZIP_NN_OPTION_MASK, 4)); - checkH5Err(H5Pset_deflate(h5pcreate, 5)); - h5var = H5Dcreate(h5File, "/is_left", H5T_STD_U32LE, h5space, - H5P_DEFAULT, h5pcreate, H5P_DEFAULT); - checkH5Err(h5var); - checkH5Err(H5Pclose(h5pcreate)); - } - - hsize_t dim = totalElements; - hid_t h5memspace = H5Screate_simple(1, &dim, 0L); - checkH5Err(h5memspace); - - checkH5Err(H5Sselect_all(h5memspace)); - - checkH5Err(H5Dwrite(h5var, H5T_NATIVE_UINT32, h5memspace, h5space, H5P_DEFAULT, isLeft)); - - checkH5Err(H5Sclose(h5memspace)); - checkH5Err(H5Dclose(h5var)); - checkH5Err(H5Sclose(h5space)); - - delete [] isLeft; - - checkH5Err(H5Fclose(h5File)); - - return 0; + return 0; } diff --git a/trace/src/CMakeLists.txt b/trace/src/CMakeLists.txt index 9c744c3..4d07fee 100644 --- a/trace/src/CMakeLists.txt +++ b/trace/src/CMakeLists.txt @@ -12,4 +12,3 @@ target_link_libraries(trace PRIVATE ${NETCDF_LIBRARIES} ${HDF5_HL_LIBRARIES} ${HDF5_C_LIBRARIES} MPI::MPI_CXX) - diff --git a/trace/src/image.cpp b/trace/src/image.cpp index 5af15fe..a8c8b03 100644 --- a/trace/src/image.cpp +++ b/trace/src/image.cpp @@ -1,4 +1,5 @@ #include "image.h" + #include #define STB_IMAGE_WRITE_IMPLEMENTATION @@ -7,75 +8,75 @@ #include #include -Image::Image(int nranks, int nsubregions, timespec const& t_min, timespec const& t_max, - long resolution) +Image::Image( + int nranks, int nsubregions, const timespec& t_min, const timespec& t_max, long resolution) : nsubregions_(nsubregions), t_min_(t_min), t_max_(t_max), resolution_(resolution) { - auto time_diff = difftime(t_min, t_max); - std::size_t time_slots = 1 + (time_diff - 1) / resolution; + auto time_diff = difftime(t_min, t_max); + std::size_t time_slots = 1 + (time_diff - 1) / resolution; - std::size_t root = sqrt(time_slots / static_cast(nranks + Offset)); - width_ = 1 + (time_slots - 1) / (root); - height_ = nranks * root + Offset * (root - 1); - stride_ = width_ * (nranks + Offset); + std::size_t root = sqrt(time_slots / static_cast(nranks + Offset)); + width_ = 1 + (time_slots - 1) / (root); + height_ = nranks * root + Offset * (root - 1); + stride_ = width_ * (nranks + Offset); - auto npixels = width_ * height_ * Comp; - px_ = new unsigned char[npixels]; - memset(px_, 0, npixels * sizeof(unsigned char)); + auto npixels = width_ * height_ * Comp; + px_ = new unsigned char[npixels]; + memset(px_, 0, npixels * sizeof(unsigned char)); - for (int r = 0; r < root - 1; ++r) { - for (int o = 0; o < nranks; ++o) { - for (int c = 0; c < Comp * width_; ++c) { - px_[Comp * (r * stride_ + o * width_) + c] = 255; - } - } + for (int r = 0; r < root - 1; ++r) { + for (int o = 0; o < nranks; ++o) { + for (int c = 0; c < Comp * width_; ++c) { + px_[Comp * (r * stride_ + o * width_) + c] = 255; + } } + } } Image::~Image() { delete[] px_; } -std::tuple Image::time_slots(timespec const& begin, - timespec const& end) { - auto d1 = difftime(t_min_, begin); - auto d2 = difftime(t_min_, end); - long long t1 = d1 / resolution_; - long long t2 = 1 + (d2 - 1) / resolution_; - double w1 = static_cast(d1) / resolution_ - t1; - double w2 = t2 - static_cast(d2) / resolution_; - return {t1, w1, t2, w2}; +std::tuple Image::time_slots(const timespec& begin, + const timespec& end) { + auto d1 = difftime(t_min_, begin); + auto d2 = difftime(t_min_, end); + long long t1 = d1 / resolution_; + long long t2 = 1 + (d2 - 1) / resolution_; + double w1 = static_cast(d1) / resolution_ - t1; + double w2 = t2 - static_cast(d2) / resolution_; + return {t1, w1, t2, w2}; } unsigned char* Image::pixel(int rank, long long time_slot) { - auto const block_row = time_slot / width_; - auto const block_begin = time_slot % width_; - return px_ + Comp * (block_row * stride_ + rank * width_ + block_begin); + const auto block_row = time_slot / width_; + const auto block_begin = time_slot % width_; + return px_ + Comp * (block_row * stride_ + rank * width_ + block_begin); } -void Image::add_pixels(int rank, timespec const& begin, timespec const& end, int subRegion) { - if (region_ >= sizeof(Colors) / sizeof(int)) { - throw std::runtime_error("OMG! Not enough colours."); - } - auto [t1, w1, t2, w2] = time_slots(begin, end); - for (auto t = t1; t <= t2; ++t) { - double w = t == t1 ? w1 : (t == t2 ? w2 : 1.0); - unsigned char* px = pixel(rank, t); - for (int c = 0; c < Comp; ++c) { - int channel = 8 * (Comp - c - 1); - unsigned char col = (Colors[region_] & (0xff << channel)) >> channel; - col *= (nsubregions_ - subRegion) / static_cast(nsubregions_); - px[c] = w * col + (1.0 - w) * px[c]; - } +void Image::add_pixels(int rank, const timespec& begin, const timespec& end, int subRegion) { + if (region_ >= sizeof(Colors) / sizeof(int)) { + throw std::runtime_error("OMG! Not enough colours."); + } + auto [t1, w1, t2, w2] = time_slots(begin, end); + for (auto t = t1; t <= t2; ++t) { + double w = t == t1 ? w1 : (t == t2 ? w2 : 1.0); + unsigned char* px = pixel(rank, t); + for (int c = 0; c < Comp; ++c) { + int channel = 8 * (Comp - c - 1); + unsigned char col = (Colors[region_] & (0xff << channel)) >> channel; + col *= (nsubregions_ - subRegion) / static_cast(nsubregions_); + px[c] = w * col + (1.0 - w) * px[c]; } + } } -void Image::add(std::vector const& offset, std::vector const& sample) { - for (std::size_t i = 0; i < offset.size() - 1; ++i) { - for (int j = offset[i]; j < offset[i + 1]; ++j) { - add_pixels(i, sample[j].begin, sample[j].end, sample[j].subRegion); - } +void Image::add(const std::vector& offset, const std::vector& sample) { + for (std::size_t i = 0; i < offset.size() - 1; ++i) { + for (int j = offset[i]; j < offset[i + 1]; ++j) { + add_pixels(i, sample[j].begin, sample[j].end, sample[j].subRegion); } - ++region_; + } + ++region_; } -void Image::write(char const* file_name) { - stbi_write_png(file_name, width_, height_, Comp, px_, width_ * Comp); +void Image::write(const char* file_name) { + stbi_write_png(file_name, width_, height_, Comp, px_, width_ * Comp); } diff --git a/trace/src/image.h b/trace/src/image.h index a07cb38..a856c1b 100644 --- a/trace/src/image.h +++ b/trace/src/image.h @@ -7,35 +7,38 @@ #include #include -inline long long difftime(timespec const& start, timespec const& end) { - return 1000000000L * (end.tv_sec - start.tv_sec) + end.tv_nsec - start.tv_nsec; +inline long long difftime(const timespec& start, const timespec& end) { + return 1000000000L * (end.tv_sec - start.tv_sec) + end.tv_nsec - start.tv_nsec; } class Image { -public: - static constexpr int Comp = 3; - static constexpr int Offset = 5; - static constexpr unsigned int Colors[] = {0xff0000, 0x00ff00, 0x0000ff}; - - Image(int nranks, int nsubregions, timespec const& t_min, timespec const& t_max, - long resolution = 1000000000L); - ~Image(); - void add(std::vector const& offset, std::vector const& sample); - void write(char const* file_name); - -private: - std::tuple time_slots(timespec const& begin, - timespec const& end); - unsigned char* pixel(int rank, long long time_slot); - void add_pixels(int rank, timespec const& begin, timespec const& end, int subRegion); - - int nsubregions_; - timespec t_min_, t_max_; - long long resolution_; - int region_ = 0; - - int width_, height_, stride_; - unsigned char* px_; + public: + static constexpr int Comp = 3; + static constexpr int Offset = 5; + static constexpr unsigned int Colors[] = {0xff0000, 0x00ff00, 0x0000ff}; + + Image(int nranks, + int nsubregions, + const timespec& t_min, + const timespec& t_max, + long resolution = 1000000000L); + ~Image(); + void add(const std::vector& offset, const std::vector& sample); + void write(const char* file_name); + + private: + std::tuple time_slots(const timespec& begin, + const timespec& end); + unsigned char* pixel(int rank, long long time_slot); + void add_pixels(int rank, const timespec& begin, const timespec& end, int subRegion); + + int nsubregions_; + timespec t_min_, t_max_; + long long resolution_; + int region_ = 0; + + int width_, height_, stride_; + unsigned char* px_; }; #endif // IMAGE_20201125_H diff --git a/trace/src/input.cpp b/trace/src/input.cpp index df19617..c7abb2c 100644 --- a/trace/src/input.cpp +++ b/trace/src/input.cpp @@ -5,52 +5,52 @@ #include void check_err(const int stat, const int line, const char* file) { - if (stat != NC_NOERR) { - (void)fprintf(stderr, "line %d of %s: %s\n", line, file, nc_strerror(stat)); - fflush(stderr); - exit(1); - } + if (stat != NC_NOERR) { + (void)fprintf(stderr, "line %d of %s: %s\n", line, file, nc_strerror(stat)); + fflush(stderr); + exit(1); + } } -std::pair, std::vector> read(char const* file_name) { - int ncid; +std::pair, std::vector> read(const char* file_name) { + int ncid; - /* dimension ids */ - int rank_dim; - int sample_dim; + /* dimension ids */ + int rank_dim; + int sample_dim; - /* dimension lengths */ - size_t rank_len; - size_t sample_len; + /* dimension lengths */ + size_t rank_len; + size_t sample_len; - /* variable ids */ - int offset_id; - int sample_id; + /* variable ids */ + int offset_id; + int sample_id; - size_t start[] = {0}; - size_t count[1]; + size_t start[] = {0}; + size_t count[1]; - check_err(nc_open(file_name, 0, &ncid), __LINE__, __FILE__); + check_err(nc_open(file_name, 0, &ncid), __LINE__, __FILE__); - /* get dimensions */ - check_err(nc_inq_dimid(ncid, "rank", &rank_dim), __LINE__, __FILE__); - check_err(nc_inq_dimlen(ncid, rank_dim, &rank_len), __LINE__, __FILE__); - check_err(nc_inq_dimid(ncid, "sample", &sample_dim), __LINE__, __FILE__); - check_err(nc_inq_dimlen(ncid, sample_dim, &sample_len), __LINE__, __FILE__); + /* get dimensions */ + check_err(nc_inq_dimid(ncid, "rank", &rank_dim), __LINE__, __FILE__); + check_err(nc_inq_dimlen(ncid, rank_dim, &rank_len), __LINE__, __FILE__); + check_err(nc_inq_dimid(ncid, "sample", &sample_dim), __LINE__, __FILE__); + check_err(nc_inq_dimlen(ncid, sample_dim, &sample_len), __LINE__, __FILE__); - auto offset = std::vector(rank_len); - auto sample = std::vector(sample_len); + auto offset = std::vector(rank_len); + auto sample = std::vector(sample_len); - /* get variable ids */ - check_err(nc_inq_varid(ncid, "offset", &offset_id), __LINE__, __FILE__); - check_err(nc_inq_varid(ncid, "sample", &sample_id), __LINE__, __FILE__); + /* get variable ids */ + check_err(nc_inq_varid(ncid, "offset", &offset_id), __LINE__, __FILE__); + check_err(nc_inq_varid(ncid, "sample", &sample_id), __LINE__, __FILE__); - count[0] = rank_len; - check_err(nc_get_vara(ncid, offset_id, start, count, offset.data()), __LINE__, __FILE__); - count[0] = sample_len; - check_err(nc_get_vara(ncid, sample_id, start, count, sample.data()), __LINE__, __FILE__); + count[0] = rank_len; + check_err(nc_get_vara(ncid, offset_id, start, count, offset.data()), __LINE__, __FILE__); + count[0] = sample_len; + check_err(nc_get_vara(ncid, sample_id, start, count, sample.data()), __LINE__, __FILE__); - check_err(nc_close(ncid), __LINE__, __FILE__); + check_err(nc_close(ncid), __LINE__, __FILE__); - return {offset, sample}; + return {offset, sample}; } diff --git a/trace/src/input.h b/trace/src/input.h index a496be2..6c39d2e 100644 --- a/trace/src/input.h +++ b/trace/src/input.h @@ -6,6 +6,6 @@ #include #include -std::pair, std::vector> read(char const* file_name); +std::pair, std::vector> read(const char* file_name); #endif // TRACE_20201125_H diff --git a/trace/src/sample.h b/trace/src/sample.h index cec06c4..90c2780 100644 --- a/trace/src/sample.h +++ b/trace/src/sample.h @@ -4,10 +4,10 @@ #include struct Sample { - timespec begin; - timespec end; - unsigned int loopLength; - unsigned int subRegion; + timespec begin; + timespec end; + unsigned int loopLength; + unsigned int subRegion; }; #endif // SAMPLE_20201125_H diff --git a/trace/src/trace.cpp b/trace/src/trace.cpp index b9f3673..a2f6ba0 100644 --- a/trace/src/trace.cpp +++ b/trace/src/trace.cpp @@ -8,54 +8,53 @@ #include #include -bool operator<(timespec const& a, timespec const& b) { - if (a.tv_sec == b.tv_sec) { - return a.tv_nsec < b.tv_nsec; - } - return a.tv_sec < b.tv_sec; +bool operator<(const timespec& a, const timespec& b) { + if (a.tv_sec == b.tv_sec) { + return a.tv_nsec < b.tv_nsec; + } + return a.tv_sec < b.tv_sec; } int main(int argc, char** argv) { - if (argc < 3) { - std::cerr << "Usage: ./trace [resolution in nsec]" << std::endl; - return -1; - } - const auto prefix = std::string(argv[1]); - char const* output = argv[2]; - long resolution = argc >= 4 ? atoi(argv[3]) : 1000000000L; - - const std::vector names{"computeLocalIntegration.nc", - "computeNeighboringIntegration.nc", - "computeDynamicRupture.nc"}; - std::vector, std::vector>> stats; - for (auto const& name : names) { - const auto f = prefix + name; - stats.emplace_back(read(f.c_str())); - } - - auto cmp_begin = [](auto const& a, auto const& b) { return a.begin < b.begin; }; - auto cmp_end = [](auto const& a, auto const& b) { return a.end < b.end; }; - auto cmp_subRegion = [](auto const& a, auto const& b) { return a.subRegion < b.subRegion; }; - - timespec t_min{std::numeric_limits::max(), std::numeric_limits::max()}; - timespec t_max{0, 0}; - unsigned subRegion_max = 0; - - for (const auto& st : stats) { - auto my_t_min = std::min_element(st.second.begin(), st.second.end(), cmp_begin)->begin; - t_min = std::min(t_min, my_t_min); - auto my_t_max = std::max_element(st.second.begin(), st.second.end(), cmp_end)->end; - t_max = std::max(t_max, my_t_max); - auto my_subRegion_max = - std::max_element(st.second.begin(), st.second.end(), cmp_subRegion)->subRegion; - subRegion_max = std::max(subRegion_max, my_subRegion_max); - } - - auto im = Image(stats[0].first.size() - 1, subRegion_max + 1, t_min, t_max, resolution); - for (const auto& st : stats) { - im.add(st.first, st.second); - } - im.write(output); - - return 0; + if (argc < 3) { + std::cerr << "Usage: ./trace [resolution in nsec]" << std::endl; + return -1; + } + const auto prefix = std::string(argv[1]); + const char* output = argv[2]; + long resolution = argc >= 4 ? atoi(argv[3]) : 1000000000L; + + const std::vector names{ + "computeLocalIntegration.nc", "computeNeighboringIntegration.nc", "computeDynamicRupture.nc"}; + std::vector, std::vector>> stats; + for (const auto& name : names) { + const auto f = prefix + name; + stats.emplace_back(read(f.c_str())); + } + + auto cmp_begin = [](const auto& a, const auto& b) { return a.begin < b.begin; }; + auto cmp_end = [](const auto& a, const auto& b) { return a.end < b.end; }; + auto cmp_subRegion = [](const auto& a, const auto& b) { return a.subRegion < b.subRegion; }; + + timespec t_min{std::numeric_limits::max(), std::numeric_limits::max()}; + timespec t_max{0, 0}; + unsigned subRegion_max = 0; + + for (const auto& st : stats) { + auto my_t_min = std::min_element(st.second.begin(), st.second.end(), cmp_begin)->begin; + t_min = std::min(t_min, my_t_min); + auto my_t_max = std::max_element(st.second.begin(), st.second.end(), cmp_end)->end; + t_max = std::max(t_max, my_t_max); + auto my_subRegion_max = + std::max_element(st.second.begin(), st.second.end(), cmp_subRegion)->subRegion; + subRegion_max = std::max(subRegion_max, my_subRegion_max); + } + + auto im = Image(stats[0].first.size() - 1, subRegion_max + 1, t_min, t_max, resolution); + for (const auto& st : stats) { + im.add(st.first, st.second); + } + im.write(output); + + return 0; }