Mercurial > octave-libtiff
view libinterp/corefcn/__contourc__.cc @ 31215:0fb5ca242d66
__tiff__.cc: added documentation to all internal functions
* libtinterp/corefcn/__tiff__.cc: Added documentation to all internal
functions to enable using help and print_usage.
author | magedrifaat <magedrifaat@gmail.com> |
---|---|
date | Wed, 07 Sep 2022 18:00:04 +0200 |
parents | 796f54d4ddbf |
children |
line wrap: on
line source
//////////////////////////////////////////////////////////////////////// // // Contour lines for function evaluated on a grid. // // Copyright (C) 2001-2022 The Octave Project Developers // // See the file COPYRIGHT.md in the top-level directory of this // distribution or <https://octave.org/copyright/>. // // Adapted to an oct file from the stand alone contourl by Victro Munoz // Copyright (C) 2004 Victor Munoz // // Based on contour plot routine (plcont.c) in PLPlot package // http://plplot.org/ // // Copyright (C) 1995, 2000, 2001 Maurice LeBrun // Copyright (C) 2000, 2002 Joao Cardoso // Copyright (C) 2000, 2001, 2002, 2004 Alan W. Irwin // Copyright (C) 2004 Andrew Ross // // // This file is part of Octave. // // Octave is free software: you can redistribute it and/or modify it // under the terms of the GNU General Public License as published by // the Free Software Foundation, either version 3 of the License, or // (at your option) any later version. // // Octave is distributed in the hope that it will be useful, but // WITHOUT ANY WARRANTY; without even the implied warranty of // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the // GNU General Public License for more details. // // You should have received a copy of the GNU General Public License // along with Octave; see the file COPYING. If not, see // <https://www.gnu.org/licenses/>. // //////////////////////////////////////////////////////////////////////// #if defined (HAVE_CONFIG_H) # include "config.h" #endif #include <limits> #include "defun.h" #include "ov.h" OCTAVE_NAMESPACE_BEGIN // FIXME: this looks like trouble... static Matrix this_contour; static Matrix contourc; static int elem; // This is the quanta in which we increase this_contour. #define CONTOUR_QUANT 50 // Add a coordinate point (x,y) to this_contour. static void add_point (double x, double y) { if (elem % CONTOUR_QUANT == 0) this_contour = this_contour.append (Matrix (2, CONTOUR_QUANT, 0)); this_contour(0, elem) = x; this_contour(1, elem) = y; elem++; } // Add contents of current contour to contourc. // this_contour.cols () - 1; static void end_contour (void) { if (elem > 2) { this_contour(1, 0) = elem - 1; contourc = contourc.append (this_contour.extract_n (0, 0, 2, elem)); } this_contour = Matrix (); elem = 0; } // Start a new contour, and add contents of current one to contourc. static void start_contour (double lvl, double x, double y) { end_contour (); this_contour.resize (2, 0); add_point (lvl, 0); add_point (x, y); } static void drawcn (const RowVector& X, const RowVector& Y, const Matrix& Z, double lvl, int r, int c, double ct_x, double ct_y, unsigned int start_edge, bool first, charMatrix& mark) { double px[4], py[4], pz[4], tmp; unsigned int stop_edge, pt[2]; // Continue while next facet is not done yet. while (r >= 0 && c >= 0 && r < mark.rows () && c < mark.cols () && mark(r, c) > 0) { //get x, y, and z - lvl for current facet px[0] = px[3] = X(c); px[1] = px[2] = X(c+1); py[0] = py[1] = Y(r); py[2] = py[3] = Y(r+1); pz[3] = Z(r+1, c) - lvl; pz[2] = Z(r+1, c + 1) - lvl; pz[1] = Z(r, c+1) - lvl; pz[0] = Z(r, c) - lvl; // Facet edge and point naming assignment. // // 0-----1 .-0-. // | | | | // | | 3 1 // | | | | // 3-----2 .-2-. // Get mark value of current facet. char id = static_cast<char> (mark(r, c)); // Check startedge s. if (start_edge == 255) { // Find start edge. for (unsigned int k = 0; k < 4; k++) if (static_cast<char> (1 << k) & id) start_edge = k; } if (start_edge == 255) break; // Decrease mark value of current facet for start edge. mark(r, c) -= static_cast<char> (1 << start_edge); // Next point (clockwise). pt[0] = start_edge; pt[1] = (pt[0] + 1) % 4; // Calculate contour segment start if first of contour. if (first) { tmp = fabs (pz[pt[1]]) / fabs (pz[pt[0]]); if (math::isnan (tmp)) ct_x = ct_y = 0.5; else { ct_x = px[pt[0]] + (px[pt[1]] - px[pt[0]])/(1 + tmp); ct_y = py[pt[0]] + (py[pt[1]] - py[pt[0]])/(1 + tmp); } start_contour (lvl, ct_x, ct_y); first = false; } // Find stop edge. // FIXME: perhaps this should use a while loop? for (unsigned int k = 1; k <= 4; k++) { if (start_edge == 0 || start_edge == 2) stop_edge = (start_edge + k) % 4; else stop_edge = (start_edge - k) % 4; if (static_cast<char> (1 << stop_edge) & id) break; } pt[0] = stop_edge; pt[1] = (pt[0] + 1) % 4; tmp = fabs (pz[pt[1]]) / fabs (pz[pt[0]]); if (math::isnan (tmp)) ct_x = ct_y = 0.5; else { ct_x = px[pt[0]] + (px[pt[1]] - px[pt[0]])/(1 + tmp); ct_y = py[pt[0]] + (py[pt[1]] - py[pt[0]])/(1 + tmp); } // Add point to contour. add_point (ct_x, ct_y); // Decrease id value of current facet for start edge. mark(r, c) -= static_cast<char> (1 << stop_edge); // Find next facet. if (stop_edge == 0) r--; else if (stop_edge == 1) c++; else if (stop_edge == 2) r++; else if (stop_edge == 3) c--; // Go to next facet. start_edge = (stop_edge + 2) % 4; } } static void mark_facets (const Matrix& Z, charMatrix& mark, double lvl) { unsigned int nr = mark.rows (); unsigned int nc = mark.cols (); double f[4]; for (unsigned int c = 0; c < nc; c++) for (unsigned int r = 0; r < nr; r++) { f[0] = Z(r, c) - lvl; f[1] = Z(r, c+1) - lvl; f[3] = Z(r+1, c) - lvl; f[2] = Z(r+1, c+1) - lvl; for (unsigned int i = 0; i < 4; i++) if (fabs(f[i]) < std::numeric_limits<double>::epsilon ()) f[i] = std::numeric_limits<double>::epsilon (); if (f[1] * f[2] < 0) mark(r, c) += 2; if (f[0] * f[3] < 0) mark(r, c) += 8; } for (unsigned int r = 0; r < nr; r++) for (unsigned int c = 0; c < nc; c++) { f[0] = Z(r, c) - lvl; f[1] = Z(r, c+1) - lvl; f[3] = Z(r+1, c) - lvl; f[2] = Z(r+1, c+1) - lvl; for (unsigned int i = 0; i < 4; i++) if (fabs(f[i]) < std::numeric_limits<double>::epsilon ()) f[i] = std::numeric_limits<double>::epsilon (); if (f[0] * f[1] < 0) mark(r, c) += 1; if (f[2] * f[3] < 0) mark(r, c) += 4; } } static void cntr (const RowVector& X, const RowVector& Y, const Matrix& Z, double lvl) { unsigned int nr = Z.rows (); unsigned int nc = Z.cols (); charMatrix mark (nr - 1, nc - 1, 0); mark_facets (Z, mark, lvl); // Find contours that start at a domain edge. for (unsigned int c = 0; c < nc - 1; c++) { // Top. if (mark(0, c) & 1) drawcn (X, Y, Z, lvl, 0, c, 0.0, 0.0, 0, true, mark); // Bottom. if (mark(nr - 2, c) & 4) drawcn (X, Y, Z, lvl, nr - 2, c, 0.0, 0.0, 2, true, mark); } for (unsigned int r = 0; r < nr - 1; r++) { // Left. if (mark(r, 0) & 8) drawcn (X, Y, Z, lvl, r, 0, 0.0, 0.0, 3, true, mark); // Right. if (mark(r, nc - 2) & 2) drawcn (X, Y, Z, lvl, r, nc - 2, 0.0, 0.0, 1, true, mark); } for (unsigned int r = 0; r < nr - 1; r++) for (unsigned int c = 0; c < nc - 1; c++) if (mark (r, c) > 0) drawcn (X, Y, Z, lvl, r, c, 0.0, 0.0, 255, true, mark); } DEFUN (__contourc__, args, , doc: /* -*- texinfo -*- @deftypefn {} {@var{c} =} __contourc__ (@var{x}, @var{y}, @var{z}, @var{levels}) Calculate Z-level contours (isolines). @end deftypefn */) { RowVector X = args(0).row_vector_value (); RowVector Y = args(1).row_vector_value (); Matrix Z = args(2).matrix_value (); RowVector L = args(3).row_vector_value (); contourc.resize (2, 0); for (int i = 0; i < L.numel (); i++) cntr (X, Y, Z, L(i)); end_contour (); return octave_value (contourc); } /* ## No test needed for internal helper function. %!assert (1) */ OCTAVE_NAMESPACE_END