src/Plugins/DevelopedSurface/IO/vtkDevelopedSurface.cxx

   1 // Copyright (C) 2017  EDF R&D
   2 //
   3 // This library is free software; you can redistribute it and/or
   4 // modify it under the terms of the GNU Lesser General Public
   5 // License as published by the Free Software Foundation; either
   6 // version 2.1 of the License, or (at your option) any later version.
   7 //
   8 // This library is distributed in the hope that it will be useful,
   9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  11 // Lesser General Public License for more details.
  12 //
  13 // You should have received a copy of the GNU Lesser General Public
  14 // License along with this library; if not, write to the Free Software
  15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19 // Author : Anthony Geay (EDF R&D)
  20
  21 #include "vtkDevelopedSurface.h"
  22 #include "VTKToMEDMem.hxx"
  23
  24 #include "vtkAdjacentVertexIterator.h"
  25 #include "vtkIntArray.h"
  26 #include "vtkCellData.h"
  27 #include "vtkPointData.h"
  28 #include "vtkCylinder.h"
  29 #include "vtkNew.h"
  30 #include "vtkCutter.h"
  31 #include "vtkTransform.h"
  32
  33 #include "vtkStreamingDemandDrivenPipeline.h"
  34 #include "vtkUnstructuredGrid.h"
  35 #include  "vtkMultiBlockDataSet.h"
  36
  37 #include "vtkInformationStringKey.h"
  38 #include "vtkAlgorithmOutput.h"
  39 #include "vtkObjectFactory.h"
  40 #include "vtkMutableDirectedGraph.h"
  41 #include "vtkMultiBlockDataSet.h"
  42 #include "vtkDataSet.h"
  43 #include "vtkInformationVector.h"
  44 #include "vtkInformation.h"
  45 #include "vtkDataArraySelection.h"
  46 #include "vtkTimeStamp.h"
  47 #include "vtkInEdgeIterator.h"
  48 #include "vtkInformationDataObjectKey.h"
  49 #include "vtkExecutive.h"
  50 #include "vtkVariantArray.h"
  51 #include "vtkStringArray.h"
  52 #include "vtkDoubleArray.h"
  53 #include "vtkFloatArray.h"
  54 #include "vtkCharArray.h"
  55 #include "vtkUnsignedCharArray.h"
  56 #include "vtkDataSetAttributes.h"
  57 #include "vtkDemandDrivenPipeline.h"
  58 #include "vtkDataObjectTreeIterator.h"
  59 #include "vtkWarpScalar.h"
  60
  61 #include "MEDCouplingMemArray.hxx"
  62
  63 #include "VTKMEDTraits.hxx"
  64
  65 #include <map>
  66 #include <deque>
  67 #include <sstream>
  68 #include <algorithm>
  69
  70 vtkStandardNewMacro(vtkDevelopedSurface);
  71
  72 ///////////////////
  73
  74 template<class T>
  75 struct VTKTraits
  76 {
  77 };
  78
  79 template<>
  80 struct VTKTraits<double>
  81 {
  82   typedef vtkDoubleArray ArrayType;
  83 };
  84
  85 template<>
  86 struct VTKTraits<float>
  87 {
  88   typedef vtkFloatArray ArrayType;
  89 };
  90
  91 void ExtractInfo(vtkInformationVector *inputVector, vtkDataSet *& usgIn)
  92 {
  93   vtkInformation *inputInfo(inputVector->GetInformationObject(0));
  94   vtkDataSet *input(0);
  95   vtkDataSet *input0(vtkDataSet::SafeDownCast(inputInfo->Get(vtkDataObject::DATA_OBJECT())));
  96   vtkMultiBlockDataSet *input1(vtkMultiBlockDataSet::SafeDownCast(inputInfo->Get(vtkDataObject::DATA_OBJECT())));
  97   if(input0)
  98     input=input0;
  99   else
 100     {
 101       if(!input1)
 102         throw MZCException("Input dataSet must be a DataSet or single elt multi block dataset expected !");
 103       if(input1->GetNumberOfBlocks()!=1)
 104         throw MZCException("Input dataSet is a multiblock dataset with not exactly one block ! Use MergeBlocks or ExtractBlocks filter before calling this filter !");
 105       vtkDataObject *input2(input1->GetBlock(0));
 106       if(!input2)
 107         throw MZCException("Input dataSet is a multiblock dataset with exactly one block but this single element is NULL !");
 108       vtkDataSet *input2c(vtkDataSet::SafeDownCast(input2));
 109       if(!input2c)
 110         throw MZCException("Input dataSet is a multiblock dataset with exactly one block but this single element is not a dataset ! Use MergeBlocks or ExtractBlocks filter before calling this filter !");
 111       input=input2c;
 112     }
 113   if(!input)
 114     throw MZCException("Input data set is NULL !");
 115   vtkPointData *att(input->GetPointData());
 116   if(!att)
 117     throw MZCException("Input dataset has no point data attribute ! Impossible to deduce a developed surface on it !");
 118   usgIn=input;
 119 }
 120
 121 class vtkDevelopedSurface::vtkInternals
 122 {
 123 public:
 124   vtkNew<vtkCutter> Cutter;
 125 };
 126
 127 ////////////////////
 128
 129 vtkDevelopedSurface::vtkDevelopedSurface():_cyl(nullptr),Internal(new vtkInternals),InvertStatus(false),OffsetInRad(0.)
 130 {
 131   //this->RegisterFilter(this->Internal->Cutter.GetPointer());
 132 }
 133
 134 vtkDevelopedSurface::~vtkDevelopedSurface()
 135 {
 136   delete this->Internal;
 137 }
 138
 139 void vtkDevelopedSurface::SetInvertWay(bool invertStatus)
 140 {
 141   this->InvertStatus=invertStatus;
 142   this->Modified();
 143 }
 144
 145 void vtkDevelopedSurface::SetThetaOffset(double offsetInDegrees)
 146 {
 147   double tmp(std::min(offsetInDegrees,180.));
 148   tmp=std::max(tmp,-180.);
 149   this->OffsetInRad=tmp/180.*M_PI;
 150   this->Modified();
 151 }
 152
 153 int vtkDevelopedSurface::RequestInformation(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
 154 {
 155   //std::cerr << "########################################## vtkDevelopedSurface::RequestInformation ##########################################" << std::endl;
 156   try
 157     {
 158       vtkDataSet *usgIn(0);
 159       ExtractInfo(inputVector[0],usgIn);
 160     }
 161   catch(MZCException& e)
 162     {
 163       std::ostringstream oss;
 164       oss << "Exception has been thrown in vtkDevelopedSurface::RequestInformation : " << e.what() << std::endl;
 165       if(this->HasObserver("ErrorEvent") )
 166         this->InvokeEvent("ErrorEvent",const_cast<char *>(oss.str().c_str()));
 167       else
 168         vtkOutputWindowDisplayErrorText(const_cast<char *>(oss.str().c_str()));
 169       vtkObject::BreakOnError();
 170       return 0;
 171     }
 172   return 1;
 173 }
 174
 175 std::vector<int> UnWrapByDuplicatingNodes(vtkCellArray *ca, vtkIdType& offset, const MEDCoupling::DataArrayDouble *thetas)
 176 {
 177   std::vector<int> ret;
 178   vtkIdType nbCells(ca->GetNumberOfCells());
 179   vtkIdType *conn(ca->GetPointer());
 180   const double *tptr(thetas->begin());
 181   for(vtkIdType i=0;i<nbCells;i++)
 182     {
 183       vtkIdType sz(*conn++);
 184       double ma(-std::numeric_limits<double>::max()),mi(std::numeric_limits<double>::max());
 185       for(vtkIdType j=0;j<sz;j++)
 186         {
 187           double tmp_theta(tptr[conn[j]]);
 188           ma=std::max(ma,tmp_theta);
 189           mi=std::min(mi,tmp_theta);
 190         }
 191       if(ma-mi>M_PI)
 192         {
 193           for(vtkIdType j=0;j<sz;j++)
 194             {
 195               double tmp_theta(tptr[conn[j]]);
 196               if(tmp_theta<=M_PI)
 197                 {
 198                   ret.push_back(conn[j]);
 199                   conn[j]=offset++;
 200                 }
 201             }
 202         }
 203       conn+=sz;
 204     }
 205   return ret;
 206 }
 207
 208 template<class T>
 209 void DealArray(vtkDataSetAttributes *pd, int pos, typename MEDFileVTKTraits<T>::VtkType *arr, std::vector<int>& nodeSel)
 210 {
 211   int nbc(arr->GetNumberOfComponents());
 212   std::size_t nbt(nodeSel.size());
 213   vtkSmartPointer< typename MEDFileVTKTraits<T>::VtkType > newArr;
 214   newArr.TakeReference(MEDFileVTKTraits<T>::VtkType::New());
 215   newArr->SetNumberOfComponents(nbc);
 216   newArr->SetNumberOfTuples(nbt);
 217   T *ptr(newArr->GetPointer(0));
 218   const T *inPtr(arr->GetPointer(0));
 219   for(std::size_t i=0;i<nbt;i++)
 220     {
 221       ptr=std::copy(inPtr+nbc*nodeSel[i],inPtr+nbc*(nodeSel[i]+1),ptr);
 222     }
 223   newArr->SetName(arr->GetName());
 224   arr->DeepCopy(newArr);
 225 }
 226
 227 void ToDouble(vtkDataArray *coords, vtkSmartPointer<vtkDoubleArray>& coordsOut)
 228 {
 229   vtkDoubleArray *coords2(vtkDoubleArray::SafeDownCast(coords));
 230   vtkFloatArray *coords3(vtkFloatArray::SafeDownCast(coords));
 231   if(!coords2 && !coords3)
 232     throw MZCException("Input coordinates are neither float64 or float32 !");
 233   //
 234   if(coords2)
 235     {
 236       coordsOut.TakeReference(coords2);
 237       coords2->Register(0);
 238     }
 239   else
 240     {
 241       coordsOut.TakeReference(vtkDoubleArray::New());
 242       coordsOut->SetNumberOfComponents(3);
 243       vtkIdType nbTuples(coords3->GetNumberOfTuples());
 244       coordsOut->SetNumberOfTuples(nbTuples);
 245       std::copy(coords3->GetPointer(0),coords3->GetPointer(0)+3*nbTuples,coordsOut->GetPointer(0));
 246     }
 247 }
 248
 249 void dealWith(vtkPolyData *outdata, const double center[3], const double axis[3], double radius, double eps, bool invertThetaInc, double offsetInRad)
 250 {
 251   vtkDataArray *coords(outdata->GetPoints()->GetData());
 252   if(coords->GetNumberOfComponents()!=3)
 253     throw MZCException("Input coordinates are expected to have 3 components !");
 254   //
 255   vtkIdType nbNodes(coords->GetNumberOfTuples());
 256   if(nbNodes==0)
 257     throw MZCException("No points -> impossible to develop anything !");
 258   //
 259   vtkSmartPointer<vtkDoubleArray> zeCoords;
 260   ToDouble(coords,zeCoords);
 261   //
 262   double axis_cross_Z[3]={axis[1],-axis[0],0.};
 263   double n_axis(sqrt(axis_cross_Z[0]*axis_cross_Z[0]+axis_cross_Z[1]*axis_cross_Z[1]));
 264   if(n_axis>eps)
 265     {
 266       double ang(asin(n_axis));
 267       if(axis[2]<0.)
 268         ang=M_PI-ang;
 269       MEDCoupling::DataArrayDouble::Rotate3DAlg(center,axis_cross_Z,ang,nbNodes,zeCoords->GetPointer(0),zeCoords->GetPointer(0));
 270     }
 271   //
 272   MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> c_cyl;
 273   {
 274     MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> cc(MEDCoupling::DataArrayDouble::New()); cc->alloc(nbNodes,3);
 275     double *ccPtr(cc->getPointer());
 276     const double *zeCoordsPtr(zeCoords->GetPointer(0));
 277     for(vtkIdType i=0;i<nbNodes;i++,ccPtr+=3,zeCoordsPtr+=3)
 278       {
 279         std::transform(zeCoordsPtr,zeCoordsPtr+3,center,ccPtr,std::minus<double>());
 280       }
 281     c_cyl=cc->fromCartToCyl();
 282   }
 283   MEDCoupling::MCAuto<MEDCoupling::MEDFileData> mfd(MEDCoupling::MEDFileData::New());
 284   WriteMEDFileFromVTKDataSet(mfd,outdata,{},0.,0);
 285   bool a;
 286   {
 287     MEDCoupling::MEDFileMeshes *ms(mfd->getMeshes());
 288     if(ms->getNumberOfMeshes()!=1)
 289       throw MZCException("Unexpected number of meshes !");
 290     MEDCoupling::MEDFileMesh *mm(ms->getMeshAtPos(0));
 291     MEDCoupling::MEDFileUMesh *mmu(dynamic_cast<MEDCoupling::MEDFileUMesh *>(mm));
 292     if(!mmu)
 293       throw MZCException("Expecting unstructured one !");
 294     MEDCoupling::MCAuto<MEDCoupling::MEDCouplingUMesh> m0(mmu->getMeshAtLevel(0));
 295     {
 296       int v(0);
 297       MEDCoupling::MCAuto<MEDCoupling::DataArrayInt> c0s(m0->getCellIdsLyingOnNodes(&v,&v+1,false));
 298       if(c0s->empty())
 299         throw MZCException("Orphan node 0 !");
 300       std::vector<int> nodes0;
 301       m0->getNodeIdsOfCell(c0s->getIJ(0,0),nodes0);
 302       MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> tmp0(c_cyl->selectByTupleIdSafe(nodes0.data(),nodes0.data()+nodes0.size()));
 303       tmp0=tmp0->keepSelectedComponents({1});
 304       double tmp(tmp0->getMaxAbsValueInArray());
 305       a=tmp>0.;
 306     }
 307   }
 308   //
 309   constexpr double EPS_FOR_RADIUS=1e-2;
 310   MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> rs(c_cyl->keepSelectedComponents({0}));
 311   if(!rs->isUniform(radius,radius*EPS_FOR_RADIUS))
 312     {
 313       double mi(rs->getMinValueInArray()),ma(rs->getMaxValueInArray());
 314       std::ostringstream oss; oss << "Looks not really a cylinder within given precision ! Range is [" << mi << "," << ma << "] expecting " << radius << " within precision of " << radius*EPS_FOR_RADIUS << " !";
 315       throw MZCException(oss.str());
 316     }
 317   double tetha0(c_cyl->getIJ(0,1));
 318   {
 319     double *ccylptr(c_cyl->getPointer()+1);
 320     double mi02(std::numeric_limits<double>::max());
 321     for(vtkIdType i=0;i<nbNodes;i++,ccylptr+=3)
 322       {
 323         *ccylptr-=tetha0+offsetInRad;
 324         mi02=std::min(mi02,*ccylptr);
 325       }
 326     if(mi02>0.)
 327       {
 328         ccylptr=c_cyl->getPointer()+1;
 329         for(vtkIdType i=0;i<nbNodes;i++,ccylptr+=3)
 330           {
 331             if(*ccylptr>=2*M_PI)
 332               {
 333                 *ccylptr+=-2*M_PI;
 334               }
 335           }
 336       }
 337   }
 338   {
 339     MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> c_cyl_2(c_cyl->keepSelectedComponents({1}));
 340     c_cyl_2->abs();
 341     MEDCoupling::MCAuto<MEDCoupling::DataArrayInt> poses(c_cyl_2->findIdsInRange(0.,eps));
 342     c_cyl->setPartOfValuesSimple3(0.,poses->begin(),poses->end(),1,2,1);
 343   }
 344   //
 345   if(a ^ (!invertThetaInc))
 346     {
 347       MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> tmp(c_cyl->keepSelectedComponents({1}));
 348       tmp=tmp->negate();
 349       std::for_each(tmp->getPointer(),tmp->getPointer()+tmp->getNumberOfTuples(),[](double& v) { if(v==-0.) v=0.; });
 350       c_cyl->setPartOfValues1(tmp,0,nbNodes,1,1,2,1);
 351     }
 352   MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> c_cyl_post(c_cyl->keepSelectedComponents({1}));
 353   {
 354     double *c_cyl_post_ptr(c_cyl_post->getPointer());
 355     for(vtkIdType i=0;i<nbNodes;i++,c_cyl_post_ptr++)
 356       {
 357         if(*c_cyl_post_ptr<0.)
 358           *c_cyl_post_ptr+=2*M_PI;
 359       }
 360   }
 361   vtkCellArray *cb(outdata->GetPolys());
 362   vtkIdType offset(nbNodes);
 363   std::vector<int> dupNodes(UnWrapByDuplicatingNodes(cb,offset,c_cyl_post));
 364   //
 365   MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> c_cyl_post2(c_cyl_post->selectByTupleId(dupNodes.data(),dupNodes.data()+dupNodes.size()));
 366   c_cyl_post2->applyLin(1.,2*M_PI);
 367   c_cyl_post=MEDCoupling::DataArrayDouble::Aggregate(c_cyl_post,c_cyl_post2);
 368   MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> z0(c_cyl->keepSelectedComponents({2}));
 369   MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> z1(z0->selectByTupleId(dupNodes.data(),dupNodes.data()+dupNodes.size()));
 370   z0=MEDCoupling::DataArrayDouble::Aggregate(z0,z1);
 371   //
 372   std::size_t outNbNodes(z0->getNumberOfTuples());
 373   vtkSmartPointer<vtkDoubleArray> zeCoords2;
 374   zeCoords2.TakeReference(vtkDoubleArray::New());
 375   zeCoords2->SetNumberOfComponents(3);
 376   zeCoords2->SetNumberOfTuples(outNbNodes);
 377   {
 378     const double *tptr(c_cyl_post->begin()),*zptr(z0->begin());
 379     double *outPtr(zeCoords2->GetPointer(0));
 380     for(std::size_t i=0;i<outNbNodes;i++,tptr++,zptr++,outPtr+=3)
 381       {
 382         outPtr[0]=radius*(*tptr);
 383         outPtr[1]=(*zptr);
 384         outPtr[2]=0.;
 385       }
 386   }
 387   //
 388   outdata->GetPoints()->SetData(zeCoords2);
 389   // now post process nodes
 390   std::vector<int> nodeSel(nbNodes+dupNodes.size());
 391   {
 392     int cnt(0);
 393     std::for_each(nodeSel.begin(),nodeSel.begin()+nbNodes,[&cnt](int& v){ v=cnt++; });
 394     std::copy(dupNodes.begin(),dupNodes.end(),nodeSel.begin()+nbNodes);
 395   }
 396   vtkDataSetAttributes *pd(outdata->GetPointData());
 397   int nba(pd->GetNumberOfArrays());
 398   for(int i=0;i<nba;i++)
 399     {
 400       vtkDataArray *arr(pd->GetArray(i));
 401       {
 402         vtkIntArray *arr0(vtkIntArray::SafeDownCast(arr));
 403         if(arr0)
 404           {
 405             DealArray<int>(pd,i,arr0,nodeSel);
 406             continue;
 407           }
 408       }
 409       {
 410         vtkFloatArray *arr0(vtkFloatArray::SafeDownCast(arr));
 411         if(arr0)
 412           {
 413             DealArray<float>(pd,i,arr0,nodeSel);
 414             continue;
 415           }
 416       }
 417       {
 418         vtkDoubleArray *arr0(vtkDoubleArray::SafeDownCast(arr));
 419         if(arr0)
 420           {
 421             DealArray<double>(pd,i,arr0,nodeSel);
 422             continue;
 423           }
 424       }
 425     }
 426 }
 427
 428 int vtkDevelopedSurface::RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
 429 {
 430   //std::cerr << "########################################## vtkDevelopedSurface::RequestData        ##########################################" << std::endl;
 431   try
 432     {
 433       if(!_cyl)
 434         throw MZCException("No cylinder object as cut function !");
 435       double center[3],axis[3],radius;
 436       vtkAbstractTransform* trf(_cyl->GetTransform());
 437       {
 438         _cyl->GetCenter(center);
 439         _cyl->GetAxis(axis[0],axis[1],axis[2]);
 440         radius=_cyl->GetRadius();
 441       }
 442       if(trf)
 443         {
 444           double axis3[3]={center[0]+0.,center[1]+1.,center[2]+0.},axis4[3];
 445           trf->TransformPoint(axis3,axis4);
 446           std::transform(axis4,axis4+3,center,axis,[](double a, double b) { return b-a; });
 447           axis[1]=-axis[1];
 448           if(std::isnan(axis[0]) && std::isnan(axis[1]) && std::isnan(axis[2]))
 449             { axis[0]=0.; axis[1]=-1.; axis[2]=0.; }
 450         }
 451       //std::cerr << trf << " jjj " << axis[0] << " " << axis[1] << " " << axis[2]  << " : " << center[0] <<  " " << center[1] << " " << center[2] <<  " " " " << " -> " << radius << std::endl;
 452       vtkDataSet *usgIn(0);
 453       ExtractInfo(inputVector[0],usgIn);
 454       vtkSmartPointer<vtkPolyData> outData;
 455       {
 456         vtkNew<vtkCutter> Cutter;
 457         Cutter->SetInputData(usgIn);
 458         Cutter->SetCutFunction(_cyl);
 459         Cutter->Update();
 460         vtkDataSet *zeComputedOutput(Cutter->GetOutput());
 461         vtkPolyData *zeComputedOutput2(vtkPolyData::SafeDownCast(zeComputedOutput));
 462         if(!zeComputedOutput2)
 463           throw MZCException("Unexpected output of cutter !");
 464         outData.TakeReference(zeComputedOutput2);
 465         zeComputedOutput2->Register(0);
 466       }
 467       if(outData->GetNumberOfCells()==0)
 468         return 1;// no cells -> nothing to do
 469       //
 470       dealWith(outData,center,axis,radius,1e-7,this->InvertStatus,this->OffsetInRad);
 471       //finish
 472       vtkInformation *outInfo(outputVector->GetInformationObject(0));
 473       vtkPolyData *output(vtkPolyData::SafeDownCast(outInfo->Get(vtkDataObject::DATA_OBJECT())));
 474       output->ShallowCopy(outData);
 475     }
 476   catch(MZCException& e)
 477     {
 478       std::ostringstream oss;
 479       oss << "Exception has been thrown in vtkDevelopedSurface::RequestInformation : " << e.what() << std::endl;
 480       if(this->HasObserver("ErrorEvent") )
 481         this->InvokeEvent("ErrorEvent",const_cast<char *>(oss.str().c_str()));
 482       else
 483         vtkOutputWindowDisplayErrorText(const_cast<char *>(oss.str().c_str()));
 484       vtkObject::BreakOnError();
 485       return 0;
 486     }
 487   return 1;
 488 }
 489
 490 int vtkDevelopedSurface::FillOutputPortInformation( int vtkNotUsed(port), vtkInformation* info)
 491 {
 492   info->Set(vtkDataObject::DATA_TYPE_NAME(), "vtkPolyData");
 493   return 1;
 494 }
 495
 496
 497 void vtkDevelopedSurface::PrintSelf(ostream& os, vtkIndent indent)
 498 {
 499   this->Superclass::PrintSelf(os, indent);
 500 }
 501
 502 void vtkDevelopedSurface::SetCutFunction(vtkImplicitFunction* func)
 503 {
 504   vtkCylinder *cyl(vtkCylinder::SafeDownCast(func));
 505   if(cyl)
 506     {
 507       _cyl=cyl;
 508       this->Modified();
 509     }
 510 }
 511
 512 vtkMTimeType vtkDevelopedSurface::GetMTime()
 513 {
 514   vtkMTimeType maxMTime = this->Superclass::GetMTime(); // My MTime
 515   if(_cyl)
 516     {
 517       maxMTime=std::max(maxMTime,_cyl->GetMTime());
 518     }
 519   return maxMTime;
 520 }