Branch data Line data Source code
1 : : // *****************************************************************************
2 : : /*!
3 : : \file src/IO/PDFWriter.cpp
4 : : \copyright 2012-2015 J. Bakosi,
5 : : 2016-2018 Los Alamos National Security, LLC.,
6 : : 2019-2021 Triad National Security, LLC.,
7 : : 2022-2024 J. Bakosi
8 : : All rights reserved. See the LICENSE file for details.
9 : : \brief Univariate PDF writer
10 : : \brief PDF writer class definition
11 : : \details This file defines a PDF writer class that facilitates outputing
12 : : probability density functions (PDFs) into files in various formats using
13 : : various configurations.
14 : : */
15 : : // *****************************************************************************
16 : :
17 : : #include <iomanip>
18 : :
19 : : #include "NoWarning/exodusII.hpp"
20 : :
21 : : #include "PDFWriter.hpp"
22 : : #include "Exception.hpp"
23 : :
24 : : using tk::PDFWriter;
25 : :
26 : 726 : PDFWriter::PDFWriter( const std::string& filename,
27 : : const std::string& format,
28 : 726 : std::streamsize precision ) :
29 : 726 : Writer( filename )
30 : : // *****************************************************************************
31 : : // Constructor
32 : : //! \param[in] filename Output filename to which output the PDF
33 : : //! \param[in] format Configure floating-point output format for ASCII output
34 : : //! \param[in] precision Configure precision for floating-point ASCII output
35 : : // *****************************************************************************
36 : : {
37 : : // Set floating-point format for output file stream
38 [ - + ]: 726 : if (format == "default")
39 : : {} //m_outFile << std::defaultfloat; GCC does not yet support this
40 [ - - ]: 0 : else if (format == "fixed")
41 [ - - ]: 0 : m_outFile << std::fixed;
42 [ - - ]: 0 : else if (format == "scientific")
43 [ - - ]: 0 : m_outFile << std::scientific;
44 [ - - ][ - - ]: 0 : else Throw( "Text floating-point format not recognized." );
[ - - ]
45 : :
46 : : // Set numeric precision for output file stream if the input makes sense
47 [ + - ][ + - ]: 726 : if (precision > 0 && precision < std::numeric_limits< tk::real >::digits10+2)
48 : 726 : m_outFile << std::setprecision( static_cast<int>(precision) );
49 : 726 : }
50 : :
51 : : void
52 : 726 : PDFWriter::writeTxt( const UniPDF& pdf, const tk::ctr::PDFInfo& info ) const
53 : : // *****************************************************************************
54 : : // Write out standardized univariate PDF to file
55 : : //! \param[in] pdf Univariate PDF
56 : : //! \param[in] info PDF metadata
57 : : // *****************************************************************************
58 : : {
59 : 726 : const auto& name = info.name;
60 : 726 : const auto& uext = info.exts;
61 : 726 : const auto& vars = info.vars;
62 : 726 : const auto& it = info.it;
63 : 726 : const auto& time = info.time;
64 : :
65 [ + - ]: 726 : assertSampleSpaceDimensions< 1 >( vars );
66 [ + - ]: 726 : assertSampleSpaceExtents< 1 >( uext );
67 : :
68 : : // Query and optionally override number of bins and minimum of sample space if
69 : : // user-specified extents were given and copy probabilities from pdf to an
70 : : // array for output
71 : : std::size_t nbi;
72 : : tk::real min, max;
73 : 726 : std::vector< tk::real > outpdf;
74 : : tk::real binsize;
75 : : std::array< long, 2*UniPDF::dim > ext;
76 [ + - ]: 726 : extents( pdf, uext, nbi, min, max, binsize, ext, outpdf );
77 : :
78 : : // Output header
79 : : m_outFile << "# vim: filetype=sh:\n#\n"
80 : 726 : << "# Univariate PDF: " << name << '(' << vars[0] << ')' << '\n'
81 : : << "# -----------------------------------------------\n"
82 [ + - ][ + - ]: 726 : << "# Numeric precision: " << m_outFile.precision() << '\n'
[ + - ][ + - ]
[ + - ][ + - ]
[ + - ][ + - ]
[ + - ][ + - ]
83 [ + - ][ + - ]: 726 : << "# Bin size: " << binsize << '\n'
84 [ + - ][ + - ]: 726 : << "# Number of bins estimated: " << ext[1] - ext[0] + 1
[ + - ][ + - ]
85 : : << '\n'
86 [ + - ][ + - ]: 726 : << "# Number of bins output: " << nbi << '\n'
87 [ + - ][ + - ]: 726 : << "# Sample space extent: [" << min << " : " << max << "]\n"
[ + - ][ + - ]
[ + - ]
88 [ + - ][ + - ]: 726 : << "# Integral: " << pdf.integral() << "\n"
[ + - ][ + - ]
[ + - ]
89 [ + - ][ + - ]: 726 : << "# Iteration: " << it << "\n"
[ + - ]
90 [ + - ][ + - ]: 726 : << "# Physical time: " << time << "\n#\n"
[ + - ]
91 : : << "# Example step-by-step visualization with gnuplot\n"
92 : : << "# -----------------------------------------------\n"
93 : : << "# gnuplot> set grid\n"
94 : : << "# gnuplot> unset key\n"
95 : 726 : << "# gnuplot> set xlabel \"" << vars[0] << "\"\n"
96 : 726 : << "# gnuplot> set ylabel \"" << name << "(" << vars[0] << ")\"\n"
97 [ + - ][ + - ]: 726 : << "# gnuplot> plot ";
[ + - ][ + - ]
[ + - ][ + - ]
[ + - ][ + - ]
[ + - ][ + - ]
[ + - ][ + - ]
[ + - ][ + - ]
98 [ - + ][ - - ]: 726 : if (!uext.empty()) m_outFile << "[" << uext[0] << ':' << uext[1] << "] ";
[ - - ][ - - ]
[ - - ][ - - ]
99 : 726 : m_outFile << "\"" << m_filename << "\" with points\n#\n"
100 : : << "# Gnuplot one-liner for quick copy-paste\n"
101 : : << "# -----------------------------------------------\n"
102 : 726 : << "# set grid; unset key; set xlabel \"" << vars[0]
103 : 726 : << "\"; set ylabel \"" << name << "(" << vars[0]
104 [ + - ][ + - ]: 726 : << ")\"; plot";
[ + - ][ + - ]
[ + - ][ + - ]
[ + - ][ + - ]
[ + - ][ + - ]
[ + - ][ + - ]
105 [ - + ][ - - ]: 726 : if (!uext.empty()) m_outFile << " [" << uext[0] << ':' << uext[1] << "]";
[ - - ][ - - ]
[ - - ][ - - ]
106 : 726 : m_outFile << " \"" << m_filename << "\" w p\n#\n"
107 : 1452 : << "# Data columns: " << vars[0] << ", " << name << "(" << vars[0]
108 [ + - ][ + - ]: 726 : << ")\n# -----------------------------------------------\n";
[ + - ][ + - ]
[ + - ][ + - ]
[ + - ][ + - ]
[ + - ][ + - ]
109 : :
110 : : // If no user-specified sample space extents, output pdf map directly
111 [ + - ]: 726 : if (uext.empty()) {
112 [ + + ]: 163071 : for (const auto& p : pdf.map())
113 [ + - ][ + - ]: 162345 : m_outFile << binsize * static_cast<tk::real>(p.first) << '\t'
114 : 162345 : << static_cast<tk::real>(p.second) / binsize /
115 [ + - ]: 162345 : static_cast<tk::real>(pdf.nsample())
116 [ + - ]: 162345 : << std::endl;
117 : : } else { // If user-specified sample space extents, output outpdf array
118 : 0 : std::size_t bin = 0;
119 [ - - ]: 0 : for (const auto& p : outpdf)
120 [ - - ][ - - ]: 0 : m_outFile << binsize * static_cast<tk::real>(bin++) + uext[0] << '\t'
121 [ - - ][ - - ]: 0 : << p << std::endl;
122 : : }
123 : 726 : }
124 : : void
125 : 726 : PDFWriter::extents( const UniPDF& pdf,
126 : : const std::vector< tk::real >& uext,
127 : : std::size_t& nbi,
128 : : tk::real& min,
129 : : tk::real& max,
130 : : tk::real& binsize,
131 : : std::array< long, 2*UniPDF::dim >& ext,
132 : : std::vector< tk::real >& outpdf ) const
133 : : // *****************************************************************************
134 : : // Query extents and other metadata of univariate PDF sample space
135 : : //! \details Query and optionally override number of bins and minimum of sample
136 : : //! space if user-specified extents were given and copy probabilities from
137 : : //! pdf to an array for output for plotting univariate PDF.
138 : : //! \param[in] pdf Univariate PDF object
139 : : //! \param[in] uext User-specified extents of sample space
140 : : //! \param[inout] nbi Number of bins
141 : : //! \param[inout] min Minimum value of sample space
142 : : //! \param[inout] max Maximum value of sample space
143 : : //! \param[inout] binsize Bin size
144 : : //! \param[inout] ext Extents of sample space
145 : : //! \param[inout] outpdf PDF ready to be written out to file
146 : : // *****************************************************************************
147 : : {
148 : 726 : assertSampleSpaceExtents< 1 >( uext );
149 : :
150 : : // Query bin size and extents of sample space from PDF
151 : 726 : binsize = pdf.binsize();
152 : 726 : ext = pdf.extents();
153 : :
154 : : // Compute number of bins of sample space (min bins: 1)
155 [ - + ][ - - ]: 726 : Assert( ext[1] >= ext[0], "Wrong extents in PDFWriter::extents" );
[ - - ][ - - ]
156 : 726 : nbi = static_cast< std::size_t >( ext[1] - ext[0] + 1 );
157 : :
158 : : // Compute minimum and maximum of sample space
159 : 726 : min = binsize * static_cast< tk::real >( ext[0] );
160 : 726 : max = binsize * static_cast< tk::real >( ext[1] );
161 : :
162 : : // Override number of bins and minimum if user-specified extents were given,
163 : : // and copy probabilities from pdf to an array for output
164 [ - + ]: 726 : if (!uext.empty()) {
165 : : // Override number of bins by that based on user-specified extents
166 [ - - ][ - - ]: 0 : Assert( uext[1] >= uext[0],
[ - - ][ - - ]
167 : : "Wrong user-defined extents in PDFWriter::extents" );
168 : 0 : nbi = static_cast< std::size_t >(
169 : 0 : std::lround( (uext[1] - uext[0]) / binsize ) );
170 : : // Override extents
171 : 0 : min = uext[0];
172 : 0 : max = uext[1];
173 : :
174 : : // Size output pdf to user-requested dimensions to overridden nbi and
175 : : // initialize output probabilities to zero
176 [ - - ]: 0 : outpdf = std::vector< tk::real >( nbi, 0.0 );
177 : :
178 : : // Fill requested region of pdf to be output from computed pdf
179 [ - - ]: 0 : for (const auto& p : pdf.map()) {
180 : : // Compute (i.e., shift) bin indices relative to user-requested extents
181 : 0 : const auto bin = p.first - std::lround( uext[0] / binsize );
182 : : // Only copy probability value if shifted bin indices fall within
183 : : // user-requested extents (lower inclusive, upper exclusive)
184 [ - - ][ - - ]: 0 : if (bin >= 0 && bin < std::lround( (uext[1] - uext[0]) / binsize )) {
[ - - ]
185 [ - - ][ - - ]: 0 : Assert( static_cast<std::size_t>(bin) < nbi,
[ - - ][ - - ]
186 : : "Bin overflow in user-specified-extent-based bin "
187 : : "calculation of univariate PDF extents." );
188 : : // Copy normalized probability to output pdf
189 : 0 : outpdf[ static_cast<std::size_t>(bin) ] =
190 : 0 : p.second / binsize / static_cast<tk::real>(pdf.nsample());
191 : : }
192 : : }
193 : : }
194 : 726 : }
|
