00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028 #include "treeNode.H"
00029 #include "octree.H"
00030 #include "treeLeaf.H"
00031 #include "treeBoundBox.H"
00032 #include <OpenFOAM/long.H>
00033 #include <OpenFOAM/linePointRef.H>
00034
00035
00036
00037 template <class Type>
00038 const Foam::label Foam::treeNode<Type>::leafOffset(100);
00039
00040
00041
00042
00043 template <class Type>
00044 void Foam::treeNode<Type>::setAsNode(const label octant)
00045 {
00046 subNodeTypes_ |= (0x1 << octant);
00047 }
00048
00049
00050 template <class Type>
00051 void Foam::treeNode<Type>::setAsLeaf(const label octant)
00052 {
00053 subNodeTypes_ &= ~(0x1 << octant);
00054 }
00055
00056
00057
00058 template <class Type>
00059 void Foam::treeNode<Type>::setNodePtr
00060 (
00061 const label octant,
00062 treeElem<Type>* treeNodePtr
00063 )
00064 {
00065 setAsNode(octant);
00066 subNodes_[octant] = treeNodePtr;
00067 }
00068
00069
00070
00071 template <class Type>
00072 void Foam::treeNode<Type>::setLeafPtr
00073 (
00074 const label octant,
00075 treeElem<Type>* treeLeafPtr
00076 )
00077 {
00078 setAsLeaf(octant);
00079 subNodes_[octant] = treeLeafPtr;
00080 }
00081
00082
00083 template <class Type>
00084 void Foam::treeNode<Type>::setVolType
00085 (
00086 const label octant,
00087 const label type
00088 )
00089 {
00090 if ((type < 0) || (type > 3))
00091 {
00092 FatalErrorIn("treeNode<Type>::setVolType(const label, const label)")
00093 << "Type " << type << " not within range 0..3" << endl;
00094 }
00095
00096
00097 volType_ &= ~(0x3 << 2*octant);
00098
00099
00100 volType_ |= (type << 2*octant);
00101 }
00102
00103
00104 template <class Type>
00105 void Foam::treeNode<Type>::space(Ostream& os, const label n)
00106 {
00107 for (label i=0; i<n; i++)
00108 {
00109 os<< ' ';
00110 }
00111 }
00112
00113
00114
00115 template <class Type>
00116 const Foam::treeLeaf<Type>* Foam::treeNode<Type>::findLeafLineOctant
00117 (
00118 const int level,
00119 const Type& shapes,
00120 const label octant,
00121 const vector& direction,
00122 point& start,
00123 const point& end
00124 ) const
00125 {
00126 static const char* functionName =
00127 "treeNode<Type>::findLeafLineOctant"
00128 "(const int, const Type&, const label, const vector&,"
00129 " point&, const point&)";
00130
00131 if (debug & 2)
00132 {
00133 space(Pout, 2*level);
00134 Pout<< "findLeafLineOctant : bb:" << this->bb()
00135 << " start:" << start
00136 << " end:" << end
00137 << " mid:" << midpoint()
00138 << " Searching octant:" << octant
00139 << endl;
00140 }
00141
00142 if (subNodes()[octant])
00143 {
00144 if (isNode(octant))
00145 {
00146
00147 const treeNode<Type>* subNodePtr = getNodePtr(octant);
00148
00149 if (subNodePtr->bb().contains(direction, start))
00150 {
00151
00152 const treeLeaf<Type>* subLeafPtr = subNodePtr->findLeafLine
00153 (
00154 level + 1,
00155 shapes,
00156 start,
00157 end
00158 );
00159
00160 if (debug & 2)
00161 {
00162 space(Pout, 2*level);
00163 Pout<< "findLeafLineOctant : bb:" << this->bb()
00164 << " returning from sub treeNode"
00165 << " with start:" << start << " subLeaf:"
00166 << long(subLeafPtr) << endl;
00167 }
00168
00169 return subLeafPtr;
00170 }
00171 else
00172 {
00173 FatalErrorIn(functionName)
00174 << "Sub node " << subNodePtr->bb()
00175 << " at octant " << octant
00176 << " does not contain start " << start
00177 << abort(FatalError);
00178 }
00179 }
00180 else
00181 {
00182
00183 const treeLeaf<Type>* subLeafPtr = getLeafPtr(octant);
00184
00185 if (subLeafPtr->bb().contains(direction, start))
00186 {
00187
00188 point tmp;
00189 if
00190 (
00191 !subLeafPtr->bb().intersects
00192 (
00193 end,
00194 start,
00195 tmp
00196 )
00197 )
00198 {
00199 FatalErrorIn(functionName)
00200 << "Sub leaf contains start " << start
00201 << " but line does not intersect its bb "
00202 << subLeafPtr->bb()
00203 << abort(FatalError);
00204 }
00205 start = tmp;
00206
00207 if (debug & 2)
00208 {
00209 space(Pout, 2*level);
00210 Pout<< "findLeafLineOctant : returning from intersecting"
00211 << " treeLeaf " << subLeafPtr->bb()
00212 << " with start:" << start << " subLeaf:"
00213 << long(subLeafPtr) << endl;
00214 }
00215
00216 return subLeafPtr;
00217 }
00218 else
00219 {
00220 FatalErrorIn(functionName)
00221 << "Sub leaf " << subLeafPtr->bb()
00222 << " at octant " << octant
00223 << " does not contain start " << start
00224 << abort(FatalError);
00225 }
00226 }
00227 }
00228 else
00229 {
00230
00231 const treeBoundBox emptyBb = this->bb().subBbox(midpoint(), octant);
00232
00233 if (emptyBb.contains(direction, start))
00234 {
00235 if (debug & 2)
00236 {
00237 space(Pout, 2*level);
00238 Pout<< "findLeafLineOctant : Empty node. Octant:" << octant
00239 << " start:" << start
00240 << " bb:" << this->bb()
00241 << " emptyBb:" << emptyBb << endl;
00242 }
00243
00244
00245 point tmp;
00246 if
00247 (
00248 !emptyBb.intersects
00249 (
00250 end,
00251 start,
00252 tmp
00253 )
00254 )
00255 {
00256 FatalErrorIn(functionName)
00257 << "Empty node contains start " << start
00258 << " but line does not intersect its (calculated)"
00259 << " bb " << emptyBb
00260 << endl << "This might be due to truncation error"
00261 << abort(FatalError);
00262 }
00263 start = tmp;
00264
00265 if (debug & 2)
00266 {
00267 space(Pout, 2*level);
00268 Pout<< "findLeafLineOctant : returning from intersecting with"
00269 << " empty " << emptyBb
00270 << " with start:" << start << " subLeaf:" << 0 << endl;
00271 }
00272
00273 return NULL;
00274 }
00275 else
00276 {
00277 FatalErrorIn(functionName)
00278 << "Empty node " << emptyBb
00279 << " at octant " << octant
00280 << " does not contain start " << start
00281 << abort(FatalError);
00282 }
00283 }
00284
00285 FatalErrorIn(functionName)
00286 << "Octant " << octant << " of cube " << this->bb()
00287 << " does not contain start " << start
00288 << abort(FatalError);
00289
00290 return NULL;
00291 }
00292
00293
00294
00295
00296
00297 template <class Type>
00298 Foam::treeNode<Type>::treeNode(const treeBoundBox& bb)
00299 :
00300 treeElem<Type>(bb),
00301 treeNodeName(),
00302 mid_(bb.midpoint()),
00303 subNodeTypes_(0),
00304 volType_(0)
00305 {
00306 for (label octantI=0; octantI<8; octantI++)
00307 {
00308 subNodes_[octantI] = NULL;
00309 setVolType(octantI, octree<Type>::UNKNOWN);
00310 }
00311 }
00312
00313
00314
00315 template <class Type>
00316 Foam::treeNode<Type>::treeNode(Istream& is)
00317 {
00318 is >> *this;
00319 }
00320
00321
00322
00323
00324 template <class Type>
00325 Foam::treeNode<Type>::~treeNode()
00326 {
00327 for (int octant=0; octant<8; octant++)
00328 {
00329 if (subNodes()[octant])
00330 {
00331 if (isNode(octant))
00332 {
00333 delete getNodePtr(octant);
00334 }
00335 else
00336 {
00337 delete getLeafPtr(octant);
00338 }
00339 }
00340 }
00341 }
00342
00343
00344
00345
00346
00347 template <class Type>
00348 void Foam::treeNode<Type>::distribute
00349 (
00350 const label level,
00351 octree<Type>& top,
00352 const Type& shapes,
00353 const labelList& indices
00354 )
00355 {
00356 if (debug & 1)
00357 {
00358 space(Pout, level);
00359 Pout<< "treeNode::distributing " << indices.size() << endl;
00360 }
00361
00362
00363 for (label octant=0; octant<8; octant++)
00364 {
00365 if (subNodes()[octant])
00366 {
00367 printNode(Pout, level);
00368 FatalErrorIn
00369 (
00370 "treeNode<Type>::distribute(const label, octree<Type>&, "
00371 "const Type&, const labelList&)"
00372 ) << "subNode already available at octant:" << octant
00373 << abort(FatalError);
00374 }
00375 else
00376 {
00377 treeLeaf<Type>* subLeafPtr = new treeLeaf<Type>
00378 (
00379 this->bb().subBbox(midpoint(), octant),
00380 indices.size()
00381 );
00382
00383 top.setLeaves(top.nLeaves() + 1);
00384 setLeafPtr(octant, subLeafPtr);
00385 }
00386 }
00387
00388
00389
00390 forAll(indices, i)
00391 {
00392 const label shapei = indices[i];
00393
00394 for (label octant=0; octant<8; octant++)
00395 {
00396 treeLeaf<Type>* leafPtr = getLeafPtr(octant);
00397
00398 if (shapes.overlaps(shapei, leafPtr->bb()))
00399 {
00400 if (debug == 1)
00401 {
00402 space(Pout, level);
00403 Pout<< "inserting " << shapei;
00404 shapes.write(Pout, shapei);
00405 Pout<< " into " << leafPtr->bb() << endl;
00406 }
00407 leafPtr->insert(shapei);
00408 top.setEntries(top.nEntries() + 1);
00409 }
00410 }
00411 }
00412
00413
00414 for (label octant=0; octant<8; octant++)
00415 {
00416 treeLeaf<Type>* subLeafPtr = getLeafPtr(octant);
00417
00418 if (subLeafPtr->size() == 0)
00419 {
00420
00421 setLeafPtr(octant, NULL);
00422 delete subLeafPtr;
00423 top.setLeaves(top.nLeaves() - 1);
00424 }
00425 else
00426 {
00427
00428 subLeafPtr->trim();
00429 }
00430 }
00431
00432 if (debug & 1)
00433 {
00434 space(Pout, level);
00435 Pout<< "end of treeNode::distribute" << endl;
00436 }
00437 }
00438
00439
00440
00441 template <class Type>
00442 void Foam::treeNode<Type>::redistribute
00443 (
00444 const label level,
00445 octree<Type>& top,
00446 const Type& shapes,
00447 const label refineLevel
00448 )
00449 {
00450 if (debug & 1)
00451 {
00452 space(Pout, level);
00453 Pout<< "treeNode::redistribute with level:" << level
00454 << " refineLevel:" << refineLevel << endl;
00455 }
00456
00457
00458 if (level < refineLevel)
00459 {
00460 for (label octant=0; octant<8; octant++)
00461 {
00462 if (subNodes()[octant])
00463 {
00464 if (isNode(octant))
00465 {
00466 getNodePtr(octant)->redistribute
00467 (
00468 level+1,
00469 top,
00470 shapes,
00471 refineLevel
00472 );
00473 }
00474 }
00475 }
00476 }
00477 else
00478 {
00479
00480 if (debug & 1)
00481 {
00482 space(Pout, level);
00483 Pout<< "treeNode::redistribute : now at correct level" << endl;
00484 }
00485
00486
00487 for (label octant=0; octant<8; octant++)
00488 {
00489 if (subNodes()[octant])
00490 {
00491 if (isNode(octant))
00492 {
00493 FatalErrorIn
00494 (
00495 "treeNode<Type>::redistribute(const int, octree& top,"
00496 "const int, const treeBoundBox&)"
00497 ) << "found treeNode instead of treeLeaf" << endl
00498 << abort(FatalError);
00499 }
00500 else
00501 {
00502 treeLeaf<Type>* leafPtr = getLeafPtr(octant);
00503
00504 treeLeaf<Type>* newSubPtr = leafPtr->redistribute
00505 (
00506 level,
00507 top,
00508 shapes
00509 );
00510
00511 if (newSubPtr && (newSubPtr != leafPtr))
00512 {
00513
00514
00515 if (debug & 1)
00516 {
00517 Pout<< "deleting "
00518 << top.nEntries() - leafPtr->size()
00519 << " entries" << endl;
00520 }
00521 top.setEntries(top.nEntries() - leafPtr->size());
00522
00523 delete leafPtr;
00524
00525 top.setLeaves(top.nLeaves() - 1);
00526
00527 setNodePtr(octant, newSubPtr);
00528 }
00529 }
00530 }
00531 }
00532 if (debug & 1)
00533 {
00534 space(Pout, level);
00535 Pout<< "end of treeNode::redistribute for correct level" << endl;
00536 }
00537 }
00538
00539 if (debug & 1)
00540 {
00541 space(Pout, level);
00542 Pout<< "return from treeNode::redistribute with bb:" << this->bb()
00543 << endl;
00544 }
00545 }
00546
00547
00548
00549 template <class Type>
00550 Foam::label Foam::treeNode<Type>::setSubNodeType
00551 (
00552 const label level,
00553 octree<Type>& top,
00554 const Type& shapes
00555 )
00556 {
00557 if (debug & 4)
00558 {
00559 space(Pout, level);
00560 Pout<< "treeNode::setSubNodeType with level:" << level
00561 << " bb:" << this->bb() << endl;
00562 }
00563
00564 label myType = -1;
00565
00566 for (label octant=0; octant<8; octant++)
00567 {
00568 label subType = -1;
00569
00570 if (subNodes()[octant])
00571 {
00572 if (isNode(octant))
00573 {
00574 subType = getNodePtr(octant)->setSubNodeType
00575 (
00576 level+1,
00577 top,
00578 shapes
00579 );
00580 }
00581 else
00582 {
00583 subType = getLeafPtr(octant)->setSubNodeType
00584 (
00585 level+1,
00586 top,
00587 shapes
00588 );
00589 }
00590 }
00591 else
00592 {
00593
00594
00595 const treeBoundBox subBb = this->bb().subBbox(midpoint(), octant);
00596
00597 subType = shapes.getSampleType(top, subBb.midpoint());
00598 }
00599
00600 if (debug & 4)
00601 {
00602 space(Pout, level);
00603 Pout<< "treeNode::setSubNodeType : setting octant with bb:"
00604 << this->bb().subBbox(midpoint(), octant)
00605 << " to type:" << octree<Type>::volType(subType) << endl;
00606 }
00607 setVolType(octant, subType);
00608
00609
00610
00611 if (myType == -1)
00612 {
00613 myType = subType;
00614 }
00615 else if (subType != myType)
00616 {
00617 myType = octree<Type>::MIXED;
00618 }
00619 }
00620
00621 if (debug & 4)
00622 {
00623 space(Pout, level);
00624 Pout<< "return from treeNode::setSubNodeType with type:"
00625 << octree<Type>::volType(myType)
00626 << " bb:" << this->bb() << endl;
00627 }
00628
00629 return myType;
00630 }
00631
00632
00633
00634 template <class Type>
00635 Foam::label Foam::treeNode<Type>::getSampleType
00636 (
00637 const label level,
00638 const octree<Type>& top,
00639 const Type& shapes,
00640 const point& sample
00641 ) const
00642 {
00643 if (debug & 4)
00644 {
00645 space(Pout, level);
00646 Pout<< "treeNode::getSampleType with level:" << level
00647 << " bb:" << this->bb() << " sample:" << sample << endl;
00648 }
00649
00650
00651 bool onEdge = false;
00652
00653 label octant = this->bb().subOctant(midpoint(), sample, onEdge);
00654
00655 label type = getVolType(octant);
00656
00657 if (type == octree<Type>::MIXED)
00658 {
00659
00660 if (subNodes()[octant])
00661 {
00662 if (isNode(octant))
00663 {
00664
00665 type = getNodePtr(octant)->getSampleType
00666 (
00667 level + 1,
00668 top,
00669 shapes,
00670 sample
00671 );
00672 }
00673 else
00674 {
00675
00676 type = getLeafPtr(octant)->getSampleType
00677 (
00678 level + 1,
00679 top,
00680 shapes,
00681 sample
00682 );
00683 }
00684 }
00685 else
00686 {
00687
00688 FatalErrorIn
00689 (
00690 "treeNode<Type>::getSampleType"
00691 "(const label, octree<Type>&, const Type&, const point&)"
00692 ) << "Empty node bb:" << this->bb().subBbox(midpoint(), octant)
00693 << " has non-mixed type:"
00694 << octree<Type>::volType(type)
00695 << abort(FatalError);
00696 }
00697 }
00698
00699 if (type == octree<Type>::MIXED)
00700 {
00701 FatalErrorIn
00702 (
00703 "treeNode<Type>::getSampleType"
00704 "(const label, octree<Type>&, const Type&, const point&)"
00705 ) << "Type is MIXED when searching for " << sample
00706 << " at level " << this->bb() << endl
00707 << "This probably is because the octree has not been constructed"
00708 << " with search facility." << exit(FatalError);
00709 }
00710
00711 if (debug & 4)
00712 {
00713 space(Pout, level);
00714 Pout<< "return from treeNode::getSampleType with type:"
00715 << octree<Type>::volType(type)
00716 << " bb:" << this->bb()
00717 << " sample:" << sample << endl;
00718 }
00719 return type;
00720 }
00721
00722
00723 template <class Type>
00724 Foam::label Foam::treeNode<Type>::find
00725 (
00726 const Type& shapes,
00727 const point& sample
00728 ) const
00729 {
00730
00731
00732 bool onEdge = false;
00733
00734 label octant = this->bb().subOctant(midpoint(), sample, onEdge);
00735
00736 if (subNodes()[octant])
00737 {
00738 if (isNode(octant))
00739 {
00740
00741 return getNodePtr(octant)->find(shapes, sample);
00742 }
00743 else
00744 {
00745
00746 return getLeafPtr(octant)->find(shapes, sample);
00747 }
00748 }
00749 return -1;
00750 }
00751
00752
00753 template <class Type>
00754 bool Foam::treeNode<Type>::findTightest
00755 (
00756 const Type& shapes,
00757 const point& sample,
00758 treeBoundBox& tightest
00759 ) const
00760 {
00761 bool changed = false;
00762 bool onEdge = false;
00763
00764 label sampleOctant = this->bb().subOctant(midpoint(), sample, onEdge);
00765
00766
00767
00768
00769 for (label octantI=0; octantI<8; octantI++)
00770 {
00771 label octant;
00772 if (octantI == 0)
00773 {
00774
00775 octant = sampleOctant;
00776 }
00777 else if (octantI == sampleOctant)
00778 {
00779 octant = 0;
00780 }
00781 else
00782 {
00783 octant = octantI;
00784 }
00785
00786 if (subNodes()[octant])
00787 {
00788 if (isNode(octant))
00789 {
00790
00791 const treeNode<Type>* subNodePtr = getNodePtr(octant);
00792
00793 if (subNodePtr->bb().overlaps(tightest))
00794 {
00795
00796 changed |= subNodePtr->findTightest
00797 (
00798 shapes,
00799 sample,
00800 tightest
00801 );
00802 }
00803 }
00804 else
00805 {
00806
00807 const treeLeaf<Type>* subLeafPtr = getLeafPtr(octant);
00808
00809 if (subLeafPtr->bb().overlaps(tightest))
00810 {
00811
00812 changed |= subLeafPtr->findTightest
00813 (
00814 shapes,
00815 sample,
00816 tightest
00817 );
00818 }
00819 }
00820 }
00821 }
00822
00823 return changed;
00824 }
00825
00826
00827 template <class Type>
00828 bool Foam::treeNode<Type>::findNearest
00829 (
00830 const Type& shapes,
00831 const point& sample,
00832 treeBoundBox& tightest,
00833 label& tightestI,
00834 scalar& tightestDist
00835 ) const
00836 {
00837 if (debug & 8)
00838 {
00839 Pout<< "In findNearest with sample:" << sample << " cube:"
00840 << this->bb() << " tightest:" << tightest << endl;
00841 }
00842
00843 bool changed = false;
00844 bool onEdge = false;
00845
00846 label sampleOctant = this->bb().subOctant(midpoint(), sample, onEdge);
00847
00848
00849
00850
00851 for (label octantI=0; octantI<8; octantI++)
00852 {
00853 label octant;
00854 if (octantI == 0)
00855 {
00856
00857 octant = sampleOctant;
00858 }
00859 else if (octantI == sampleOctant)
00860 {
00861 octant = 0;
00862 }
00863 else
00864 {
00865 octant = octantI;
00866 }
00867
00868 if (subNodes()[octant])
00869 {
00870 if (isNode(octant))
00871 {
00872
00873 const treeNode<Type>* subNodePtr = getNodePtr(octant);
00874
00875 if (subNodePtr->bb().overlaps(tightest))
00876 {
00877
00878 changed |= subNodePtr->findNearest
00879 (
00880 shapes,
00881 sample,
00882 tightest,
00883 tightestI,
00884 tightestDist
00885 );
00886 }
00887 }
00888 else
00889 {
00890
00891 const treeLeaf<Type>* subLeafPtr = getLeafPtr(octant);
00892
00893 if (subLeafPtr->bb().overlaps(tightest))
00894 {
00895
00896 changed |= subLeafPtr->findNearest
00897 (
00898 shapes,
00899 sample,
00900 tightest,
00901 tightestI,
00902 tightestDist
00903 );
00904 }
00905 }
00906 }
00907 }
00908
00909 if (debug & 8)
00910 {
00911 Pout<< "Exiting findNearest for sample:" << sample << " cube:"
00912 << this->bb() << " tightestI:" << tightestI << endl;
00913 }
00914
00915 return changed;
00916 }
00917
00918
00919 template <class Type>
00920 bool Foam::treeNode<Type>::findNearest
00921 (
00922 const Type& shapes,
00923 const linePointRef& ln,
00924 treeBoundBox& tightest,
00925 label& tightestI,
00926 point& linePoint,
00927 point& shapePoint
00928 ) const
00929 {
00930 bool changed = false;
00931 bool onEdge = false;
00932
00933 label sampleOctant = this->bb().subOctant(midpoint(), ln.centre(), onEdge);
00934
00935
00936
00937
00938 for (label octantI=0; octantI<8; octantI++)
00939 {
00940 label octant;
00941 if (octantI == 0)
00942 {
00943
00944 octant = sampleOctant;
00945 }
00946 else if (octantI == sampleOctant)
00947 {
00948 octant = 0;
00949 }
00950 else
00951 {
00952 octant = octantI;
00953 }
00954
00955 if (subNodes()[octant])
00956 {
00957 if (isNode(octant))
00958 {
00959
00960 const treeNode<Type>* subNodePtr = getNodePtr(octant);
00961
00962 if (subNodePtr->bb().overlaps(tightest))
00963 {
00964
00965 changed |= subNodePtr->findNearest
00966 (
00967 shapes,
00968 ln,
00969 tightest,
00970 tightestI,
00971 linePoint,
00972 shapePoint
00973 );
00974 }
00975 }
00976 else
00977 {
00978
00979 const treeLeaf<Type>* subLeafPtr = getLeafPtr(octant);
00980
00981 if (subLeafPtr->bb().overlaps(tightest))
00982 {
00983
00984 changed |= subLeafPtr->findNearest
00985 (
00986 shapes,
00987 ln,
00988 tightest,
00989 tightestI,
00990 linePoint,
00991 shapePoint
00992 );
00993 }
00994 }
00995 }
00996 }
00997
00998 return changed;
00999 }
01000
01001
01002 template <class Type>
01003 bool Foam::treeNode<Type>::findBox
01004 (
01005 const Type& shapes,
01006 const boundBox& box,
01007 labelHashSet& elements
01008 ) const
01009 {
01010 bool changed = false;
01011 bool onEdge = false;
01012
01013 label sampleOctant = this->bb().subOctant
01014 (
01015 midpoint(),
01016 box.midpoint(),
01017 onEdge
01018 );
01019
01020
01021
01022
01023 for (label octantI=0; octantI<8; octantI++)
01024 {
01025 label octant;
01026 if (octantI == 0)
01027 {
01028
01029 octant = sampleOctant;
01030 }
01031 else if (octantI == sampleOctant)
01032 {
01033 octant = 0;
01034 }
01035 else
01036 {
01037 octant = octantI;
01038 }
01039
01040 if (subNodes()[octant])
01041 {
01042 if (isNode(octant))
01043 {
01044
01045 const treeNode<Type>* subNodePtr = getNodePtr(octant);
01046
01047 if (subNodePtr->bb().overlaps(box))
01048 {
01049
01050 changed |= subNodePtr->findBox(shapes, box, elements);
01051 }
01052 }
01053 else
01054 {
01055
01056 const treeLeaf<Type>* subLeafPtr = getLeafPtr(octant);
01057
01058 if (subLeafPtr->bb().overlaps(box))
01059 {
01060
01061 changed |= subLeafPtr->findBox(shapes, box, elements);
01062 }
01063 }
01064 }
01065 }
01066
01067 return changed;
01068 }
01069
01070
01071
01072 template <class Type>
01073 const Foam::treeLeaf<Type>* Foam::treeNode<Type>::findLeafLine
01074 (
01075 const int level,
01076 const Type& shapes,
01077 point& start,
01078 const point& end
01079 ) const
01080 {
01081 if (debug & 2)
01082 {
01083 space(Pout, 2*level);
01084 Pout<< "findLeafLine : bb:" << this->bb() << " mid:" << midpoint()
01085 << " start:" << start << endl;
01086 }
01087
01088 scalar typDim = this->bb().avgDim();
01089
01090 const vector direction = end - start;
01091
01092
01093
01094
01095
01096 label iter = 0;
01097
01098 while (true)
01099 {
01100 if (!this->bb().contains(direction, start))
01101 {
01102 if (debug & 2)
01103 {
01104 space(Pout, 2*level);
01105 Pout<< "findLeafLine : Start not inside bb " << this->bb()
01106 << ". Returning with start:" << start << " subLeaf:"
01107 << 0 << endl;
01108 }
01109 return NULL;
01110 }
01111
01112
01113 if ((mag(start - end)/typDim) < SMALL)
01114 {
01115 if (debug & 2)
01116 {
01117 space(Pout, 2*level);
01118 Pout<< "findLeafLine : start equals end"
01119 << ". Returning with start:" << start << " subLeaf:"
01120 << 0 << endl;
01121 }
01122 return NULL;
01123 }
01124
01125 if (iter >= 4)
01126 {
01127
01128 break;
01129 }
01130
01131 bool onEdge = false;
01132 label octant = this->bb().subOctant
01133 (
01134 midpoint(), direction, start, onEdge
01135 );
01136
01137
01138 const treeLeaf<Type>* leafPtr = findLeafLineOctant
01139 (
01140 level,
01141 shapes,
01142 octant,
01143 direction,
01144 start,
01145 end
01146 );
01147
01148 if (leafPtr)
01149 {
01150
01151 if (debug & 2)
01152 {
01153 space(Pout, 2*level);
01154 Pout<< "findLeafLine : Found treeLeaf"
01155 << ". Returning with start:" << start << " subLeaf:"
01156 << long(leafPtr) << endl;
01157 }
01158
01159 return leafPtr;
01160 }
01161
01162 iter++;
01163 }
01164
01165
01166 FatalErrorIn
01167 (
01168 "treeNode<Type>::findLeafLine"
01169 "(const label, octree<Type>&, point&,"
01170 " const point&)"
01171 ) << "Did not leave bb " << this->bb()
01172 << " after " << iter
01173 << " iterations of updating starting point."
01174 << "start:" << start << " end:" << end
01175 << abort(FatalError);
01176
01177 return NULL;
01178 }
01179
01180
01181 template <class Type>
01182 void Foam::treeNode<Type>::findLeaves
01183 (
01184 List<treeLeaf<Type>*>& leafArray,
01185 label& leafIndex
01186 ) const
01187 {
01188
01189 for (label octant=0; octant<8; octant++)
01190 {
01191 if (subNodes()[octant])
01192 {
01193 if (isNode(octant))
01194 {
01195
01196 const treeNode<Type>* subNodePtr = getNodePtr(octant);
01197 subNodePtr->findLeaves(leafArray, leafIndex);
01198 }
01199 else
01200 {
01201
01202 treeLeaf<Type>* subLeafPtr = getLeafPtr(octant);
01203 leafArray[leafIndex++] = subLeafPtr;
01204 }
01205 }
01206 }
01207 }
01208
01209
01210 template <class Type>
01211 void Foam::treeNode<Type>::findLeaves
01212 (
01213 List<const treeLeaf<Type>*>& leafArray,
01214 label& leafIndex
01215 ) const
01216 {
01217
01218 for (label octant=0; octant<8; octant++)
01219 {
01220 if (subNodes()[octant])
01221 {
01222 if (isNode(octant))
01223 {
01224
01225 const treeNode<Type>* subNodePtr = getNodePtr(octant);
01226 subNodePtr->findLeaves(leafArray, leafIndex);
01227 }
01228 else
01229 {
01230
01231 treeLeaf<Type>* subLeafPtr = getLeafPtr(octant);
01232 leafArray[leafIndex++] = subLeafPtr;
01233 }
01234 }
01235 }
01236 }
01237
01238
01239 template <class Type>
01240 void Foam::treeNode<Type>::printNode
01241 (
01242 Ostream& os,
01243 const label level
01244 ) const
01245 {
01246 space(os, 2*level);
01247
01248 os << "node:" << this->bb() << endl;
01249
01250 for (label octant=0; octant<8; octant++)
01251 {
01252 label type = getVolType(octant);
01253
01254 string typeString = octree<Type>::volType(type);
01255
01256 if (!subNodes_[octant])
01257 {
01258 space(os, level);
01259 os << octant << ":" << typeString << " : null" << endl;
01260 }
01261 else if (isNode(octant))
01262 {
01263 space(os, level);
01264 os << octant << ":" << typeString << " : node" << endl;
01265 getNodePtr(octant)->printNode(os, level+1);
01266 }
01267 else
01268 {
01269 space(os, level);
01270 os << octant << ":" << typeString << " : leaf" << endl;
01271
01272 treeLeaf<Type>* leafPtr = getLeafPtr(octant);
01273 leafPtr->printLeaf(os, level+1);
01274 }
01275 }
01276 }
01277
01278
01279 template <class Type>
01280 void Foam::treeNode<Type>::writeOBJ
01281 (
01282 Ostream& os,
01283 const label level,
01284 label& vertNo
01285 ) const
01286 {
01287 point midPoint(this->bb().midpoint());
01288
01289 label midVertNo = vertNo;
01290 os << "v " << midPoint.x() << " " << midPoint.y() << " "
01291 << midPoint.z() << endl;
01292 vertNo++;
01293
01294 for (label octant=0; octant<8; octant++)
01295 {
01296 if (subNodes_[octant])
01297 {
01298 if (isNode(octant))
01299 {
01300 treeNode<Type>* nodePtr = getNodePtr(octant);
01301
01302 point subMidPoint(nodePtr->bb().midpoint());
01303 os << "v " << subMidPoint.x() << " " << subMidPoint.y() << " "
01304 << subMidPoint.z() << endl;
01305 os << "l " << midVertNo + 1<< " " << vertNo + 1 << endl;
01306 vertNo++;
01307
01308 nodePtr->writeOBJ(os, level+1, vertNo);
01309 }
01310 else
01311 {
01312 treeLeaf<Type>* leafPtr = getLeafPtr(octant);
01313
01314 point subMidPoint(leafPtr->bb().midpoint());
01315 os << "v " << subMidPoint.x() << " " << subMidPoint.y() << " "
01316 << subMidPoint.z() << endl;
01317 os << "l " << midVertNo + 1<< " " << vertNo + 1 << endl;
01318 vertNo++;
01319
01320
01321 }
01322 }
01323 }
01324 }
01325
01326
01327
01328
01329 template <class Type>
01330 Foam::Istream& Foam::operator>>(Istream& is, treeNode<Type>& oc)
01331 {
01332 for (label octant = 0; octant < 8; octant++)
01333 {
01334 oc.subNodes_[octant] = NULL;
01335 }
01336
01337 is >> oc.bb();
01338
01339 label nPtrs;
01340
01341
01342 is >> nPtrs;
01343
01344 is.readBegin("treeNode");
01345 for (label octant = 0; octant < nPtrs; octant++)
01346 {
01347 label index;
01348 is >> index;
01349
01350 if (index >= treeNode<Type>::leafOffset)
01351 {
01352
01353 treeLeaf<Type>* leafPtr = new treeLeaf<Type>(is);
01354 oc.setLeafPtr(index - treeNode<Type>::leafOffset, leafPtr);
01355 }
01356 else
01357 {
01358 oc.setNodePtr(index, new treeNode<Type>(is));
01359 }
01360 }
01361
01362
01363 is.readEnd("treeNode");
01364
01365
01366 is.check("Istream& operator>>(Istream&, treeNode&)");
01367
01368 return is;
01369 }
01370
01371
01372 template <class Type>
01373 Foam::Ostream& Foam::operator<<(Ostream& os, const treeNode<Type>& tn)
01374 {
01375
01376
01377
01378 label nPtrs = 0;
01379 for (label octant = 0; octant < 8; octant++)
01380 {
01381 if (tn.subNodes_[octant])
01382 {
01383 if (tn.isNode(octant) || tn.getLeafPtr(octant)->indices().size())
01384 {
01385 nPtrs++;
01386 }
01387 }
01388 }
01389
01390
01391
01392 os << token::SPACE << tn.bb() << token::SPACE << nPtrs
01393 << token::SPACE << token::BEGIN_LIST;
01394
01395 for (label octant = 0; octant < 8; octant++)
01396 {
01397 if (tn.subNodes_[octant])
01398 {
01399 if (tn.isNode(octant))
01400 {
01401 const treeNode<Type>* subNodePtr = tn.getNodePtr(octant);
01402
01403
01404 os << token::SPACE << octant << token::SPACE << *subNodePtr
01405 << token::NL;
01406 }
01407 else if (tn.getLeafPtr(octant)->indices().size())
01408 {
01409
01410 const treeLeaf<Type>* subLeafPtr = tn.getLeafPtr(octant);
01411
01412 os << token::SPACE << octant + treeNode<Type>::leafOffset
01413 << token::SPACE << *subLeafPtr
01414 << token::NL;
01415 }
01416 }
01417 }
01418
01419 os << token::SPACE << token::END_LIST;
01420
01421 return os;
01422 }
01423
01424
01425