| 1 | /* | 
| 2 | * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved. | 
| 3 | * | 
| 4 | * The University of Notre Dame grants you ("Licensee") a | 
| 5 | * non-exclusive, royalty free, license to use, modify and | 
| 6 | * redistribute this software in source and binary code form, provided | 
| 7 | * that the following conditions are met: | 
| 8 | * | 
| 9 | * 1. Redistributions of source code must retain the above copyright | 
| 10 | *    notice, this list of conditions and the following disclaimer. | 
| 11 | * | 
| 12 | * 2. Redistributions in binary form must reproduce the above copyright | 
| 13 | *    notice, this list of conditions and the following disclaimer in the | 
| 14 | *    documentation and/or other materials provided with the | 
| 15 | *    distribution. | 
| 16 | * | 
| 17 | * This software is provided "AS IS," without a warranty of any | 
| 18 | * kind. All express or implied conditions, representations and | 
| 19 | * warranties, including any implied warranty of merchantability, | 
| 20 | * fitness for a particular purpose or non-infringement, are hereby | 
| 21 | * excluded.  The University of Notre Dame and its licensors shall not | 
| 22 | * be liable for any damages suffered by licensee as a result of | 
| 23 | * using, modifying or distributing the software or its | 
| 24 | * derivatives. In no event will the University of Notre Dame or its | 
| 25 | * licensors be liable for any lost revenue, profit or data, or for | 
| 26 | * direct, indirect, special, consequential, incidental or punitive | 
| 27 | * damages, however caused and regardless of the theory of liability, | 
| 28 | * arising out of the use of or inability to use software, even if the | 
| 29 | * University of Notre Dame has been advised of the possibility of | 
| 30 | * such damages. | 
| 31 | * | 
| 32 | * SUPPORT OPEN SCIENCE!  If you use OpenMD or its source code in your | 
| 33 | * research, please cite the appropriate papers when you publish your | 
| 34 | * work.  Good starting points are: | 
| 35 | * | 
| 36 | * [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005). | 
| 37 | * [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006). | 
| 38 | * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008). | 
| 39 | * [4] Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010). | 
| 40 | */ | 
| 41 |  | 
| 42 |  | 
| 43 | #include <algorithm> | 
| 44 | #include <fstream> | 
| 45 | #include "applications/staticProps/RNEMDStats.hpp" | 
| 46 | #include "primitives/Molecule.hpp" | 
| 47 | #include "utils/PhysicalConstants.hpp" | 
| 48 |  | 
| 49 | namespace OpenMD { | 
| 50 |  | 
| 51 | RNEMDZ::RNEMDZ(SimInfo* info, const std::string& filename, | 
| 52 | const std::string& sele, int nzbins) | 
| 53 | : SlabStatistics(info, filename, sele, nzbins) { | 
| 54 |  | 
| 55 | setOutputName(getPrefix(filename) + ".rnemdZ"); | 
| 56 |  | 
| 57 | temperature = new OutputData; | 
| 58 | temperature->units =  "K"; | 
| 59 | temperature->title =  "Temperature"; | 
| 60 | temperature->dataType = odtReal; | 
| 61 | temperature->dataHandling = odhAverage; | 
| 62 | temperature->accumulator.reserve(nBins_); | 
| 63 | for (int i = 0; i < nBins_; i++) | 
| 64 | temperature->accumulator.push_back( new Accumulator() ); | 
| 65 | data_.push_back(temperature); | 
| 66 |  | 
| 67 | velocity = new OutputData; | 
| 68 | velocity->units = "angstroms/fs"; | 
| 69 | velocity->title =  "Velocity"; | 
| 70 | velocity->dataType = odtVector3; | 
| 71 | velocity->dataHandling = odhAverage; | 
| 72 | velocity->accumulator.reserve(nBins_); | 
| 73 | for (int i = 0; i < nBins_; i++) | 
| 74 | velocity->accumulator.push_back( new VectorAccumulator() ); | 
| 75 | data_.push_back(velocity); | 
| 76 |  | 
| 77 | density = new OutputData; | 
| 78 | density->units =  "g cm^-3"; | 
| 79 | density->title =  "Density"; | 
| 80 | density->dataType = odtReal; | 
| 81 | density->dataHandling = odhAverage; | 
| 82 | density->accumulator.reserve(nBins_); | 
| 83 | for (int i = 0; i < nBins_; i++) | 
| 84 | density->accumulator.push_back( new Accumulator() ); | 
| 85 | data_.push_back(density); | 
| 86 | } | 
| 87 |  | 
| 88 | void RNEMDZ::processFrame(int istep) { | 
| 89 | RealType z; | 
| 90 |  | 
| 91 | hmat_ = currentSnapshot_->getHmat(); | 
| 92 | for (int i = 0; i < nBins_; i++) { | 
| 93 | z = (((RealType)i + 0.5) / (RealType)nBins_) * hmat_(2,2); | 
| 94 | dynamic_cast<Accumulator*>(z_->accumulator[i])->add(z); | 
| 95 | } | 
| 96 | volume_ = currentSnapshot_->getVolume(); | 
| 97 |  | 
| 98 |  | 
| 99 | Molecule* mol; | 
| 100 | RigidBody* rb; | 
| 101 | StuntDouble* sd; | 
| 102 | SimInfo::MoleculeIterator mi; | 
| 103 | Molecule::RigidBodyIterator rbIter; | 
| 104 | int i; | 
| 105 |  | 
| 106 | vector<RealType> binMass(nBins_, 0.0); | 
| 107 | vector<Vector3d> binP(nBins_, V3Zero); | 
| 108 | vector<RealType> binKE(nBins_, 0.0); | 
| 109 | vector<unsigned int> binDof(nBins_, 0); | 
| 110 |  | 
| 111 | for (mol = info_->beginMolecule(mi); mol != NULL; | 
| 112 | mol = info_->nextMolecule(mi)) { | 
| 113 |  | 
| 114 | // change the positions of atoms which belong to the rigidbodies | 
| 115 |  | 
| 116 | for (rb = mol->beginRigidBody(rbIter); rb != NULL; | 
| 117 | rb = mol->nextRigidBody(rbIter)) { | 
| 118 | rb->updateAtomVel(); | 
| 119 | } | 
| 120 | } | 
| 121 |  | 
| 122 | if (evaluator_.isDynamic()) { | 
| 123 | seleMan_.setSelectionSet(evaluator_.evaluate()); | 
| 124 | } | 
| 125 |  | 
| 126 | // loop over the selected atoms: | 
| 127 |  | 
| 128 | for (sd = seleMan_.beginSelected(i); sd != NULL; | 
| 129 | sd = seleMan_.nextSelected(i)) { | 
| 130 |  | 
| 131 | // figure out where that object is: | 
| 132 | Vector3d pos = sd->getPos(); | 
| 133 | Vector3d vel = sd->getVel(); | 
| 134 | RealType m = sd->getMass(); | 
| 135 |  | 
| 136 | int bin = getBin(pos); | 
| 137 |  | 
| 138 | binMass[bin] += m; | 
| 139 | binP[bin] += m * vel; | 
| 140 | binKE[bin] += 0.5 * (m * vel.lengthSquare()); | 
| 141 | binDof[bin] += 3; | 
| 142 |  | 
| 143 | if (sd->isDirectional()) { | 
| 144 | Vector3d angMom = sd->getJ(); | 
| 145 | Mat3x3d I = sd->getI(); | 
| 146 | if (sd->isLinear()) { | 
| 147 | int i = sd->linearAxis(); | 
| 148 | int j = (i + 1) % 3; | 
| 149 | int k = (i + 2) % 3; | 
| 150 | binKE[bin] += 0.5 * (angMom[j] * angMom[j] / I(j, j) + | 
| 151 | angMom[k] * angMom[k] / I(k, k)); | 
| 152 | binDof[bin] += 2; | 
| 153 | } else { | 
| 154 | binKE[bin] += 0.5 * (angMom[0] * angMom[0] / I(0, 0) + | 
| 155 | angMom[1] * angMom[1] / I(1, 1) + | 
| 156 | angMom[2] * angMom[2] / I(2, 2)); | 
| 157 | binDof[bin] += 3; | 
| 158 | } | 
| 159 | } | 
| 160 | } | 
| 161 |  | 
| 162 | for (int i = 0; i < nBins_; i++) { | 
| 163 |  | 
| 164 | if (binDof[i] > 0) { | 
| 165 | RealType temp = 2.0 * binKE[i] / (binDof[i] * PhysicalConstants::kb * | 
| 166 | PhysicalConstants::energyConvert); | 
| 167 | RealType den = binMass[i] * nBins_ * PhysicalConstants::densityConvert | 
| 168 | / volume_; | 
| 169 | Vector3d vel = binP[i] / binMass[i]; | 
| 170 |  | 
| 171 | dynamic_cast<Accumulator *>(temperature->accumulator[i])->add(temp); | 
| 172 | dynamic_cast<VectorAccumulator *>(velocity->accumulator[i])->add(vel); | 
| 173 | dynamic_cast<Accumulator *>(density->accumulator[i])->add(den); | 
| 174 | dynamic_cast<Accumulator *>(counts_->accumulator[i])->add(1); | 
| 175 | } | 
| 176 | } | 
| 177 | } | 
| 178 |  | 
| 179 | void RNEMDZ::processStuntDouble(StuntDouble* sd, int bin) { | 
| 180 | } | 
| 181 |  | 
| 182 | RNEMDR::RNEMDR(SimInfo* info, const std::string& filename, | 
| 183 | const std::string& sele, int nrbins) | 
| 184 | : ShellStatistics(info, filename, sele, nrbins) { | 
| 185 |  | 
| 186 |  | 
| 187 | setOutputName(getPrefix(filename) + ".rnemdR"); | 
| 188 |  | 
| 189 | temperature = new OutputData; | 
| 190 | temperature->units =  "K"; | 
| 191 | temperature->title =  "Temperature"; | 
| 192 | temperature->dataType = odtReal; | 
| 193 | temperature->dataHandling = odhAverage; | 
| 194 | temperature->accumulator.reserve(nBins_); | 
| 195 | for (int i = 0; i < nBins_; i++) | 
| 196 | temperature->accumulator.push_back( new Accumulator() ); | 
| 197 | data_.push_back(temperature); | 
| 198 |  | 
| 199 | angularVelocity = new OutputData; | 
| 200 | angularVelocity->units = "angstroms^2/fs"; | 
| 201 | angularVelocity->title =  "Angular Velocity"; | 
| 202 | angularVelocity->dataType = odtVector3; | 
| 203 | angularVelocity->dataHandling = odhAverage; | 
| 204 | angularVelocity->accumulator.reserve(nBins_); | 
| 205 | for (int i = 0; i < nBins_; i++) | 
| 206 | angularVelocity->accumulator.push_back( new VectorAccumulator() ); | 
| 207 | data_.push_back(angularVelocity); | 
| 208 |  | 
| 209 | density = new OutputData; | 
| 210 | density->units =  "g cm^-3"; | 
| 211 | density->title =  "Density"; | 
| 212 | density->dataType = odtReal; | 
| 213 | density->dataHandling = odhAverage; | 
| 214 | density->accumulator.reserve(nBins_); | 
| 215 | for (int i = 0; i < nBins_; i++) | 
| 216 | density->accumulator.push_back( new Accumulator() ); | 
| 217 | data_.push_back(density); | 
| 218 | } | 
| 219 |  | 
| 220 |  | 
| 221 | void RNEMDR::processFrame(int istep) { | 
| 222 |  | 
| 223 | Molecule* mol; | 
| 224 | RigidBody* rb; | 
| 225 | StuntDouble* sd; | 
| 226 | SimInfo::MoleculeIterator mi; | 
| 227 | Molecule::RigidBodyIterator rbIter; | 
| 228 | int i; | 
| 229 |  | 
| 230 | vector<RealType> binMass(nBins_, 0.0); | 
| 231 | vector<Mat3x3d>  binI(nBins_); | 
| 232 | vector<Vector3d> binL(nBins_, V3Zero); | 
| 233 | vector<RealType> binKE(nBins_, 0.0); | 
| 234 | vector<int> binDof(nBins_, 0); | 
| 235 |  | 
| 236 | for (mol = info_->beginMolecule(mi); mol != NULL; | 
| 237 | mol = info_->nextMolecule(mi)) { | 
| 238 |  | 
| 239 | // change the positions of atoms which belong to the rigidbodies | 
| 240 |  | 
| 241 | for (rb = mol->beginRigidBody(rbIter); rb != NULL; | 
| 242 | rb = mol->nextRigidBody(rbIter)) { | 
| 243 | rb->updateAtomVel(); | 
| 244 | } | 
| 245 | } | 
| 246 |  | 
| 247 | if (evaluator_.isDynamic()) { | 
| 248 | seleMan_.setSelectionSet(evaluator_.evaluate()); | 
| 249 | } | 
| 250 |  | 
| 251 | // loop over the selected atoms: | 
| 252 |  | 
| 253 | for (sd = seleMan_.beginSelected(i); sd != NULL; | 
| 254 | sd = seleMan_.nextSelected(i)) { | 
| 255 |  | 
| 256 | // figure out where that object is: | 
| 257 | int bin = getBin( sd->getPos() ); | 
| 258 |  | 
| 259 | if (bin >= 0 && bin < nBins_)  { | 
| 260 |  | 
| 261 | Vector3d rPos = sd->getPos() - coordinateOrigin_; | 
| 262 | Vector3d vel = sd->getVel(); | 
| 263 | RealType m = sd->getMass(); | 
| 264 | Vector3d L = m * cross(rPos, vel); | 
| 265 | Mat3x3d I(0.0); | 
| 266 | I = outProduct(rPos, rPos) * m; | 
| 267 | RealType r2 = rPos.lengthSquare(); | 
| 268 | I(0, 0) += m * r2; | 
| 269 | I(1, 1) += m * r2; | 
| 270 | I(2, 2) += m * r2; | 
| 271 |  | 
| 272 | binMass[bin] += m; | 
| 273 | binI[bin] += I; | 
| 274 | binL[bin] += L; | 
| 275 | binKE[bin] += 0.5 * (m * vel.lengthSquare()); | 
| 276 | binDof[bin] += 3; | 
| 277 |  | 
| 278 | if (sd->isDirectional()) { | 
| 279 | Vector3d angMom = sd->getJ(); | 
| 280 | Mat3x3d Ia = sd->getI(); | 
| 281 | if (sd->isLinear()) { | 
| 282 | int i = sd->linearAxis(); | 
| 283 | int j = (i + 1) % 3; | 
| 284 | int k = (i + 2) % 3; | 
| 285 | binKE[bin] += 0.5 * (angMom[j] * angMom[j] / Ia(j, j) + | 
| 286 | angMom[k] * angMom[k] / Ia(k, k)); | 
| 287 | binDof[bin] += 2; | 
| 288 | } else { | 
| 289 | binKE[bin] += 0.5 * (angMom[0] * angMom[0] / Ia(0, 0) + | 
| 290 | angMom[1] * angMom[1] / Ia(1, 1) + | 
| 291 | angMom[2] * angMom[2] / Ia(2, 2)); | 
| 292 | binDof[bin] += 3; | 
| 293 | } | 
| 294 | } | 
| 295 | } | 
| 296 | } | 
| 297 |  | 
| 298 | for (int i = 0; i < nBins_; i++) { | 
| 299 | RealType rinner = (RealType)i * binWidth_; | 
| 300 | RealType router = (RealType)(i+1) * binWidth_; | 
| 301 | if (binDof[i] > 0) { | 
| 302 | RealType temp = 2.0 * binKE[i] / (binDof[i] * PhysicalConstants::kb * | 
| 303 | PhysicalConstants::energyConvert); | 
| 304 | RealType den = binMass[i] * 3.0 * PhysicalConstants::densityConvert | 
| 305 | / (4.0 * M_PI * (pow(router,3) - pow(rinner,3))); | 
| 306 |  | 
| 307 | Vector3d omega = binI[i].inverse() * binL[i]; | 
| 308 |  | 
| 309 | dynamic_cast<Accumulator *>(temperature->accumulator[i])->add(temp); | 
| 310 | dynamic_cast<VectorAccumulator *>(angularVelocity->accumulator[i])->add(omega); | 
| 311 | dynamic_cast<Accumulator *>(density->accumulator[i])->add(den); | 
| 312 | dynamic_cast<Accumulator *>(counts_->accumulator[i])->add(1); | 
| 313 | } | 
| 314 | } | 
| 315 | } | 
| 316 |  | 
| 317 |  | 
| 318 | void RNEMDR::processStuntDouble(StuntDouble* sd, int bin) { | 
| 319 | } | 
| 320 |  | 
| 321 | RNEMDRTheta::RNEMDRTheta(SimInfo* info, const std::string& filename, | 
| 322 | const std::string& sele, int nrbins, int nangleBins) | 
| 323 | : ShellStatistics(info, filename, sele, nrbins), nAngleBins_(nangleBins) { | 
| 324 |  | 
| 325 | Globals* simParams = info->getSimParams(); | 
| 326 | RNEMDParameters* rnemdParams = simParams->getRNEMDParameters(); | 
| 327 | bool hasAngularMomentumFluxVector = rnemdParams->haveAngularMomentumFluxVector(); | 
| 328 |  | 
| 329 | if (hasAngularMomentumFluxVector) { | 
| 330 | std::vector<RealType> amf = rnemdParams->getAngularMomentumFluxVector(); | 
| 331 | if (amf.size() != 3) { | 
| 332 | sprintf(painCave.errMsg, | 
| 333 | "RNEMDRTheta: Incorrect number of parameters specified for angularMomentumFluxVector.\n" | 
| 334 | "\tthere should be 3 parameters, but %lu were specified.\n", | 
| 335 | amf.size()); | 
| 336 | painCave.isFatal = 1; | 
| 337 | simError(); | 
| 338 | } | 
| 339 | fluxVector_.x() = amf[0]; | 
| 340 | fluxVector_.y() = amf[1]; | 
| 341 | fluxVector_.z() = amf[2]; | 
| 342 | } else { | 
| 343 |  | 
| 344 | std::string fluxStr = rnemdParams->getFluxType(); | 
| 345 | if (fluxStr.find("Lx") != std::string::npos) { | 
| 346 | fluxVector_ = V3X; | 
| 347 | } else if (fluxStr.find("Ly") != std::string::npos) { | 
| 348 | fluxVector_ = V3Y; | 
| 349 | } else { | 
| 350 | fluxVector_ = V3Z; | 
| 351 | } | 
| 352 | } | 
| 353 |  | 
| 354 | fluxVector_.normalize(); | 
| 355 |  | 
| 356 | setOutputName(getPrefix(filename) + ".rnemdRTheta"); | 
| 357 |  | 
| 358 | angularVelocity = new OutputData; | 
| 359 | angularVelocity->units = "angstroms^2/fs"; | 
| 360 | angularVelocity->title =  "Angular Velocity"; | 
| 361 | angularVelocity->dataType = odtArray2d; | 
| 362 | angularVelocity->dataHandling = odhAverage; | 
| 363 | angularVelocity->accumulatorArray2d.reserve(nBins_); | 
| 364 | for (int i = 0; i < nBins_; i++) { | 
| 365 | angularVelocity->accumulatorArray2d[i].reserve(nAngleBins_); | 
| 366 | for (int j = 0 ; j < nAngleBins_; j++) { | 
| 367 | angularVelocity->accumulatorArray2d[i][j] = new Accumulator(); | 
| 368 | } | 
| 369 | } | 
| 370 | data_.push_back(angularVelocity); | 
| 371 |  | 
| 372 | } | 
| 373 |  | 
| 374 |  | 
| 375 | std::pair<int,int> RNEMDRTheta::getBins(Vector3d pos) { | 
| 376 | std::pair<int,int> result; | 
| 377 |  | 
| 378 | Vector3d rPos = pos - coordinateOrigin_; | 
| 379 | RealType cosAngle= dot(rPos, fluxVector_) / rPos.length(); | 
| 380 |  | 
| 381 | result.first = int(rPos.length() / binWidth_); | 
| 382 | result.second = int( (nAngleBins_ - 1) * 0.5 * (cosAngle + 1.0) ); | 
| 383 | return result; | 
| 384 | } | 
| 385 |  | 
| 386 | void RNEMDRTheta::processStuntDouble(StuntDouble* sd, int bin) { | 
| 387 | } | 
| 388 |  | 
| 389 | void RNEMDRTheta::processFrame(int istep) { | 
| 390 |  | 
| 391 | Molecule* mol; | 
| 392 | RigidBody* rb; | 
| 393 | StuntDouble* sd; | 
| 394 | SimInfo::MoleculeIterator mi; | 
| 395 | Molecule::RigidBodyIterator rbIter; | 
| 396 | int i; | 
| 397 |  | 
| 398 | vector<vector<Mat3x3d> >  binI; | 
| 399 | vector<vector<Vector3d> > binL; | 
| 400 | vector<vector<int> > binCount; | 
| 401 |  | 
| 402 | for (mol = info_->beginMolecule(mi); mol != NULL; | 
| 403 | mol = info_->nextMolecule(mi)) { | 
| 404 |  | 
| 405 | // change the positions of atoms which belong to the rigidbodies | 
| 406 |  | 
| 407 | for (rb = mol->beginRigidBody(rbIter); rb != NULL; | 
| 408 | rb = mol->nextRigidBody(rbIter)) { | 
| 409 | rb->updateAtomVel(); | 
| 410 | } | 
| 411 | } | 
| 412 |  | 
| 413 | if (evaluator_.isDynamic()) { | 
| 414 | seleMan_.setSelectionSet(evaluator_.evaluate()); | 
| 415 | } | 
| 416 |  | 
| 417 | // loop over the selected atoms: | 
| 418 |  | 
| 419 | for (sd = seleMan_.beginSelected(i); sd != NULL; | 
| 420 | sd = seleMan_.nextSelected(i)) { | 
| 421 |  | 
| 422 | // figure out where that object is: | 
| 423 | std::pair<int,int> bins = getBins( sd->getPos() ); | 
| 424 |  | 
| 425 | if (bins.first >= 0 && bins.first < nBins_)  { | 
| 426 | if (bins.second >= 0 && bins.second < nAngleBins_) { | 
| 427 |  | 
| 428 | Vector3d rPos = sd->getPos() - coordinateOrigin_; | 
| 429 | Vector3d vel = sd->getVel(); | 
| 430 | RealType m = sd->getMass(); | 
| 431 | Vector3d L = m * cross(rPos, vel); | 
| 432 | Mat3x3d I(0.0); | 
| 433 | I = outProduct(rPos, rPos) * m; | 
| 434 | RealType r2 = rPos.lengthSquare(); | 
| 435 | I(0, 0) += m * r2; | 
| 436 | I(1, 1) += m * r2; | 
| 437 | I(2, 2) += m * r2; | 
| 438 |  | 
| 439 | binI[bins.first][bins.second] += I; | 
| 440 | binL[bins.first][bins.second] += L; | 
| 441 | binCount[bins.first][bins.second]++; | 
| 442 | } | 
| 443 | } | 
| 444 | } | 
| 445 |  | 
| 446 |  | 
| 447 | for (int i = 0; i < nBins_; i++) { | 
| 448 | for (int j = 0; j < nAngleBins_; j++) { | 
| 449 |  | 
| 450 | if (binCount[i][j] > 0) { | 
| 451 | Vector3d omega = binI[i][j].inverse() * binL[i][j]; | 
| 452 | RealType omegaProj = dot(omega, fluxVector_); | 
| 453 |  | 
| 454 | dynamic_cast<Accumulator *>(angularVelocity->accumulatorArray2d[i][j])->add(omegaProj); | 
| 455 | } | 
| 456 | } | 
| 457 | } | 
| 458 | } | 
| 459 |  | 
| 460 | void RNEMDRTheta::writeOutput() { | 
| 461 |  | 
| 462 | vector<OutputData*>::iterator i; | 
| 463 | OutputData* outputData; | 
| 464 |  | 
| 465 | ofstream outStream(outputFilename_.c_str()); | 
| 466 | if (outStream.is_open()) { | 
| 467 |  | 
| 468 | //write title | 
| 469 | outStream << "# SPATIAL STATISTICS\n"; | 
| 470 | outStream << "#"; | 
| 471 |  | 
| 472 | for(outputData = beginOutputData(i); outputData; | 
| 473 | outputData = nextOutputData(i)) { | 
| 474 | outStream << "\t" << outputData->title << | 
| 475 | "(" << outputData->units << ")"; | 
| 476 | // add some extra tabs for column alignment | 
| 477 | if (outputData->dataType == odtVector3) outStream << "\t\t"; | 
| 478 | } | 
| 479 |  | 
| 480 | outStream << std::endl; | 
| 481 |  | 
| 482 | outStream.precision(8); | 
| 483 |  | 
| 484 | for (int j = 0; j < nBins_; j++) { | 
| 485 |  | 
| 486 | int counts = counts_->accumulator[j]->count(); | 
| 487 |  | 
| 488 | if (counts > 0) { | 
| 489 | for(outputData = beginOutputData(i); outputData; | 
| 490 | outputData = nextOutputData(i)) { | 
| 491 |  | 
| 492 | int n = outputData->accumulator[j]->count(); | 
| 493 | if (n != 0) { | 
| 494 | writeData( outStream, outputData, j ); | 
| 495 | } | 
| 496 | } | 
| 497 | outStream << std::endl; | 
| 498 | } | 
| 499 | } | 
| 500 | } | 
| 501 | } | 
| 502 | } | 
| 503 |  |