00001 #ifndef LOFT_H
00002 #define LOFT_H
00003 
00004 #include <vector>
00005 #include <list>
00006 #include <stdlib.h>
00007 #include <stdio.h>
00008 #include <algorithm>
00009 
00010 #include "Surface.h"
00011 #include "SplineCurve.h"
00012 
00013 bool mypredicate (int i, int j) {
00014   return (i>j);
00015 }
00016 
00017 template <class DataType>
00018 class Loft : public Surface<DataType> {
00019   typedef Part<DataType,3> PartType;
00020   typedef Connection ConBase;
00021   typedef Segment<DataType,3> SegType;
00022   typedef Facet<DataType> FacetType;
00023   typedef SplineCurve<DataType,3> SplineCurveType;
00024   typedef Surface<DataType> SurfaceType;
00025 
00026   enum interpolation { linear, overhauser };
00027 
00028 public:
00029 
00030   Loft() {
00031     SurfaceType();
00032     SurfaceType::loft = 1;
00033   }
00034   ~Loft() {}
00035 
00036   virtual void read(const int type, std::string file) {
00037     
00038 
00039     ( mklog() << "\nREADING "<<file<<" type = "<<type<<" dim = "<<3<<std::endl).flush();
00040     switch(type) {
00041     case PartType::pointList: {
00042 
00043       break;
00044     }
00045     case PartType::VTK_curve: {
00046 
00047       break;
00048     }
00049     case PartType::Brep2D: {
00050 
00051       break;
00052     }
00053     case PartType::Brep3D: {
00054       int numPoints, numFacets, dim, rdim, idtmp0, idtmp1, idtmp2;
00055       DataType vtmp;
00056       ( mklog() << "Loft reading filetype = 3D Brep .iss\n").flush();
00057       int i=0;
00058       std::string tmp;
00059 
00060       std::ifstream infile;
00061       infile.open(file.c_str(), std::ifstream::in);
00062       if (!infile.is_open()) {
00063         ( mklog() << std::system("pwd") << "\n\n\n\t\t\tCANNOT OPEN "<<file<<"\n\n\n").flush();
00064         std::cerr << std::system("pwd") << "\n\n\n\t\t\tCANNOT OPEN "<<file<<"\n\n\n";
00065         exit (-1);
00066       }
00067       setName(file);
00068       infile >> dim; infile >> rdim;
00069       infile >> numPoints;
00070       while (infile.good() && i < numPoints) {
00071 
00072         if (!infile.good()) break;
00073 #ifdef LOFT_DEBUG_ON
00074         ( mklog() << " reading line "<<i<<"\t").flush();
00075 #endif
00076         PType p(0.);
00077 
00078         infile >> vtmp;         p(0) = vtmp;
00079         infile >> vtmp;         p(1) = vtmp;
00080         infile >> vtmp;         p(2) = vtmp;
00081 
00082 #ifdef LOFT_DEBUG_ON
00083         ( mklog() << " read " << p << std::endl ).flush();
00084 #endif
00085         AddPoint(p);
00086 #ifdef LOFT_DEBUG_ON
00087         ( mklog() << Name() <<" # "<< i <<"\t"<< GetNthPoint(i)(0)
00088                   <<"\t"<< GetNthPoint(i)(1) <<"\t"<< GetNthPoint(i)(2) << std::endl ).flush();
00089 #endif
00090         ( mklog() << Name() <<" # "<< i <<"\t"<< GetNthPoint(i)(0)
00091                   <<"\t"<< GetNthPoint(i)(1) <<"\t"<< GetNthPoint(i)(2) << std::endl ).flush();
00092         i++;
00093       }
00094       i = 0;
00095       infile >> numFacets;
00096       while (infile.good() && i < numFacets) {
00097 
00098         if (!infile.good()) break;
00099 
00100         FacetType c;
00101 
00102         infile >> idtmp0;               c.setNthCon(0,idtmp0);
00103         infile >> idtmp1;               c.setNthCon(1,idtmp1);
00104         infile >> idtmp2;               c.setNthCon(2,idtmp2);
00105         AddConnection(c);
00106 #ifdef LOFT_DEBUG_ON
00107         ( mklog()<< Name()<<" # "<<i<<"\t"<<GetNthFacet(i)->getNthCon(0)<<"\t"
00108                  <<GetNthFacet(i)->getNthCon(1)<<"\t"<<GetNthFacet(i)->getNthCon(2)<< std::endl
00109                 <<i<<"\t"<<SurfaceType::connectivity[i].getNthCon(0)<<"\t"<<SurfaceType::connectivity[i].getNthCon(1)<<"\t"
00110                <<SurfaceType::connectivity[i].getNthCon(2) ).flush();
00111 #endif
00112         i++;
00113       }
00114 
00115       if ( norm_2(GetNthPoint(0)-GetNthPoint(GetNumPoints()-1) ) < TOL ) {
00116         ( mklog() << "closed\n").flush();
00117         setClosure(1);
00118       }
00119       else
00120       {
00121         ( mklog() << "open\n").flush();
00122         setClosure(0);
00123       }
00124 
00125       infile.close();
00126 
00127       SurfaceType::measure();
00128 
00129       break;
00130     }
00131     case PartType::STL: {
00132 
00133       break;
00134     }
00135     case PartType::loft: {
00136       int itmp;
00137       DataType vtmp;
00138       std::string stmp;
00139 
00140       std::ifstream infile;
00141       infile.open(file.c_str(), std::ifstream::in);
00142       if (!infile.is_open()) {
00143         ( mklog() <<std::system("pwd")<< "\n\n\n\t\t\tLoft read() CANNOT OPEN " << file << "\n\n\n").flush();
00144         std::cerr <<std::system("pwd")<< "\n\n\n\t\t\tLoft read() CANNOT OPEN " << file << "\n\n\n";
00145         exit (-1);
00146       }
00147 
00148       setName(file);
00149       infile >> itmp;
00150       setNProfiles(itmp);
00151       infile >> itmp;
00152       setNChords(itmp);
00153       infile >> itmp;
00154       setNSpans(itmp);
00155       infile >> itmp; spline_intType = itmp;
00156       for (register int i = 0; i < getNProfiles(); i++) {
00157         PType position(0.);
00158         PType orientation(0.);
00159         PType scale(0.);
00160         PType reflection(0.);
00161         MType transtmp(4,4);
00162         transtmp.assign(imatrix(4));
00163         infile >> itmp;
00164         AddProfileType(itmp);
00165         infile >> stmp;
00166         AddProfileName(stmp);
00167         infile >> vtmp; position(0) = vtmp;
00168         infile >> vtmp; position(1) = vtmp;
00169         infile >> vtmp; position(2) = vtmp;
00170         AddPosition(position);
00171         infile >> vtmp; orientation(0) = (vtmp);
00172         infile >> vtmp; orientation(1) = (vtmp);
00173         infile >> vtmp; orientation(2) = (vtmp);
00174         AddOrientation(orientation);
00175         infile >> vtmp; reflection(0) = (vtmp);
00176         infile >> vtmp; reflection(1) = (vtmp);
00177         infile >> vtmp; reflection(2) = (vtmp);
00178         AddReflection(reflection);
00179         infile >> vtmp; scale(0) = vtmp;
00180         infile >> vtmp; scale(1) = vtmp;
00181         infile >> vtmp; scale(2) = vtmp;
00182         AddScale(scale);
00183         AddTransform(transtmp);
00184 
00185         infile >> itmp;  addIntType(itmp);
00186 
00187       }
00188 
00189       infile.close();
00190 
00191       
00192       SplineCurveType * profileSpline;
00193       MType rTmp(4,4);
00194       for (register int i = 0; i < getNProfiles(); i++) {
00195         profileSpline = new SplineCurveType;
00196         AddProfileSpline(profileSpline);
00197         GetNthProfileSpline(i)->read(profileTypes[i],profileNames[i] );
00198 
00199         if ( norm_2(reflections[i]) > TOL  ) {
00200           rTmp = reflMat3D(reflections[i]);
00201           ( mklog() << Name() << " profile " << i << " " << reflections[i] << " " << rTmp << std::endl ).flush();
00202           GetNthProfileSpline(i)->GetCPoints_org().clear();
00203           prodMVP( rTmp, GetNthProfileSpline(i)->GetNodes(), GetNthProfileSpline(i)->GetCPoints_org() );
00204 
00205           GetNthProfileSpline(i)->GetNodes().clear();
00206           std::reverse_copy( GetNthProfileSpline(i)->GetCPoints_org().begin(), GetNthProfileSpline(i)->GetCPoints_org().end(), std::back_inserter (GetNthProfileSpline(i)->GetNodes() ) );
00207 
00208           GetNthProfileSpline(i)->GetCPoints().clear();
00209 
00210         }
00211 
00212         if (intType[i] == overhauser) {
00213           GetNthProfileSpline(i)->interpolate(getNSpans());
00214         }
00215         else if (intType[i] == linear) {
00216           GetNthProfileSpline(i)->interpolate_linear(getNSpans());
00217         }
00218         GetNthProfileSpline(i)->measure();
00219       }
00220 #ifdef LOFT_DEBUG_ON
00221       ( mklog() << "Done reading and interpolating Prolifes\n").flush();
00222 #endif
00223       GenLoft();
00224       makeFacets(0);
00225 
00226       break;
00227     }
00228     }
00229   }
00230 
00231   void GenLoft() {
00232 #ifdef LOFT_DEBUG_ON
00233     ( mklog() << " Genloft() Begun @@@@@@@@@@@@\n" ).flush();
00234 #endif
00235     PType fp;
00236     int c=0;
00237 
00238     SurfaceType::cPoints.clear();
00239     SurfaceType::cPoints_org.clear();
00240     spanSplines.clear();
00241 
00242     PType zero(0.), rCenScal(0.);
00243     MType pTran(4,4), tTmp(4,4), rTmp(4,4), sTmp(4,4);
00244     pTran.assign(imatrix(4));
00245     tTmp.assign(imatrix(4));
00246     rTmp.assign(imatrix(4));
00247     sTmp.assign(imatrix(4));
00248 
00249     for (register int i = 0; i < getNProfiles(); i++) {
00250 
00251       pTran = tranMat3D(-GetNthProfileSpline(i)->Centroid());
00252       
00253       rTmp = rotMat3D(orientations[i]);
00254       pTran = boost::numeric::ublas::prod(rTmp, pTran);
00255       
00256       sTmp = scalMat3D(scales[i]);
00257       pTran = boost::numeric::ublas::prod(sTmp, pTran);
00258       
00259       rCenScal = prodMP(sTmp, GetNthProfileSpline(i)->Centroid());
00260       tTmp = tranMat3D(rCenScal);
00261       pTran = boost::numeric::ublas::prod(tTmp,  pTran);
00262       
00263       tTmp = tranMat3D(positions[i] - GetNthProfileSpline(i)->Centroid());
00264       pTran = boost::numeric::ublas::prod(tTmp,  pTran);
00265 
00266       SetNthTransform(i,pTran);
00267       prodMVP( GetNthTransform(i), GetNthProfileSpline(i)->GetCPoints_org(), GetNthProfileSpline(i)->GetCPoints() );
00268 
00269 #ifdef LOFT_DEBUG_ON
00270       ( mklog() << " " << Name() << "    positions="<< positions[i] << " scales="<< scales[i]
00271                    <<" orientations="<< orientations[i] << std::endl
00272                   << " MAT " << GetNthTransform(i) << std::endl ).flush();
00273 #endif
00274 
00275     }
00276 
00277 #ifdef LOFT_DEBUG_ON
00278     ( mklog() << " Genloft() Done @@@@@@@@@@@@\n" ).flush();
00279 #endif
00280 
00281     spanSplines.clear();
00282     SplineCurveType * spantmp;
00283 
00284     for (c=0; c < getNSpans(); c++) {
00285       spantmp = new SplineCurveType;
00286 
00287       spantmp->setClosure(1);
00288 
00289       for (register int i = 0; i < getNProfiles(); i++) { 
00290         PType p(0.);
00291 
00292         p = GetNthProfileSpline(i)->GetNthPoint(c);
00293 
00294 #ifdef LOFT_DEBUG_ON
00295         ( mklog() << GetNthProfileSpline(i)->Name() << " " << p << std::endl ).flush();
00296 #endif
00297         spantmp->AddSplineNode(p);
00298       }
00299 
00300       if ( spline_intType == overhauser ) {
00301         if (getNChords() > 1) spantmp->interpolate(getNChords());
00302       }
00303       else if ( spline_intType == linear ) {
00304         if (getNChords() > 1) spantmp->interpolate_linear(getNChords());
00305       }
00306       for (register int j = 0; j < spantmp->GetNumPoints(); j++ ) {
00307         AddPoint(spantmp->GetNthPoint(j));
00308 #ifdef LOFT_DEBUG_ON
00309         ( mklog() << Name() << " " << GetNthPoint(GetNumPoints()-1) << std::endl ).flush();
00310 #endif
00311       }
00312 
00313       AddSpanSpline(spantmp);
00314 
00315     }
00316 
00317     SurfaceType::cPoints_org.clear();
00318     std::copy (SurfaceType::cPoints.begin (), SurfaceType::cPoints.end (), std::back_inserter (SurfaceType::cPoints_org));
00319     SurfaceType::velocities_loc.resize(GetNumPoints());
00320 #ifdef LOFT_DEBUG_ON
00321     ( mklog() << "GenLoft DONE adding points "<<GetNumPoints()<<std::endl).flush();
00322 #endif
00323   }
00324 
00325   virtual void makeFacets( int flip) {
00326     int aid, bid, cid;
00327 
00328 #ifdef LOFT_DEBUG_ON
00329     int tmpfileflag = 0;
00330     FILE * fout;
00331     if (tmpfileflag == 1) fout = fopen ("tmp.stl","w");
00332     if (tmpfileflag == 1) fprintf(fout,"solid %s\n",Name().c_str());
00333 #endif
00334     PType a(0.); PType b(0.); PType c(0.);
00335 
00336     for (register int ii = 0; ii < (GetNumPoints()-getNChords()-1); ii++ ) {
00337 
00338       a = GetNthPoint(ii);
00339       aid = ii;
00340       b = GetNthPoint(ii+getNChords());
00341       bid = ii+getNChords();
00342       c = GetNthPoint(ii+getNChords()+1);
00343       cid = ii+getNChords()+1;
00344 
00345       if ( (ii+1) % getNChords() != 0 ) {
00346         if (flip==0) {
00347 #ifdef LOFT_DEBUG_ON
00348           if (tmpfileflag == 1) {
00349             fprintf(fout,"top ii=%d flip=%d : %d %d %d\n",flip,ii,aid,bid,cid);
00350             fprintf(fout,"vertex %g\t%g\t%g\n",a(0),a(1),a(2));
00351             fprintf(fout,"vertex %g\t%g\t%g\n",b(0),b(1),b(2));
00352             fprintf(fout,"vertex %g\t%g\t%g\n",c(0),c(1),c(2));
00353           }
00354 #endif
00355           FacetType facetemp(aid,bid,cid);
00356           AddConnection(facetemp);
00357 
00358         }
00359         else {
00360 #ifdef LOFT_DEBUG_ON
00361           if (tmpfileflag == 1) {
00362             fprintf(fout,"bottom ii=%d flip=%d : %d %d %d\n",flip,ii,aid,cid,bid);
00363             fprintf(fout,"vertex %g\t%g\t%g\n",a(0),a(1),a(2));
00364             fprintf(fout,"vertex %g\t%g\t%g\n",c(0),c(1),c(2));
00365             fprintf(fout,"vertex %g\t%g\t%g\n",b(0),b(1),b(2));
00366           }
00367 #endif
00368           FacetType facetemp(aid,cid,bid);
00369           AddConnection(facetemp);
00370         }
00371 #ifdef LOFT_DEBUG_ON
00372         if (tmpfileflag == 1) fprintf(fout,"endloop\nendfacet\n");
00373 #endif
00374 
00375         b = GetNthPoint(ii+1);
00376         bid = ii+1;
00377 
00378         if (flip==0) {
00379 #ifdef LOFT_DEBUG_ON
00380           if (tmpfileflag == 1) {
00381             fprintf(fout,"top ii=%d flip=%d : %d %d %d\n",flip,ii,aid,cid,bid);
00382             fprintf(fout,"vertex %g\t%g\t%g\n",a(0),a(1),a(2));
00383             fprintf(fout,"vertex %g\t%g\t%g\n",b(0),b(1),b(2));
00384             fprintf(fout,"vertex %g\t%g\t%g\n",c(0),c(1),c(2));
00385           }
00386 #endif
00387           FacetType facetemp(aid,cid,bid);
00388           AddConnection(facetemp);
00389 
00390         }
00391         else {
00392 #ifdef LOFT_DEBUG_ON
00393           if (tmpfileflag == 1) {
00394             fprintf(fout,"bottom ii=%d flip=%d : %d %d %d\n",flip,ii,aid,bid,cid);
00395             fprintf(fout,"vertex %g\t%g\t%g\n",a(0),a(1),a(2));
00396             fprintf(fout,"vertex %g\t%g\t%g\n",b(0),b(1),b(2));
00397             fprintf(fout,"vertex %g\t%g\t%g\n",c(0),c(1),c(2));
00398           }
00399 #endif
00400           FacetType facetemp(aid,bid,cid);
00401           AddConnection(facetemp);
00402         }
00403 #ifdef LOFT_DEBUG_ON
00404         if (tmpfileflag == 1) fprintf(fout,"endloop\nendfacet\n");
00405 #endif
00406       }
00407 
00408     }
00409 
00410     int an=0, bn=0, cn=0, tmp=0;
00411 #ifdef LOFT_DEBUG_ON
00412     if (tmpfileflag == 1) fprintf(fout,"CLOSING LOFT ends\n");
00413 #endif
00414     int cnc = getNSpans();
00415     int cnf = getNChords();
00416     int n = cnf*(cnc-1);
00417 
00418     for (int m=0; m < cnc-1; m++) {
00419       an = cnf*m;
00420       bn = n-cnf*(m);
00421       cn = cnf*(m+1);
00422       
00423       if (an != bn && an != cn && bn != cn) {
00424         a =GetNthPoint(an);
00425         b =GetNthPoint(bn);
00426         c =GetNthPoint(cn);
00427         if ( SurfaceType::validFacet(a,b,c) ) {
00428           FacetType facetemp(an,bn,cn);
00429           AddConnection(facetemp);
00430         }
00431       }
00432     }
00433 
00434     int nn = cnf-1;
00435 #ifdef LOFT_DEBUG_ON
00436     fprintf(fout,"\nClosing far end nn=%d\n",nn);
00437 #endif
00438     for (int m=0; m < cnc-1; m++) {
00439       an = nn+cnf*m;
00440       bn = nn+cnf*(m+1);
00441       cn = nn+n-cnf*m;
00442 
00443       if (an != bn && an != cn && bn != cn) {
00444         if (bn < an) {
00445           tmp = bn;
00446           bn = cn;
00447           cn = tmp;
00448         }
00449 
00450         a =GetNthPoint(an);
00451         b =GetNthPoint(bn);
00452         c =GetNthPoint(cn);
00453         if ( SurfaceType::validFacet(a,b,c) ) {
00454           FacetType facetemp(an,bn,cn);
00455           AddConnection(facetemp);
00456         }
00457       }
00458     }
00459 
00460 #ifdef LOFT_DEBUG_ON
00461     if (tmpfileflag == 1) {
00462       fprintf(fout,"endsolid %s\n",Name().c_str());
00463       fclose (fout);
00464     }
00465 #endif
00466 
00467     makeUnique();
00468     SurfaceType::measure();
00469 
00470   }
00471 
00472   virtual void makeUnique() {
00473 #ifdef LOFT_DEBUG_ON
00474     ( mklog() << "UNIQUE POINTS\n     GetNumPoints()="<<GetNumPoints()<<std::endl).flush();
00475 #endif
00476     std::list<PType> CPointList;
00477     std::list<PType>::iterator it2;
00478     std::vector<FacetType> conTemp;
00479     std::copy (SurfaceType::cPoints.begin (), SurfaceType::cPoints.end (), std::back_inserter (CPointList));
00480     std::copy (SurfaceType::connectivity.begin (), SurfaceType::connectivity.end (), std::back_inserter (conTemp));
00481     PType me(0.); PType other(0.);
00482     int rmcount = 0, j=0, kk=0, jrem;
00483     std::vector<int> removedPoints;
00484 
00485 
00486     for ( register int k = 0; k < (GetNumPoints()/2); k++) {
00487       me = GetNthPoint(k);
00488 
00489       kk = k+1;
00490       it2 = CPointList.begin();
00491       advance (it2,(kk));
00492       j=kk;
00493       while ( it2 != CPointList.end() ) {
00494         other = *it2;
00495         if (norm_2(me-other) < TOL ) {
00496           jrem = j + rmcount;
00497           for ( register int f = 0; f < GetNumConnections(); f++) {
00498 
00499             if (GetNthFacet(f)->getNthCon(0) == jrem || GetNthFacet(f)->getNthCon(1) == jrem || GetNthFacet(f)->getNthCon(2) == jrem) {
00500               if (GetNthFacet(f)->getNthCon(0) == jrem) {
00501                 conTemp[f].cPid[0] = k;
00502               }
00503               if (GetNthFacet(f)->getNthCon(1) == jrem) {
00504                 conTemp[f].cPid[1] = k;
00505               }
00506               if (GetNthFacet(f)->getNthCon(2) == jrem) {
00507                 conTemp[f].cPid[2] = k;
00508               }
00509 
00510             } 
00511           } 
00512 #ifdef LOFT_DEBUG_ON
00513          ( mklog() << "end search of Facets\n"
00514                    <<"          ______removing cPoint["<<j<<"]\n"
00515                   <<"                  about to remove it2 " << *it2 << std::endl ).flush();
00516 #endif
00517           removedPoints.push_back(j+removedPoints.size());
00518           it2 = CPointList.erase(it2);
00519           rmcount++;
00520 #ifdef LOFT_DEBUG_ON
00521           ( mklog() << "          ______removed cPoint["<<j<<"].                "<< rmcount
00522                     <<" points have been removed. CPointList.size()="<<CPointList.size()<<std::endl ).flush();
00523 #endif
00524           if ( j == int(CPointList.size()) ) break;
00525 
00526         } 
00527         ++it2;
00528         j++;
00529       } 
00530       if (k==int(CPointList.size())-1 ) {
00531         break;
00532       }
00533     }
00534 
00535 #ifdef LOFT_DEBUG_ON
00536     ( mklog() << "\nRemoved "<<rmcount<<" redundant points. CPointList.size()="<<CPointList.size()<<std::endl).flush();
00537     for (register int i = 0; i < removedPoints.size(); i++) {
00538       ( mklog() << "removedPoints["<<i<<"] = "<<removedPoints[i]<<std::endl).flush();
00539     }
00540     ( mklog() << "\n+*+*+*+*+*+*+ Search and update of facet cPIDs\n").flush();
00541 #endif
00542     std::vector<int>::iterator it;
00543     int needle2[1];
00544 
00545     for ( register int ff = 0; ff < GetNumConnections(); ff++) {
00546 
00547 #ifdef LOFT_DEBUG_ON
00548       ( mklog() << "          Remaping Facet(0) "<<ff<<" : {"<<GetNthFacet(ff)->getNthCon(0)<<","
00549                 <<GetNthFacet(ff)->getNthCon(1)<<","<<GetNthFacet(ff)->getNthCon(2)<<"} to ").flush();
00550 #endif
00551       for (register int m = 0; m < 3; m++) {
00552         needle2[0] = conTemp[ff].cPid[m];
00553         it = std::search (removedPoints.begin(), removedPoints.end(), needle2, needle2+1, mypredicate);
00554         conTemp[ff].cPid[m] -= (it-removedPoints.begin());
00555       }
00556 #ifdef LOFT_DEBUG_ON
00557       ( mklog() <<"{"<<GetNthFacet(ff)->getNthCon(0)<<","
00558                   <<GetNthFacet(ff)->getNthCon(1)<<","<<GetNthFacet(ff)->getNthCon(2)<<"}\n").flush();
00559 #endif
00560 
00561     }
00562 
00563 
00564     SurfaceType::cPoints.clear();
00565     std::copy (CPointList.begin (), CPointList.end (), std::back_inserter (SurfaceType::cPoints));
00566     SurfaceType::cPoints_org.clear();
00567     std::copy (CPointList.begin (), CPointList.end (), std::back_inserter (SurfaceType::cPoints_org));
00568     SurfaceType::connectivity.clear();
00569     std::copy (conTemp.begin (), conTemp.end (), std::back_inserter (SurfaceType::connectivity));
00570   }
00571 
00572   virtual DataType Length() const { return 0; }
00573   virtual DataType Area() const { return SurfaceType::area; }
00574   virtual DataType Volume() const { return SurfaceType::volume; }
00575   virtual void setLength(DataType val) {}
00576   virtual void setArea(DataType val) { SurfaceType::area = val; }
00577   virtual void setVolume(DataType val) { SurfaceType::volume = val; }
00578 
00579   virtual std::string Name() const { return SurfaceType::name; }
00580   virtual void setName(std::string val) { SurfaceType::name = val; }
00581 
00582   virtual DataType conArea() const { return SurfaceType::area; }
00583   void setConArea(DataType val) { SurfaceType::area = val; }
00584 
00585   virtual bool closure() const { return SurfaceType::closed; }
00586   virtual void setClosure(bool val) { SurfaceType::closed = val; }
00587 
00588   virtual PType Centroid() const { return SurfaceType::centroid; }
00589   virtual void setCentroid(PType val) { SurfaceType::centroid = val; }
00590 
00591   virtual PType operator [] (int i) const {
00592     return SurfaceType::cPoints[i];
00593   }
00594 
00595   virtual void AddPoint(const PType v) {
00596     SurfaceType::cPoints.push_back(v);
00597   }
00598 
00599   virtual int GetNumPoints() const {
00600     return ((int) SurfaceType::cPoints.size());
00601   }
00602 
00603   virtual PType GetNthPoint(int n) const {
00604     return (SurfaceType::cPoints[n]);
00605   }
00606 
00607   virtual PType* GetNthPointA(int n) {
00608     return (&SurfaceType::cPoints[n]);
00609   }
00610 
00611   virtual std::vector<PType>& GetCPoints() {
00612     return SurfaceType::cPoints;
00613   }
00614 
00615   virtual std::vector<PType>& GetCPoints_org() {
00616     return SurfaceType::cPoints_org;
00617   }
00618 
00619   virtual int GetNumCons() const {
00620     return ((int) SurfaceType::connectivity.size());
00621   }
00622 
00623   virtual int GetConnectivity() { return 3;}
00624   virtual void AddConnection(ConBase * v) {}
00625   virtual void AddConnection(SegType * v) {}
00626   virtual void AddConnection(FacetType v) {
00627     SurfaceType::connectivity.push_back(v);
00628   }
00629   virtual int GetNumConnections() const { return SurfaceType::connectivity.size(); }
00630   virtual ConBase * GetNthConnection(int n) const { return (ConBase*) &SurfaceType::connectivity[n]; }
00631   virtual SegType * GetNthSegment(int n) const { return (SegType *) 0; }
00632   virtual const FacetType * GetNthFacet(int n) const { return &SurfaceType::connectivity[n]; }
00633   virtual void setNthCon(int i, int val) {}
00634   virtual int getNthCon(int i) const { return -1; }
00635 
00636   virtual PType Normal() const { return (PType) 0.; }
00637   virtual void setNormal(PType val) {}
00638 
00639   
00640   virtual std::string getNthProfileName(int i) { return profileNames[i]; }
00641   virtual void AddProfileName(std::string val) { profileNames.push_back(val); }
00642   virtual int getNProfiles() { return NProfiles; }
00643   virtual void setNProfiles(int val) { NProfiles = val; }
00644   virtual int getNChords() { return NChords; }
00645   virtual void setNChords(int val) { NChords = val; }
00646   virtual int getNSpans() { return NSpans; }
00647   virtual void setNSpans(int val) { NSpans = val; }
00648 
00649 
00650   virtual void AddProfileType(int val) { profileTypes.push_back(val); }
00651   virtual int getNthProfileType(int i) { return profileTypes[i]; }
00652   virtual void AddPosition(PType val) { positions.push_back(val); }
00653   virtual void AddOrientation(PType val) { orientations.push_back(val); }
00654   virtual void AddReflection(PType val) { reflections.push_back(val); }
00655   virtual void AddScale(PType val) { scales.push_back(val); }
00656 
00657   virtual void AddProfileSpline(SplineCurveType * val) { profileSplines.push_back(val); }
00658   virtual SplineCurveType * GetNthProfileSpline(int i) { return profileSplines[i]; }
00659 
00660   virtual void AddSpanSpline(SplineCurveType * val) { spanSplines.push_back(val); }
00661   virtual SplineCurveType * GetNthSpanSpline(int i) { return spanSplines[i]; }
00662 
00663   virtual void AddTransform(MType val) { transformList.push_back(val); }
00664   virtual MType GetNthTransform(int i) { return transformList[i]; }
00665   virtual void SetNthTransform(int i, MType val) { transformList[i] = val; }
00666 
00667   virtual void updatePart(DataType dt, DataType time) {
00668 #ifdef LOFT_DEBUG_ON
00669     ( mklog() << "Start Loft update "<<Name()<<" dt="<<dt<<" time="<<time<<std::endl).flush();
00670 #endif
00671     if (PartType::GetMobility() == 1 || time < dt ) {
00672       MType tmp = PartType::GetDHMat();
00673 #ifdef MOTION_DEBUG_ON
00674       ( mklog() << "\tlOFT " << Name() << " DH_Mat " << tmp << std::endl ).flush();
00675       ( mklog() << "\t\ttime " << time << " " << Name() << " Velocity DH_MAT\n\t\t" << (PartType::DHMat-PartType::DHMat_old)*(DataType(1.)/dt) << std::endl ).flush();
00676 #endif
00677       if ( SurfaceType::GetDeformation() == 1 ) {
00678         
00679         GenLoft();
00680         makeUnique();
00681       }
00682       prodMVP(tmp,GetCPoints_org(),GetCPoints());
00683       if (time < dt) {
00684         PartType::setDHMat_old(tmp);
00685       }
00686       numVelMVP(PartType::GetDHMat_old(),PartType::GetDHMat(),GetCPoints_org(),SurfaceType::GetVelocities(),dt);
00687       PartType::setDHMat_old(tmp);
00688       SurfaceType::measure();
00689     }
00690     else {
00691 
00692       if ( SurfaceType::GetDeformation() == 1 ) {
00693         
00694         GenLoft();
00695         makeUnique();
00696         SurfaceType::measure();
00697       }
00698     }
00699   }
00700 
00701   virtual void deformPart(DataType dt,DataType time) {
00702     
00703   }
00704 
00705   virtual void setLoftProfileTransform(int n, PType ptemp, PType translate1, PType scale1, PType rotate1, int invert ) {
00706     ( mklog() << "setLoftDeformation Profile "<< n << "    positions="<< positions[n]
00707                  << " scales="<<scales[n] <<" orientations="<< orientations[n] << std::endl ).flush();
00708 
00709     positions[n] = translate1;
00710     scales[n] = scale1;
00711     orientations[n] = rotate1;
00712   }
00713 
00714   virtual void Restart(std::ifstream& ifs, int& pos, double& t, double& dt) {
00715 #ifdef CHECKRESTART_DEBUG_ON
00716     ( mklog() << "Restarting Loft " << Name() << " pos = " << pos << std::endl ).flush();
00717 #endif
00718     SurfaceType::Restart(ifs,pos,t,dt);
00719     pos = ifs.tellg();
00720 #ifdef CHECKRESTART_DEBUG_ON
00721     ( mklog() << "\t_____continuing Restarting Loft " << Name() << t << " dt = " << dt
00722               << " pos = " << pos << std::endl ).flush();
00723 #endif
00724     ifs.seekg(pos);
00725     ifs.read((char*)&NProfiles,sizeof(int));
00726     for (register int i=0; i < NProfiles; i++) {
00727       ifs.read((char*)&profileTypes[i],sizeof(int));
00728     }
00729 #ifdef CHECKRESTART_DEBUG_ON
00730     ( mklog() << " profileTypes read " << NProfiles <<"\n" ).flush();
00731 #endif
00732     for (register int i=0; i < NProfiles; i++) {
00733       ifs.read((char*)&positions[i](0),sizeof(DataType));
00734           ifs.read((char*)&positions[i](1),sizeof(DataType));
00735           ifs.read((char*)&positions[i](2),sizeof(DataType));
00736           ifs.read((char*)&orientations[i](0),sizeof(DataType));
00737           ifs.read((char*)&orientations[i](1),sizeof(DataType));
00738           ifs.read((char*)&orientations[i](2),sizeof(DataType));
00739           ifs.read((char*)&scales[i](0),sizeof(DataType));
00740           ifs.read((char*)&scales[i](1),sizeof(DataType));
00741           ifs.read((char*)&scales[i](2),sizeof(DataType));
00742     }
00743     ifs.read((char*)&num_transforms,sizeof(int));
00744     for (register int i=0; i < num_transforms; i++) {
00745       ifs.read((char*)&transformList[i](0,0),sizeof(DataType));
00746           ifs.read((char*)&transformList[i](0,1),sizeof(DataType));
00747           ifs.read((char*)&transformList[i](0,2),sizeof(DataType));
00748           ifs.read((char*)&transformList[i](0,3),sizeof(DataType));
00749           ifs.read((char*)&transformList[i](1,0),sizeof(DataType));
00750           ifs.read((char*)&transformList[i](1,1),sizeof(DataType));
00751           ifs.read((char*)&transformList[i](1,2),sizeof(DataType));
00752           ifs.read((char*)&transformList[i](1,3),sizeof(DataType));
00753           ifs.read((char*)&transformList[i](2,0),sizeof(DataType));
00754           ifs.read((char*)&transformList[i](2,1),sizeof(DataType));
00755           ifs.read((char*)&transformList[i](2,2),sizeof(DataType));
00756           ifs.read((char*)&transformList[i](2,3),sizeof(DataType));
00757           ifs.read((char*)&transformList[i](3,0),sizeof(DataType));
00758           ifs.read((char*)&transformList[i](3,1),sizeof(DataType));
00759           ifs.read((char*)&transformList[i](3,2),sizeof(DataType));
00760           ifs.read((char*)&transformList[i](3,3),sizeof(DataType));
00761     }
00762     pos = ifs.tellg();
00763 #ifdef CHECKRESTART_DEBUG_ON
00764     ( mklog() << " transforms read pos=" << pos << std::endl ).flush();
00765 #endif
00766 
00767     for (register int i=0; i < NProfiles; i++) {
00768       profileSplines[i]->Restart(ifs,pos,t,dt);
00769     }
00770     ifs.read((char*)&num_spanSplines,sizeof(int));
00771     pos = ifs.tellg();
00772     for (register int i=0; i < num_spanSplines; i++) {
00773       spanSplines[i]->Restart(ifs,pos,t,dt);
00774     }
00775     pos = ifs.tellg();
00776 #ifdef CHECKRESTART_DEBUG_ON
00777     ( mklog()<< " spanSplines read\n" ).flush();
00778 #endif
00779   }
00780 
00781   virtual void Checkpointing(std::ofstream& ofs) {
00782     SurfaceType::Checkpointing(ofs);
00783 #ifdef CHECKRESTART_DEBUG_ON
00784     ( mklog() << "\nCHECKPOINTING Loft " << Name()  << " NProfiles " << NProfiles << std::endl ).flush();
00785 #endif
00786     ofs.write((char*)&NProfiles,sizeof(int));
00787     for (register int i=0; i < NProfiles; i++) {
00788       ofs.write((char*)&profileTypes[i],sizeof(int));
00789     }
00790     for (register int i=0; i < NProfiles; i++) {
00791       ofs.write((char*)&positions[i](0),sizeof(DataType));
00792           ofs.write((char*)&positions[i](1),sizeof(DataType));
00793           ofs.write((char*)&positions[i](2),sizeof(DataType));
00794           ofs.write((char*)&orientations[i](0),sizeof(DataType));
00795           ofs.write((char*)&orientations[i](1),sizeof(DataType));
00796           ofs.write((char*)&orientations[i](2),sizeof(DataType));
00797           ofs.write((char*)&scales[i](0),sizeof(DataType));
00798           ofs.write((char*)&scales[i](1),sizeof(DataType));
00799           ofs.write((char*)&scales[i](2),sizeof(DataType));
00800     }
00801     num_transforms=transformList.size();
00802     ofs.write((char*)&num_transforms,sizeof(int));
00803     for (register int i=0; i < num_transforms; i++) {
00804       ofs.write((char*)&transformList[i](0,0),sizeof(DataType));
00805           ofs.write((char*)&transformList[i](0,1),sizeof(DataType));
00806           ofs.write((char*)&transformList[i](0,2),sizeof(DataType));
00807           ofs.write((char*)&transformList[i](0,3),sizeof(DataType));
00808           ofs.write((char*)&transformList[i](1,0),sizeof(DataType));
00809           ofs.write((char*)&transformList[i](1,1),sizeof(DataType));
00810           ofs.write((char*)&transformList[i](1,2),sizeof(DataType));
00811           ofs.write((char*)&transformList[i](1,3),sizeof(DataType));
00812           ofs.write((char*)&transformList[i](2,0),sizeof(DataType));
00813           ofs.write((char*)&transformList[i](2,1),sizeof(DataType));
00814           ofs.write((char*)&transformList[i](2,2),sizeof(DataType));
00815           ofs.write((char*)&transformList[i](2,3),sizeof(DataType));
00816           ofs.write((char*)&transformList[i](3,0),sizeof(DataType));
00817           ofs.write((char*)&transformList[i](3,1),sizeof(DataType));
00818           ofs.write((char*)&transformList[i](3,2),sizeof(DataType));
00819           ofs.write((char*)&transformList[i](3,3),sizeof(DataType));
00820     }
00821     for (register int i=0; i < NProfiles; i++) {
00822       profileSplines[i]->Checkpointing(ofs);
00823     }
00824     num_spanSplines=spanSplines.size();
00825     ofs.write((char*)&num_spanSplines,sizeof(int));
00826     for (register int i=0; i < int(spanSplines.size()); i++) {
00827       spanSplines[i]->Checkpointing(ofs);
00828     }
00829 
00830   }
00831 
00832   virtual PType GetNetP() { return SurfaceType::GetNetP(); }
00833 
00834   virtual void addIntType(int val) { intType.push_back(val); }
00835 
00836 protected:
00837 
00838   int NProfiles, NChords, NSpans;
00839   std::vector<int> profileTypes, intType;
00840   std::vector<std::string> profileNames;
00841   std::vector<PType> positions,orientations,scales, reflections;
00842   std::vector<SplineCurveType *> profileSplines;
00843   std::vector<SplineCurveType *> spanSplines;
00844   std::vector<MType> transformList;
00845 
00846   int profileName_size, num_transforms, num_spanSplines, spline_intType;
00847 };
00848 
00849 #endif // LOFT_H