Panda3D

pfmFile.cxx

00001 // Filename: pfmFile.cxx
00002 // Created by:  drose (23Dec10)
00003 //
00004 ////////////////////////////////////////////////////////////////////
00005 //
00006 // PANDA 3D SOFTWARE
00007 // Copyright (c) Carnegie Mellon University.  All rights reserved.
00008 //
00009 // All use of this software is subject to the terms of the revised BSD
00010 // license.  You should have received a copy of this license along
00011 // with this source code in a file named "LICENSE."
00012 //
00013 ////////////////////////////////////////////////////////////////////
00014 
00015 #include "config_grutil.h"
00016 #include "pfmFile.h"
00017 #include "virtualFileSystem.h"
00018 #include "pandaFileStream.h"
00019 #include "littleEndian.h"
00020 #include "bigEndian.h"
00021 #include "cmath.h"
00022 #include "geomNode.h"
00023 #include "geom.h"
00024 #include "geomVertexData.h"
00025 #include "geomVertexFormat.h"
00026 #include "geomPoints.h"
00027 #include "geomTriangles.h"
00028 #include "geomVertexWriter.h"
00029 #include "lens.h"
00030 #include "look_at.h"
00031 
00032 ////////////////////////////////////////////////////////////////////
00033 //     Function: PfmFile::Constructor
00034 //       Access: Published
00035 //  Description:
00036 ////////////////////////////////////////////////////////////////////
00037 PfmFile::
00038 PfmFile() {
00039   _zero_special = false;
00040   _vis_inverse = false;
00041   _vis_2d = false;
00042   clear();
00043 }
00044 
00045 ////////////////////////////////////////////////////////////////////
00046 //     Function: PfmFile::Copy Constructor
00047 //       Access: Published
00048 //  Description:
00049 ////////////////////////////////////////////////////////////////////
00050 PfmFile::
00051 PfmFile(const PfmFile &copy) :
00052   _table(copy._table),
00053   _x_size(copy._x_size),
00054   _y_size(copy._y_size),
00055   _scale(copy._scale),
00056   _num_channels(copy._num_channels),
00057   _zero_special(copy._zero_special),
00058   _vis_inverse(copy._vis_inverse),
00059   _vis_2d(copy._vis_2d)
00060 {
00061 }
00062 
00063 ////////////////////////////////////////////////////////////////////
00064 //     Function: PfmFile::Copy Assignment
00065 //       Access: Published
00066 //  Description:
00067 ////////////////////////////////////////////////////////////////////
00068 void PfmFile::
00069 operator = (const PfmFile &copy) {
00070   _table = copy._table;
00071   _x_size = copy._x_size;
00072   _y_size = copy._y_size;
00073   _scale = copy._scale;
00074   _num_channels = copy._num_channels;
00075   _zero_special = copy._zero_special;
00076   _vis_inverse = copy._vis_inverse;
00077   _vis_2d = copy._vis_2d;
00078 }
00079 
00080 ////////////////////////////////////////////////////////////////////
00081 //     Function: PfmFile::clear
00082 //       Access: Published
00083 //  Description: Eliminates all data in the file.
00084 ////////////////////////////////////////////////////////////////////
00085 void PfmFile::
00086 clear() {
00087   _x_size = 0;
00088   _y_size = 0;
00089   _scale = 1.0;
00090   _num_channels = 3;
00091   _table.clear();
00092 }
00093 
00094 ////////////////////////////////////////////////////////////////////
00095 //     Function: PfmFile::clear
00096 //       Access: Published
00097 //  Description: Resets to an empty table with a specific size.
00098 ////////////////////////////////////////////////////////////////////
00099 void PfmFile::
00100 clear(int x_size, int y_size, int num_channels) {
00101   nassertv(num_channels == 1 || num_channels == 3);
00102   nassertv(x_size >= 0 && y_size >= 0);
00103   _x_size = x_size;
00104   _y_size = y_size;
00105   _scale = 1.0;
00106   _num_channels = _num_channels;
00107 
00108   _table.clear();
00109   int size = _x_size * _y_size;
00110   _table.insert(_table.end(), size, LPoint3::zero());
00111 }
00112 
00113 ////////////////////////////////////////////////////////////////////
00114 //     Function: PfmFile::read
00115 //       Access: Published
00116 //  Description: Reads the PFM data from the indicated file, returning
00117 //               true on success, false on failure.
00118 ////////////////////////////////////////////////////////////////////
00119 bool PfmFile::
00120 read(const Filename &fullpath) {
00121   VirtualFileSystem *vfs = VirtualFileSystem::get_global_ptr();
00122 
00123   Filename filename = Filename::binary_filename(fullpath);
00124   PT(VirtualFile) file = vfs->get_file(filename);
00125   if (file == (VirtualFile *)NULL) {
00126     // No such file.
00127     grutil_cat.error()
00128       << "Could not find " << fullpath << "\n";
00129     return false;
00130   }
00131 
00132   if (grutil_cat.is_debug()) {
00133     grutil_cat.debug()
00134       << "Reading PFM file " << filename << "\n";
00135   }
00136 
00137   istream *in = file->open_read_file(true);
00138   bool success = read(*in);
00139   vfs->close_read_file(in);
00140 
00141   return success;
00142 }
00143 
00144 ////////////////////////////////////////////////////////////////////
00145 //     Function: PfmFile::read
00146 //       Access: Published
00147 //  Description: Reads the PFM data from the indicated stream,
00148 //               returning true on success, false on failure.
00149 ////////////////////////////////////////////////////////////////////
00150 bool PfmFile::
00151 read(istream &in) {
00152   clear();
00153 
00154   string identifier;
00155   in >> identifier;
00156 
00157   if (identifier == "PF") {
00158     _num_channels = 3;
00159   } else if (identifier == "Pf") {
00160     _num_channels = 1;
00161   } else {
00162     grutil_cat.error()
00163       << "Not a pfm file.\n";
00164     return false;
00165   }
00166 
00167   int width, height;
00168   PN_stdfloat scale;
00169   in >> width >> height >> scale;
00170   if (!in) {
00171     grutil_cat.error()
00172       << "Error parsing pfm header.\n";
00173     return false;
00174   }
00175 
00176   // Skip the last newline/whitespace character before the raw data
00177   // begins.
00178   in.get();
00179 
00180   bool little_endian = false;
00181   if (scale < 0) {
00182     scale = -scale;
00183     little_endian = true;
00184   }
00185   if (pfm_force_littleendian) {
00186     little_endian = true;
00187   }
00188   if (pfm_reverse_dimensions) {
00189     int t = width;
00190     width = height;
00191     height = t;
00192   }
00193 
00194   _x_size = width;
00195   _y_size = height;
00196   _scale = scale;
00197 
00198   // So far, so good.  Now read the data.
00199   int size = _x_size * _y_size;
00200   _table.reserve(size);
00201 
00202   if (little_endian) {
00203     for (int i = 0; i < size; ++i) {
00204       LPoint3 point = LPoint3::zero();
00205       for (int ci = 0; ci < _num_channels; ++ci) {
00206         PN_float32 data;
00207         in.read((char *)&data, sizeof(data));
00208         LittleEndian value(&data, sizeof(data));
00209         PN_float32 result;
00210         value.store_value(&result, sizeof(result));
00211         if (!cnan(result)) {
00212           point[ci] = result;
00213         }
00214       }
00215       _table.push_back(point);
00216     }
00217   } else {
00218     for (int i = 0; i < size; ++i) {
00219       LPoint3 point = LPoint3::zero();
00220       for (int ci = 0; ci < _num_channels; ++ci) {
00221         PN_float32 data;
00222         in.read((char *)&data, sizeof(data));
00223         BigEndian value(&data, sizeof(data));
00224         PN_float32 result;
00225         value.store_value(&result, sizeof(result));
00226         if (!cnan(result)) {
00227           point[ci] = result;
00228         }
00229       }
00230       _table.push_back(point);
00231     }
00232   }
00233 
00234   if (in.fail() && !in.eof()) {
00235     return false;
00236   }
00237 
00238   nassertr(sizeof(PN_float32) == 4, false);
00239   return true;
00240 }
00241 
00242 ////////////////////////////////////////////////////////////////////
00243 //     Function: PfmFile::write
00244 //       Access: Published
00245 //  Description: Writes the PFM data to the indicated file, returning
00246 //               true on success, false on failure.
00247 ////////////////////////////////////////////////////////////////////
00248 bool PfmFile::
00249 write(const Filename &fullpath) {
00250   Filename filename = Filename::binary_filename(fullpath);
00251   pofstream out;
00252   if (!filename.open_write(out)) {
00253     grutil_cat.error()
00254       << "Unable to open " << filename << "\n";
00255     return false;
00256   }
00257 
00258   if (grutil_cat.is_debug()) {
00259     grutil_cat.debug()
00260       << "Writing PFM file " << filename << "\n";
00261   }
00262 
00263   return write(out);
00264 }
00265 
00266 ////////////////////////////////////////////////////////////////////
00267 //     Function: PfmFile::write
00268 //       Access: Published
00269 //  Description: Writes the PFM data to the indicated stream,
00270 //               returning true on success, false on failure.
00271 ////////////////////////////////////////////////////////////////////
00272 bool PfmFile::
00273 write(ostream &out) {
00274   nassertr(is_valid(), false);
00275 
00276   if (_num_channels == 1) {
00277     out << "Pf\n";
00278   } else {
00279     out << "PF\n";
00280   }
00281   out << _x_size << " " << _y_size << "\n";
00282 
00283   PN_stdfloat scale = cabs(_scale);
00284   if (scale == 0.0f) {
00285     scale = 1.0f;
00286   }
00287 #ifndef WORDS_BIGENDIAN
00288   // Little-endian computers must write a negative scale to indicate
00289   // the little-endian nature of the output.
00290   scale = -scale;
00291 #endif
00292   out << scale << "\n";
00293 
00294   int size = _x_size * _y_size;
00295   for (int i = 0; i < size; ++i) {
00296     const LPoint3 &point = _table[i];
00297     for (int ci = 0; ci < _num_channels; ++ci) {
00298       PN_float32 data = point[ci];
00299       out.write((const char *)&data, sizeof(data));
00300     }
00301   }
00302 
00303   if (out.fail()) {
00304     return false;
00305   }
00306   nassertr(sizeof(PN_float32) == 4, false);
00307   return true;
00308 }
00309 
00310 ////////////////////////////////////////////////////////////////////
00311 //     Function: PfmFile::calc_average_point
00312 //       Access: Published
00313 //  Description: Computes the unweighted average point of all points
00314 //               within the box centered at (x, y) with the indicated
00315 //               Manhattan-distance radius.  Missing points are
00316 //               assigned the value of their nearest neighbor.
00317 //               Returns true if successful, or false if the point
00318 //               value cannot be determined.
00319 ////////////////////////////////////////////////////////////////////
00320 bool PfmFile::
00321 calc_average_point(LPoint3 &result, PN_stdfloat x, PN_stdfloat y, PN_stdfloat radius) const {
00322   result = LPoint3::zero();
00323 
00324   int min_x = int(ceil(x - radius));
00325   int min_y = int(ceil(y - radius));
00326   int max_x = int(floor(x + radius));
00327   int max_y = int(floor(y + radius));
00328 
00329   // We first construct a mini-grid of x_size by y_size integer values
00330   // to index into the main table.  This indirection allows us to fill
00331   // in the holes in the mini-grid with the nearest known values
00332   // before we compute the average.
00333   int x_size = max_x - min_x + 1;
00334   int y_size = max_y - min_y + 1;
00335   int size = x_size * y_size;
00336   if (size == 0) {
00337     return false;
00338   }
00339 
00340   pvector<MiniGridCell> mini_grid;
00341   mini_grid.insert(mini_grid.end(), size, MiniGridCell());
00342 
00343   // Now collect the known data points and apply them to the
00344   // mini-grid.
00345   min_x = max(min_x, 0);
00346   min_y = max(min_y, 0);
00347   max_x = min(max_x, _x_size - 1);
00348   max_y = min(max_y, _y_size - 1);
00349 
00350   bool got_any = false;
00351   int xi, yi;
00352   for (yi = min_y; yi <= max_y; ++yi) {
00353     for (xi = min_x; xi <= max_x; ++xi) {
00354       const LPoint3 &p = _table[yi * _x_size + xi];
00355       if (_zero_special && p == LPoint3::zero()) {
00356         continue;
00357       }
00358 
00359       int gi = (yi - min_y) * y_size + (xi - min_x);
00360       nassertr(gi >= 0 && gi < size, false);
00361       mini_grid[gi]._ti = yi * _x_size + xi;
00362       mini_grid[gi]._dist = 0;
00363       got_any = true;
00364     }
00365   }
00366 
00367   if (!got_any) {
00368     return false;
00369   }
00370 
00371   // Now recursively fill in any missing holes in the mini-grid.
00372   for (yi = 0; yi < y_size; ++yi) {
00373     for (xi = 0; xi < x_size; ++xi) {
00374       int gi = yi * x_size + xi;
00375       if (mini_grid[gi]._dist == 0) {
00376         int ti = mini_grid[gi]._ti;
00377         fill_mini_grid(&mini_grid[0], x_size, y_size, xi + 1, yi, 1, ti);
00378         fill_mini_grid(&mini_grid[0], x_size, y_size, xi - 1, yi, 1, ti);
00379         fill_mini_grid(&mini_grid[0], x_size, y_size, xi, yi + 1, 1, ti);
00380         fill_mini_grid(&mini_grid[0], x_size, y_size, xi, yi - 1, 1, ti);
00381       }
00382     }
00383   }
00384 
00385   // Now the mini-grid is completely filled, so we can compute the
00386   // average.
00387   for (int gi = 0; gi < size; ++gi) {
00388     int ti = mini_grid[gi]._ti;
00389     nassertr(ti >= 0 && ti < (int)_table.size(), false);
00390     result += _table[ti];
00391   }
00392 
00393   result /= PN_stdfloat(size);
00394   return true;
00395 }
00396 
00397 ////////////////////////////////////////////////////////////////////
00398 //     Function: PfmFile::calc_min_max
00399 //       Access: Published
00400 //  Description: Calculates the minimum and maximum x, y, and z depth
00401 //               component values, representing the bounding box of
00402 //               depth values, and places them in the indicated
00403 //               vectors.  Returns true if successful, false if the
00404 //               mesh contains no points.
00405 ////////////////////////////////////////////////////////////////////
00406 bool PfmFile::
00407 calc_min_max(LVecBase3 &min_depth, LVecBase3 &max_depth) const {
00408   bool any_points = false;
00409 
00410   min_depth = LVecBase3::zero();
00411   max_depth = LVecBase3::zero();
00412 
00413   Table::const_iterator ti;
00414   for (ti = _table.begin(); ti != _table.end(); ++ti) {
00415     const LPoint3 &p = (*ti);
00416     if (_zero_special && p == LPoint3::zero()) {
00417       continue;
00418     }
00419     
00420     if (!any_points) {
00421       min_depth = p;
00422       max_depth = p;
00423       any_points = true;
00424     } else {
00425       min_depth[0] = min(min_depth[0], p[0]);
00426       min_depth[1] = min(min_depth[1], p[1]);
00427       min_depth[2] = min(min_depth[2], p[2]);
00428       max_depth[0] = max(max_depth[0], p[0]);
00429       max_depth[1] = max(max_depth[1], p[1]);
00430       max_depth[2] = max(max_depth[2], p[2]);
00431     }
00432   }
00433 
00434   return any_points;
00435 }
00436 
00437 ////////////////////////////////////////////////////////////////////
00438 //     Function: PfmFile::resize
00439 //       Access: Published
00440 //  Description: Applies a simple filter to resample the pfm file
00441 //               in-place to the indicated size.  Don't confuse this
00442 //               with applying a scale to all of the points via
00443 //               xform().
00444 ////////////////////////////////////////////////////////////////////
00445 void PfmFile::
00446 resize(int new_x_size, int new_y_size) {
00447   if (_x_size == 0 || _y_size == 0 || new_x_size == 0 || new_y_size == 0) {
00448     clear(new_x_size, new_y_size, _num_channels);
00449     return;
00450   }
00451 
00452   if (new_x_size == _x_size && new_y_size == _y_size) {
00453     return;
00454   }
00455 
00456   Table new_data;
00457   new_data.reserve(new_x_size * new_y_size);
00458 
00459   PN_stdfloat from_x0, from_x1, from_y0, from_y1;
00460 
00461   PN_stdfloat x_scale = 1.0;
00462   PN_stdfloat y_scale = 1.0;
00463 
00464   if (new_x_size > 1) {
00465     x_scale = (PN_stdfloat)(_x_size - 1) / (PN_stdfloat)(new_x_size - 1);
00466   }
00467   if (new_y_size > 1) {
00468     y_scale = (PN_stdfloat)(_y_size - 1) / (PN_stdfloat)(new_y_size - 1);
00469   }
00470 
00471   from_y0 = 0.0;
00472   for (int to_y = 0; to_y < new_y_size; ++to_y) {
00473     from_y1 = (to_y + 0.5) * y_scale;
00474     from_y1 = min(from_y1, (PN_stdfloat) _y_size);
00475 
00476     from_x0 = 0.0;
00477     for (int to_x = 0; to_x < new_x_size; ++to_x) {
00478       from_x1 = (to_x + 0.5) * x_scale;
00479       from_x1 = min(from_x1, (PN_stdfloat) _x_size);
00480 
00481       // Now the box from (from_x0, from_y0) - (from_x1, from_y1)
00482       // but not including (from_x1, from_y1) maps to the pixel (to_x, to_y).
00483       LPoint3 result;
00484       box_filter_region(result, from_x0, from_y0, from_x1, from_y1);
00485       new_data.push_back(result);
00486 
00487       from_x0 = from_x1;
00488     }
00489     from_y0 = from_y1;
00490   }
00491 
00492   _table.swap(new_data);
00493   _x_size = new_x_size;
00494   _y_size = new_y_size;
00495 }
00496 
00497 ////////////////////////////////////////////////////////////////////
00498 //     Function: PfmFile::reverse_rows
00499 //       Access: Published
00500 //  Description: Performs an in-place reversal of the row (y) data.
00501 ////////////////////////////////////////////////////////////////////
00502 void PfmFile::
00503 reverse_rows() {
00504   nassertv(is_valid());
00505 
00506   Table reversed;
00507   reversed.reserve(_table.size());
00508   for (int yi = 0; yi < _y_size; ++yi) {
00509     int source_yi = _y_size - 1 - yi;
00510     int start = source_yi * _x_size;
00511     reversed.insert(reversed.end(), 
00512                     _table.begin() + start, _table.begin() + start + _x_size);
00513   }
00514 
00515   nassertv(reversed.size() == _table.size());
00516   _table.swap(reversed);
00517 }
00518 
00519 ////////////////////////////////////////////////////////////////////
00520 //     Function: PfmFile::xform
00521 //       Access: Published
00522 //  Description: Applies the indicated transform matrix to all points
00523 //               in-place.
00524 ////////////////////////////////////////////////////////////////////
00525 void PfmFile::
00526 xform(const LMatrix4 &transform) {
00527   nassertv(is_valid());
00528 
00529   Table::iterator ti;
00530   for (ti = _table.begin(); ti != _table.end(); ++ti) {
00531     if (_zero_special && (*ti) == LPoint3::zero()) {
00532       continue;
00533     }
00534 
00535     (*ti) = (*ti) * transform;
00536   }
00537 }
00538 
00539 ////////////////////////////////////////////////////////////////////
00540 //     Function: PfmFile::project
00541 //       Access: Published
00542 //  Description: Adjusts each (x, y, z) point of the Pfm file by
00543 //               projecting it through the indicated lens, converting
00544 //               each point to a (u, v, 0) texture coordinate.  The
00545 //               resulting file can be generated to a mesh (with
00546 //               set_vis_inverse(true) and generate_vis_mesh())
00547 //               that will apply the lens distortion to an arbitrary
00548 //               texture image.
00549 ////////////////////////////////////////////////////////////////////
00550 void PfmFile::
00551 project(const Lens *lens) {
00552   nassertv(is_valid());
00553 
00554   static LMatrix4 to_uv(0.5, 0.0, 0.0, 0.0,
00555                         0.0, 0.5, 0.0, 0.0, 
00556                         0.0, 0.0, 0.0, 0.0, 
00557                         0.5, 0.5, 0.0, 1.0);
00558   
00559   Table::iterator ti;
00560   for (ti = _table.begin(); ti != _table.end(); ++ti) {
00561     if (_zero_special && (*ti) == LPoint3::zero()) {
00562       continue;
00563     }
00564 
00565     LPoint3 &p = (*ti);
00566     LPoint3 film;
00567     lens->project(p, film);
00568     p = to_uv.xform_point(film);
00569   }
00570 }
00571 
00572 ////////////////////////////////////////////////////////////////////
00573 //     Function: PfmFile::merge
00574 //       Access: Published
00575 //  Description: Wherever there is missing data in this PfmFile (that
00576 //               is, wherever has_point() returns false), copy data
00577 //               from the other PfmFile, which must be exactly the
00578 //               same dimensions as this one.
00579 ////////////////////////////////////////////////////////////////////
00580 void PfmFile::
00581 merge(const PfmFile &other) {
00582   nassertv(is_valid() && other.is_valid());
00583   nassertv(other._x_size == _x_size && other._y_size == _y_size);
00584 
00585   if (!_zero_special) {
00586     // Trivial no-op.
00587     return;
00588   }
00589 
00590   for (size_t i = 0; i < _table.size(); ++i) {
00591     if (_table[i] == LPoint3::zero()) {
00592       _table[i] = other._table[i];
00593     }
00594   }
00595 }
00596 
00597 ////////////////////////////////////////////////////////////////////
00598 //     Function: PfmFile::compute_planar_bounds
00599 //       Access: Published
00600 //  Description: This version of this method exists for temporary
00601 //               backward compatibility only.
00602 ////////////////////////////////////////////////////////////////////
00603 PT(BoundingHexahedron) PfmFile::
00604 compute_planar_bounds(PN_stdfloat point_dist, PN_stdfloat sample_radius) const {
00605   return compute_planar_bounds(LPoint2(0.5, 0.5), point_dist, sample_radius, false);
00606 }
00607 
00608 
00609 ////////////////////////////////////////////////////////////////////
00610 //     Function: PfmFile::compute_planar_bounds
00611 //       Access: Published
00612 //  Description: Computes the minmax bounding volume of the points in
00613 //               3-D space, assuming the points represent a
00614 //               mostly-planar surface.
00615 //
00616 //               This algorithm works by sampling the (square)
00617 //               sample_radius pixels at the four point_dist corners
00618 //               around the center (cx - pd, cx + pd) and so on, to
00619 //               approximate the plane of the surface.  Then all of
00620 //               the points are projected into that plane and the
00621 //               bounding volume of the entire mesh within that plane
00622 //               is determined.  If points_only is true, the bounding
00623 //               volume of only those four points is determined.
00624 //
00625 //               center, point_dist and sample_radius are in UV space,
00626 //               i.e. in the range 0..1.
00627 ////////////////////////////////////////////////////////////////////
00628 PT(BoundingHexahedron) PfmFile::
00629 compute_planar_bounds(const LPoint2 &center, PN_stdfloat point_dist, PN_stdfloat sample_radius, bool points_only) const {
00630   LPoint3 p0, p1, p2, p3;
00631   compute_sample_point(p0, center[0] + point_dist, center[1] - point_dist, sample_radius);
00632   compute_sample_point(p1, center[0] + point_dist, center[1] + point_dist, sample_radius);
00633   compute_sample_point(p2, center[0] - point_dist, center[1] + point_dist, sample_radius);
00634   compute_sample_point(p3, center[0] - point_dist, center[1] - point_dist, sample_radius);
00635 
00636   LPoint3 normal;
00637 
00638   normal[0] = p0[1] * p1[2] - p0[2] * p1[1];
00639   normal[1] = p0[2] * p1[0] - p0[0] * p1[2];
00640   normal[2] = p0[0] * p1[1] - p0[1] * p1[0];
00641 
00642   normal[0] += p1[1] * p2[2] - p1[2] * p2[1];
00643   normal[1] += p1[2] * p2[0] - p1[0] * p2[2];
00644   normal[2] += p1[0] * p2[1] - p1[1] * p2[0];
00645 
00646   normal[0] += p2[1] * p3[2] - p2[2] * p3[1];
00647   normal[1] += p2[2] * p3[0] - p2[0] * p3[2];
00648   normal[2] += p2[0] * p3[1] - p2[1] * p3[0];
00649 
00650   normal[0] += p3[1] * p0[2] - p3[2] * p0[1];
00651   normal[1] += p3[2] * p0[0] - p3[0] * p0[2];
00652   normal[2] += p3[0] * p0[1] - p3[1] * p0[0];
00653 
00654   normal.normalize();
00655 
00656   LVector3 up = (p1 - p0) + (p2 - p3);
00657   LPoint3 pcenter = ((p0 + p1 + p2 + p3) * 0.25);
00658 
00659   // Compute the transform necessary to rotate all of the points into
00660   // the Y = 0 plane.
00661   LMatrix4 rotate;
00662   look_at(rotate, normal, up);
00663 
00664   LMatrix4 rinv;
00665   rinv.invert_from(rotate);
00666 
00667   LPoint3 trans = pcenter * rinv;
00668   rinv.set_row(3, -trans);
00669   rotate.invert_from(rinv);
00670 
00671   // Now determine the minmax.
00672   PN_stdfloat min_x, min_y, min_z, max_x, max_y, max_z;
00673   bool got_point = false;
00674   if (points_only) {
00675     LPoint3 points[4] = {
00676       p0 * rinv,
00677       p1 * rinv,
00678       p2 * rinv,
00679       p3 * rinv,
00680     };
00681     for (int i = 0; i < 4; ++i) {
00682       const LPoint3 &point = points[i];
00683       if (!got_point) {
00684         min_x = point[0];
00685         min_y = point[1];
00686         min_z = point[2];
00687         max_x = point[0];
00688         max_y = point[1];
00689         max_z = point[2];
00690         got_point = true;
00691       } else {
00692         min_x = min(min_x, point[0]);
00693         min_y = min(min_y, point[1]);
00694         min_z = min(min_z, point[2]);
00695         max_x = max(max_x, point[0]);
00696         max_y = max(max_y, point[1]);
00697         max_z = max(max_z, point[2]);
00698       }
00699     }
00700 
00701   } else {
00702     Table::const_iterator ti;
00703     for (ti = _table.begin(); ti != _table.end(); ++ti) {
00704       if (_zero_special && (*ti) == LPoint3::zero()) {
00705         continue;
00706       }
00707       LPoint3 point = (*ti) * rinv;
00708       if (!got_point) {
00709         min_x = point[0];
00710         min_y = point[1];
00711         min_z = point[2];
00712         max_x = point[0];
00713         max_y = point[1];
00714         max_z = point[2];
00715         got_point = true;
00716       } else {
00717         min_x = min(min_x, point[0]);
00718         min_y = min(min_y, point[1]);
00719         min_z = min(min_z, point[2]);
00720         max_x = max(max_x, point[0]);
00721         max_y = max(max_y, point[1]);
00722         max_z = max(max_z, point[2]);
00723       }
00724     }
00725   }
00726     
00727   PT(BoundingHexahedron) bounds = new BoundingHexahedron
00728     (LPoint3(min_x, min_y, min_z), LPoint3(max_x, min_y, min_z),
00729      LPoint3(min_x, min_y, max_z), LPoint3(max_x, min_y, max_z),
00730      LPoint3(min_x, max_y, min_z), LPoint3(max_x, max_y, min_z),
00731      LPoint3(min_x, max_y, max_z), LPoint3(max_x, max_y, max_z));
00732 
00733   // Rotate the bounding volume back into the original space of the
00734   // screen.
00735   bounds->xform(rotate);
00736 
00737   return bounds;
00738 }
00739 
00740 ////////////////////////////////////////////////////////////////////
00741 //     Function: PfmFile::compute_sample_point
00742 //       Access: Published
00743 //  Description: Computes the average of all the point within
00744 //               sample_radius (manhattan distance) and the indicated
00745 //               point.
00746 //
00747 //               The point coordinates are given in UV space, in the
00748 //               range 0..1.
00749 ////////////////////////////////////////////////////////////////////
00750 void PfmFile::
00751 compute_sample_point(LPoint3 &result,
00752                      PN_stdfloat x, PN_stdfloat y, PN_stdfloat sample_radius) const {
00753   x *= _x_size;
00754   y *= _y_size;
00755   PN_stdfloat xr = sample_radius * _x_size;
00756   PN_stdfloat yr = sample_radius * _y_size;
00757   box_filter_region(result, x - xr, y - yr, x + xr, y + yr);
00758 }
00759 
00760 
00761 ////////////////////////////////////////////////////////////////////
00762 //     Function: PfmFile::generate_vis_points
00763 //       Access: Published
00764 //  Description: Creates a point cloud with the points of the pfm as
00765 //               3-d coordinates in space, and texture coordinates
00766 //               ranging from 0 .. 1 based on the position within the
00767 //               pfm grid.
00768 ////////////////////////////////////////////////////////////////////
00769 NodePath PfmFile::
00770 generate_vis_points() const {
00771   nassertr(is_valid(), NodePath());
00772 
00773   CPT(GeomVertexFormat) format;
00774   if (_vis_inverse) {
00775     if (_vis_2d) {
00776       format = GeomVertexFormat::get_v3t2();
00777     } else {
00778       // We need a 3-d texture coordinate if we're inverting the vis
00779       // and it's 3-d.
00780       GeomVertexArrayFormat *v3t3 = new GeomVertexArrayFormat
00781         (InternalName::get_vertex(), 3, 
00782          Geom::NT_stdfloat, Geom::C_point,
00783          InternalName::get_texcoord(), 3, 
00784          Geom::NT_stdfloat, Geom::C_texcoord);
00785       format = GeomVertexFormat::register_format(v3t3);
00786     }
00787   } else {
00788     format = GeomVertexFormat::get_v3t2();
00789   }
00790 
00791   PT(GeomVertexData) vdata = new GeomVertexData
00792     ("points", format, Geom::UH_static);
00793   vdata->set_num_rows(_x_size * _y_size);
00794   GeomVertexWriter vertex(vdata, InternalName::get_vertex());
00795   GeomVertexWriter texcoord(vdata, InternalName::get_texcoord());
00796 
00797   LPoint2 uv_scale(1.0, 1.0);
00798   if (_x_size > 1) {
00799     uv_scale[0] = 1.0f / PN_stdfloat(_x_size - 1);
00800   }
00801   if (_y_size > 1) {
00802     uv_scale[1] = 1.0f / PN_stdfloat(_y_size - 1);
00803   }
00804   
00805   for (int yi = 0; yi < _y_size; ++yi) {
00806     for (int xi = 0; xi < _x_size; ++xi) {
00807       const LPoint3 &point = get_point(xi, yi);
00808       LPoint2 uv(PN_stdfloat(xi) * uv_scale[0],
00809                  PN_stdfloat(yi) * uv_scale[1]);
00810       if (_vis_inverse) {
00811         vertex.add_data2(uv);
00812         texcoord.add_data3(point);
00813       } else if (_vis_2d) {
00814         vertex.add_data2(point[0], point[1]);
00815         texcoord.add_data2(uv);
00816       } else {
00817         vertex.add_data3(point);
00818         texcoord.add_data2(uv);
00819       }
00820     }
00821   }
00822   
00823   PT(Geom) geom = new Geom(vdata);
00824   PT(GeomPoints) points = new GeomPoints(Geom::UH_static);
00825   points->add_next_vertices(_x_size * _y_size);
00826   geom->add_primitive(points);
00827   
00828   PT(GeomNode) gnode = new GeomNode("");
00829   gnode->add_geom(geom);
00830   return NodePath(gnode);
00831 }
00832 
00833 ////////////////////////////////////////////////////////////////////
00834 //     Function: PfmFile::generate_vis_mesh
00835 //       Access: Published
00836 //  Description: Creates a triangle mesh with the points of the pfm as
00837 //               3-d coordinates in space, and texture coordinates
00838 //               ranging from 0 .. 1 based on the position within the
00839 //               pfm grid.
00840 ////////////////////////////////////////////////////////////////////
00841 NodePath PfmFile::
00842 generate_vis_mesh(MeshFace face) const {
00843   nassertr(is_valid(), NodePath());
00844   nassertr(face != 0, NodePath());
00845 
00846   if (_x_size == 1 || _y_size == 1) {
00847     // Can't generate a 1-d mesh, so generate points in this case.
00848     return generate_vis_points();
00849   }
00850   
00851   PT(GeomNode) gnode = new GeomNode("");
00852 
00853   if (face & MF_front) {
00854     make_vis_mesh_geom(gnode, false);
00855   }
00856 
00857   if (face & MF_back) {
00858     make_vis_mesh_geom(gnode, true);
00859   }
00860 
00861   return NodePath(gnode);
00862 }
00863 
00864 ////////////////////////////////////////////////////////////////////
00865 //     Function: PfmFile::make_vis_mesh_geom
00866 //       Access: Private
00867 //  Description: Returns a triangle mesh for the pfm.  If inverted is
00868 //               true, the mesh is facing the opposite direction.
00869 ////////////////////////////////////////////////////////////////////
00870 void PfmFile::
00871 make_vis_mesh_geom(GeomNode *gnode, bool inverted) const {
00872   int num_x_cells = 1;
00873   int num_y_cells = 1;
00874 
00875   int x_size = _x_size;
00876   int y_size = _y_size;
00877 
00878   // This is the number of independent vertices we will require.
00879   int num_vertices = x_size * y_size;
00880   if (num_vertices == 0) {
00881     // Trivial no-op.
00882     return;
00883   }
00884 
00885   // This is the max number of vertex indices we might add to the
00886   // GeomTriangles.  (We might actually add fewer than this due to
00887   // omitting the occasional missing data point.)
00888   int max_indices = (x_size - 1) * (y_size - 1) * 6;
00889 
00890   while (num_vertices > pfm_vis_max_vertices || max_indices > pfm_vis_max_indices) {
00891     // Too many vertices in one mesh.  Subdivide the mesh into smaller
00892     // pieces.
00893     if (num_x_cells > num_y_cells) {
00894       ++num_y_cells;
00895     } else {
00896       ++num_x_cells;
00897     }
00898 
00899     x_size = (_x_size + num_x_cells - 1) / num_x_cells + 1;
00900     y_size = (_y_size + num_y_cells - 1) / num_y_cells + 1;
00901 
00902     num_vertices = x_size * y_size;
00903     max_indices = (x_size - 1) * (y_size - 1) * 6;
00904   }
00905 
00906   // OK, now we know how many cells we need.
00907   if (grutil_cat.is_debug()) {
00908     grutil_cat.debug()
00909       << "Generating mesh with " << num_x_cells << " x " << num_y_cells
00910       << " pieces.\n";
00911   }
00912 
00913   PT(GeomVertexArrayFormat) array_format;
00914 
00915   if (_vis_2d) {
00916     // No normals needed if we're just generating a 2-d mesh.
00917     array_format = new GeomVertexArrayFormat
00918       (InternalName::get_vertex(), 3, Geom::NT_stdfloat, Geom::C_point,
00919        InternalName::get_texcoord(), 2, Geom::NT_stdfloat, Geom::C_texcoord);
00920 
00921   } else {
00922     if (_vis_inverse) {
00923       // We need a 3-d texture coordinate if we're inverting the vis
00924       // and it's 3-d.  But we still don't need normals in that case.
00925       array_format = new GeomVertexArrayFormat
00926         (InternalName::get_vertex(), 3, Geom::NT_stdfloat, Geom::C_point,
00927          InternalName::get_texcoord(), 3, Geom::NT_stdfloat, Geom::C_texcoord);
00928     } else {
00929       // Otherwise, we only need a 2-d texture coordinate, and we do
00930       // want normals.
00931       array_format = new GeomVertexArrayFormat
00932         (InternalName::get_vertex(), 3, Geom::NT_stdfloat, Geom::C_point,
00933          InternalName::get_normal(), 3, Geom::NT_stdfloat, Geom::C_vector,
00934          InternalName::get_texcoord(), 2, Geom::NT_stdfloat, Geom::C_texcoord);
00935     }
00936   }
00937 
00938   if (_flat_texcoord_name != (InternalName *)NULL) {
00939     // We need an additional texcoord column for the flat texcoords.
00940     array_format->add_column(_flat_texcoord_name, 2, 
00941                              Geom::NT_stdfloat, Geom::C_texcoord);
00942   }
00943 
00944   CPT(GeomVertexFormat) format = GeomVertexFormat::register_format(array_format);
00945 
00946   for (int yci = 0; yci < num_y_cells; ++yci) {
00947     int y_begin = (yci * _y_size) / num_y_cells;
00948     int y_end = ((yci + 1) * _y_size) / num_y_cells;
00949 
00950     // Include the first vertex from the next strip in this strip's
00951     // vertices, so we are connected.
00952     y_end = min(y_end + 1, _y_size);
00953 
00954     y_size = y_end - y_begin;
00955     if (y_size == 0) {
00956       continue;
00957     }
00958 
00959     for (int xci = 0; xci < num_x_cells; ++xci) {
00960       int x_begin = (xci * _x_size) / num_x_cells;
00961       int x_end = ((xci + 1) * _x_size) / num_x_cells;
00962       x_end = min(x_end + 1, _x_size);
00963       x_size = x_end - x_begin;
00964       if (x_size == 0) {
00965         continue;
00966       }
00967 
00968       num_vertices = x_size * y_size;
00969       max_indices = (x_size - 1) * (y_size - 1) * 6;
00970 
00971       ostringstream mesh_name;
00972       mesh_name << "mesh_" << xci << "_" << yci;
00973       PT(GeomVertexData) vdata = new GeomVertexData
00974         (mesh_name.str(), format, Geom::UH_static);
00975 
00976       vdata->set_num_rows(num_vertices);
00977       GeomVertexWriter vertex(vdata, InternalName::get_vertex());
00978       GeomVertexWriter normal(vdata, InternalName::get_normal());
00979       GeomVertexWriter texcoord(vdata, InternalName::get_texcoord());
00980       GeomVertexWriter texcoord2(vdata, _flat_texcoord_name);
00981 
00982       for (int yi = y_begin; yi < y_end; ++yi) {
00983         for (int xi = x_begin; xi < x_end; ++xi) {
00984           const LPoint3 &point = get_point(xi, yi);
00985           LPoint2 uv(PN_stdfloat(xi) / PN_stdfloat(_x_size - 1),
00986                      PN_stdfloat(yi) / PN_stdfloat(_y_size - 1));
00987 
00988           if (_vis_inverse) {
00989             vertex.add_data2(uv);
00990             texcoord.add_data3(point);
00991           } else if (_vis_2d) {
00992             vertex.add_data2(point[0], point[1]);
00993             texcoord.add_data2(uv);
00994           } else {
00995             vertex.add_data3(point);
00996             texcoord.add_data2(uv);
00997             
00998             // Calculate the normal based on two neighboring vertices.
00999             LPoint3 v[3];
01000             v[0] = get_point(xi, yi);
01001             if (xi + 1 < _x_size) {
01002               v[1] = get_point(xi + 1, yi);
01003             } else {
01004               v[1] = v[0];
01005               v[0] = get_point(xi - 1, yi);
01006             }
01007             
01008             if (yi + 1 < _y_size) {
01009               v[2] = get_point(xi, yi + 1);
01010             } else {
01011               v[2] = v[0];
01012               v[0] = get_point(xi, yi - 1);
01013             }
01014         
01015             LVector3 n = LVector3::zero();
01016             for (int i = 0; i < 3; ++i) {
01017               const LPoint3 &v0 = v[i];
01018               const LPoint3 &v1 = v[(i + 1) % 3];
01019               n[0] += v0[1] * v1[2] - v0[2] * v1[1];
01020               n[1] += v0[2] * v1[0] - v0[0] * v1[2];
01021               n[2] += v0[0] * v1[1] - v0[1] * v1[0];
01022             }
01023             n.normalize();
01024             nassertv(!n.is_nan());
01025             if (inverted) {
01026               n = -n;
01027             }
01028             normal.add_data3(n);
01029           }
01030 
01031           if (_flat_texcoord_name != (InternalName *)NULL) {
01032             texcoord2.add_data2(uv);
01033           }
01034         }
01035       }
01036   
01037       PT(Geom) geom = new Geom(vdata);
01038       PT(GeomTriangles) tris = new GeomTriangles(Geom::UH_static);
01039 
01040       tris->reserve_num_vertices(max_indices);
01041 
01042       for (int yi = y_begin; yi < y_end - 1; ++yi) {
01043         for (int xi = x_begin; xi < x_end - 1; ++xi) {
01044 
01045           if (_zero_special) {
01046             if (get_point(xi, yi) == LPoint3::zero() ||
01047                 get_point(xi, yi + 1) == LPoint3::zero() ||
01048                 get_point(xi + 1, yi + 1) == LPoint3::zero() ||
01049                 get_point(xi + 1, yi) == LPoint3::zero()) {
01050               continue;
01051             }
01052           }
01053 
01054           int xi0 = xi - x_begin;
01055           int yi0 = yi - y_begin;
01056 
01057           int vi0 = ((xi0) + (yi0) * x_size);
01058           int vi1 = ((xi0) + (yi0 + 1) * x_size);
01059           int vi2 = ((xi0 + 1) + (yi0 + 1) * x_size);
01060           int vi3 = ((xi0 + 1) + (yi0) * x_size);
01061           
01062           if (inverted) {
01063             tris->add_vertices(vi2, vi0, vi1);
01064             tris->close_primitive();
01065             
01066             tris->add_vertices(vi3, vi0, vi2);
01067             tris->close_primitive();
01068           } else {
01069             tris->add_vertices(vi2, vi1, vi0);
01070             tris->close_primitive();
01071             
01072             tris->add_vertices(vi3, vi2, vi0);
01073             tris->close_primitive();
01074           }
01075         }
01076       }
01077       geom->add_primitive(tris);
01078       gnode->add_geom(geom);
01079     }
01080   }
01081 }
01082 
01083 
01084 ////////////////////////////////////////////////////////////////////
01085 //     Function: PfmFile::box_filter_region
01086 //       Access: Private
01087 //  Description: Averages all the points in the rectangle from x0
01088 //               .. y0 to x1 .. y1 into result.  The region may be
01089 //               defined by floating-point boundaries; the result will
01090 //               be weighted by the degree of coverage of each
01091 //               included point.
01092 ////////////////////////////////////////////////////////////////////
01093 void PfmFile::
01094 box_filter_region(LPoint3 &result,
01095                   PN_stdfloat x0, PN_stdfloat y0, PN_stdfloat x1, PN_stdfloat y1) const {
01096   result = LPoint3::zero();
01097   PN_stdfloat coverage = 0.0;
01098 
01099   if (x1 < x0 || y1 < y0) {
01100     return;
01101   }
01102   nassertv(y0 >= 0.0 && y1 >= 0.0);
01103 
01104   int y = (int)y0;
01105   // Get the first (partial) row
01106   box_filter_line(result, coverage, x0, y, x1, (PN_stdfloat)(y+1)-y0);
01107 
01108   int y_last = (int)y1;
01109   if (y < y_last) {
01110     y++;
01111     while (y < y_last) {
01112       // Get each consecutive (complete) row
01113       box_filter_line(result, coverage, x0, y, x1, 1.0);
01114       y++;
01115     }
01116 
01117     // Get the final (partial) row
01118     PN_stdfloat y_contrib = y1 - (PN_stdfloat)y_last;
01119     if (y_contrib > 0.0001) {
01120       box_filter_line(result, coverage, x0, y, x1, y_contrib);
01121     }
01122   }
01123 
01124   if (coverage != 0.0) {
01125     result /= coverage;
01126   }
01127 }
01128 
01129 ////////////////////////////////////////////////////////////////////
01130 //     Function: PfmFile::box_filter_line
01131 //       Access: Private
01132 //  Description: 
01133 ////////////////////////////////////////////////////////////////////
01134 void PfmFile::
01135 box_filter_line(LPoint3 &result, PN_stdfloat &coverage,
01136                 PN_stdfloat x0, int y, PN_stdfloat x1, PN_stdfloat y_contrib) const {
01137   int x = (int)x0;
01138   // Get the first (partial) xel
01139   box_filter_point(result, coverage, x, y, (PN_stdfloat)(x+1)-x0, y_contrib);
01140 
01141   int x_last = (int)x1;
01142   if (x < x_last) {
01143     x++;
01144     while (x < x_last) {
01145       // Get each consecutive (complete) xel
01146       box_filter_point(result, coverage, x, y, 1.0, y_contrib);
01147       x++;
01148     }
01149 
01150     // Get the final (partial) xel
01151     PN_stdfloat x_contrib = x1 - (PN_stdfloat)x_last;
01152     if (x_contrib > 0.0001) {
01153       box_filter_point(result, coverage, x, y, x_contrib, y_contrib);
01154     }
01155   }
01156 }
01157 
01158 ////////////////////////////////////////////////////////////////////
01159 //     Function: PfmFile::box_filter_point
01160 //       Access: Private
01161 //  Description: 
01162 ////////////////////////////////////////////////////////////////////
01163 void PfmFile::
01164 box_filter_point(LPoint3 &result, PN_stdfloat &coverage,
01165                  int x, int y, PN_stdfloat x_contrib, PN_stdfloat y_contrib) const {
01166   const LPoint3 &point = get_point(x, y);
01167   if (_zero_special && point == LPoint3::zero()) {
01168     return;
01169   }
01170 
01171   PN_stdfloat contrib = x_contrib * y_contrib;
01172   result += point * contrib;
01173   coverage += contrib;
01174 }
01175 
01176 ////////////////////////////////////////////////////////////////////
01177 //     Function: PfmFile::fill_mini_grid
01178 //       Access: Private
01179 //  Description: A support function for calc_average_point(), this
01180 //               recursively fills in the holes in the mini_grid data
01181 //               with the index to the nearest value.
01182 ////////////////////////////////////////////////////////////////////
01183 void PfmFile::
01184 fill_mini_grid(MiniGridCell *mini_grid, int x_size, int y_size, 
01185                int xi, int yi, int dist, int ti) const {
01186   if (xi < 0 || xi >= x_size || yi < 0 || yi >= y_size) {
01187     // Out of bounds.
01188     return;
01189   }
01190 
01191   int gi = yi * x_size + xi;
01192   if (mini_grid[gi]._dist == -1 || mini_grid[gi]._dist > dist) {
01193     // Here's an undefined value that we need to populate.
01194     mini_grid[gi]._dist = dist;
01195     mini_grid[gi]._ti = ti;
01196     fill_mini_grid(mini_grid, x_size, y_size, xi + 1, yi, dist + 1, ti);
01197     fill_mini_grid(mini_grid, x_size, y_size, xi - 1, yi, dist + 1, ti);
01198     fill_mini_grid(mini_grid, x_size, y_size, xi, yi + 1, dist + 1, ti);
01199     fill_mini_grid(mini_grid, x_size, y_size, xi, yi - 1, dist + 1, ti);
01200   }
01201 }
 All Classes Functions Variables Enumerations