Field3D
|
00001 //----------------------------------------------------------------------------// 00002 00003 /* 00004 * Copyright (c) 2009 Sony Pictures Imageworks Inc 00005 * 00006 * All rights reserved. 00007 * 00008 * Redistribution and use in source and binary forms, with or without 00009 * modification, are permitted provided that the following conditions 00010 * are met: 00011 * 00012 * Redistributions of source code must retain the above copyright 00013 * notice, this list of conditions and the following disclaimer. 00014 * Redistributions in binary form must reproduce the above copyright 00015 * notice, this list of conditions and the following disclaimer in the 00016 * documentation and/or other materials provided with the 00017 * distribution. Neither the name of Sony Pictures Imageworks nor the 00018 * names of its contributors may be used to endorse or promote 00019 * products derived from this software without specific prior written 00020 * permission. 00021 * 00022 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS 00023 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT 00024 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS 00025 * FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE 00026 * COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, 00027 * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES 00028 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR 00029 * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 00030 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, 00031 * STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) 00032 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED 00033 * OF THE POSSIBILITY OF SUCH DAMAGE. 00034 */ 00035 00036 //----------------------------------------------------------------------------// 00037 00043 //----------------------------------------------------------------------------// 00044 00045 #ifndef _INCLUDED_Field3D_MACFieldUtil_H_ 00046 #define _INCLUDED_Field3D_MACFieldUtil_H_ 00047 00048 #include "MACField.h" 00049 00050 //----------------------------------------------------------------------------// 00051 00052 #include "ns.h" 00053 00054 FIELD3D_NAMESPACE_OPEN 00055 00056 //----------------------------------------------------------------------------// 00057 // Utility functions 00058 //----------------------------------------------------------------------------// 00059 00061 template <class Data_T, class Field_T> 00062 void convertMACToCellCentered(typename MACField<Data_T>::Ptr mac, 00063 typename Field_T::Ptr cc); 00064 00065 //----------------------------------------------------------------------------// 00066 00068 template <class Field_T, class Data_T> 00069 void convertCellCenteredToMAC(typename Field_T::Ptr cc, 00070 typename MACField<Data_T>::Ptr mac); 00071 00072 //----------------------------------------------------------------------------// 00073 // Implementations 00074 //----------------------------------------------------------------------------// 00075 00077 template <class Data_T, class Field_T> 00078 void convertMACToCellCentered(typename MACField<Data_T>::Ptr mac, 00079 typename Field_T::Ptr cc) 00080 { 00081 // Make sure the extents and data window match 00082 if (cc->extents().min != mac->extents().min || 00083 cc->extents().max != mac->extents().max || 00084 cc->dataWindow().min != mac->dataWindow().min || 00085 cc->dataWindow().max != mac->dataWindow().max ) { 00086 cc->setSize(mac->extents(), mac->dataWindow()); 00087 } 00088 00089 // Make sure mapping matches 00090 if (!cc->mapping()->isIdentical(mac->mapping())) { 00091 cc->setMapping(mac->mapping()); 00092 } 00093 00094 // MAC velocities are in simulation space (axis-aligned to the 00095 // mapping) because the values are stored on the faces, so rotate 00096 // vectors from simulation-space to world-space when transferring 00097 // from MAC to cell-centered 00098 00099 bool rotateVector = false; 00100 M44d ssToWsMtx; 00101 MatrixFieldMapping::Ptr mapping = 00102 field_dynamic_cast<MatrixFieldMapping>(mac->mapping()); 00103 if (mapping) { 00104 M44d localToWorldMtx = mapping->localToWorld(); 00105 V3d scale, rot, trans, shear; 00106 if (extractSHRT(localToWorldMtx, scale, shear, rot, trans, false)) { 00107 ssToWsMtx.rotate(rot); 00108 if (rot.length2() > FLT_EPSILON) 00109 rotateVector = true; 00110 } 00111 } 00112 00113 // Loop over all the voxels in the output field --- 00114 00115 typename Field_T::iterator i = cc->begin(); 00116 typename Field_T::iterator end = cc->end(); 00117 00118 if (rotateVector) { 00119 for (; i != end; ++i) { 00120 *i = mac->value(i.x, i.y, i.z) * ssToWsMtx; 00121 } 00122 } else { 00123 for (; i != end; ++i) { 00124 *i = mac->value(i.x, i.y, i.z); 00125 } 00126 } 00127 } 00128 00129 //----------------------------------------------------------------------------// 00130 00131 template <class Field_T, class Data_T> 00132 void convertCellCenteredToMAC(typename Field_T::Ptr cc, 00133 typename MACField<Data_T>::Ptr mac) 00134 { 00135 // Make sure the extents and data window match 00136 if (mac->extents().min != cc->extents().min || 00137 mac->extents().max != cc->extents().max || 00138 mac->dataWindow().min != cc->dataWindow().min || 00139 mac->dataWindow().max != cc->dataWindow().max ) { 00140 mac->setSize(cc->extents(), cc->dataWindow()); 00141 } 00142 00143 // Make sure mapping matches 00144 if (!mac->mapping()->isIdentical(cc->mapping())) { 00145 mac->setMapping(cc->mapping()); 00146 } 00147 00148 Box3i data = mac->dataWindow(); 00149 00150 // MAC velocities are in simulation space (axis-aligned to the 00151 // mapping) because the values are stored on the faces, so rotate 00152 // vectors from world-space to simulation-space when transferring 00153 // from cell-centered to MAC 00154 00155 bool rotateVector = false; 00156 M44d wsToSsMtx; 00157 MatrixFieldMapping::Ptr mapping = 00158 field_dynamic_cast<MatrixFieldMapping>(mac->mapping()); 00159 if (mapping) { 00160 M44d localToWorld = mapping->localToWorld(); 00161 V3d scale, rot, trans, shear; 00162 if (FIELD3D_EXTRACT_SHRT(localToWorld, scale, shear, rot, trans, false)) { 00163 wsToSsMtx.rotate(-rot); 00164 rotateVector = true; 00165 } 00166 } 00167 00168 // Use a pointer to a field below so we can substitute it out for 00169 // our intermediate, rotated field if necessary, without needing to 00170 // duplicate the loops. This should be more efficient CPU-wise so 00171 // we don't need to do 3 matrix multiplies per cell-centered voxel 00172 // because it's used in 3 separate loops (1 per MAC face). 00173 typename Field_T::Ptr src = cc; 00174 00175 typename Field_T::Ptr ccRotated; 00176 if (rotateVector) { 00177 ccRotated = 00178 typename Field_T::Ptr(new Field_T); 00179 ccRotated->matchDefinition(cc); 00180 00181 typename Field_T::const_iterator iIn = cc->cbegin(); 00182 typename Field_T::const_iterator endIn = cc->cend(); 00183 typename Field_T::iterator iOut = ccRotated->begin(); 00184 00185 for (; iIn != endIn; ++iIn, ++iOut) { 00186 *iOut = *iIn * wsToSsMtx; 00187 } 00188 src = ccRotated; 00189 } 00190 00191 // Set the u edge value to their closest voxel 00192 for (int k = data.min.z; k <= data.max.z; k++) { 00193 for (int j = data.min.y; j <= data.max.y; j++) { 00194 mac->u(data.min.x, j, k) = src->value(data.min.x, j, k).x; 00195 mac->u(data.max.x + 1, j, k) = src->value(data.max.x, j, k).x; 00196 } 00197 } 00198 00199 // Set the v edge value to their closest voxel 00200 for (int k = data.min.z; k <= data.max.z; k++) { 00201 for (int i = data.min.x; i <= data.max.x; i++) { 00202 mac->v(i, data.min.y, k) = src->value(i, data.min.y, k).y; 00203 mac->v(i, data.max.y + 1, k) = src->value(i, data.max.y, k).y; 00204 } 00205 } 00206 00207 // Set the w edge value to their closest voxel 00208 for (int j = data.min.y; j <= data.max.y; j++) { 00209 for (int i = data.min.x; i <= data.max.x; i++) { 00210 mac->w(i, j, data.min.z) = src->value(i, j, data.min.z).z; 00211 mac->w(i, j, data.max.z + 1) = src->value(i, j, data.max.z).z; 00212 } 00213 } 00214 00215 // Loop over internal u values 00216 for (int k = data.min.z; k <= data.max.z; ++k) { 00217 for (int j = data.min.y; j <= data.max.y; ++j) { 00218 for (int i = data.min.x + 1; i <= data.max.x; ++i) { 00219 mac->u(i, j, k) = 00220 (src->value(i, j, k).x + src->value(i - 1, j, k).x) * 0.5; 00221 } 00222 } 00223 } 00224 00225 // Loop over internal v values 00226 for (int k = data.min.z; k <= data.max.z; ++k) { 00227 for (int j = data.min.y + 1; j <= data.max.y; ++j) { 00228 for (int i = data.min.x; i <= data.max.x; ++i) { 00229 mac->v(i, j, k) = 00230 (src->value(i, j, k).y + src->value(i, j - 1, k).y) * 0.5; 00231 } 00232 } 00233 } 00234 00235 // Loop over internal w values 00236 for (int k = data.min.z + 1; k <= data.max.z; ++k) { 00237 for (int j = data.min.y; j <= data.max.y; ++j) { 00238 for (int i = data.min.x; i <= data.max.x; ++i) { 00239 mac->w(i, j, k) = 00240 (src->value(i, j, k).z + src->value(i, j, k - 1).z) * 0.5; 00241 } 00242 } 00243 } 00244 } 00245 00246 //----------------------------------------------------------------------------// 00247 00248 FIELD3D_NAMESPACE_HEADER_CLOSE 00249 00250 //----------------------------------------------------------------------------// 00251 00252 #endif // Include guard