1 //============================================================================ 2 // Copyright (c) Kitware, Inc. 3 // All rights reserved. 4 // See LICENSE.txt for details. 5 // 6 // This software is distributed WITHOUT ANY WARRANTY; without even 7 // the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR 8 // PURPOSE. See the above copyright notice for more information. 9 //============================================================================ 10 11 #include <vtkm/cont/DataSet.h> 12 #include <vtkm/cont/Initialize.h> 13 #include <vtkm/filter/Streamline.h> 14 #include <vtkm/io/VTKDataSetReader.h> 15 #include <vtkm/io/VTKDataSetWriter.h> 16 17 // Example computing streamlines. 18 // An example vector field is available in the vtk-m data directory: magField.vtk 19 // Example usage: 20 // this will advect 200 particles 50 steps using a step size of 0.01 21 // 22 // Particle_Advection <path-to-data-dir>/magField.vtk vec 200 50 0.01 output.vtk 23 // 24 main(int argc,char ** argv)25int main(int argc, char** argv) 26 { 27 auto opts = vtkm::cont::InitializeOptions::DefaultAnyDevice; 28 auto config = vtkm::cont::Initialize(argc, argv, opts); 29 30 if (argc < 8) 31 { 32 std::cerr << "Usage: " << argv[0] 33 << "dataFile varName numSeeds numSteps stepSize outputFile [options]" << std::endl; 34 std::cerr << "where options are: " << std::endl << config.Usage << std::endl; 35 return -1; 36 } 37 38 std::string dataFile = argv[1]; 39 std::string varName = argv[2]; 40 vtkm::Id numSeeds = std::stoi(argv[3]); 41 vtkm::Id numSteps = std::stoi(argv[4]); 42 vtkm::FloatDefault stepSize = std::stof(argv[5]); 43 std::string outputFile = argv[6]; 44 45 vtkm::cont::DataSet ds; 46 47 if (dataFile.find(".vtk") != std::string::npos) 48 { 49 vtkm::io::VTKDataSetReader rdr(dataFile); 50 ds = rdr.ReadDataSet(); 51 } 52 else 53 { 54 std::cerr << "Unsupported data file: " << dataFile << std::endl; 55 return -1; 56 } 57 58 //create seeds randomly placed withing the bounding box of the data. 59 vtkm::Bounds bounds = ds.GetCoordinateSystem().GetBounds(); 60 std::vector<vtkm::Particle> seeds; 61 62 for (vtkm::Id i = 0; i < numSeeds; i++) 63 { 64 vtkm::Particle p; 65 vtkm::FloatDefault rx = (vtkm::FloatDefault)rand() / (vtkm::FloatDefault)RAND_MAX; 66 vtkm::FloatDefault ry = (vtkm::FloatDefault)rand() / (vtkm::FloatDefault)RAND_MAX; 67 vtkm::FloatDefault rz = (vtkm::FloatDefault)rand() / (vtkm::FloatDefault)RAND_MAX; 68 p.Pos[0] = static_cast<vtkm::FloatDefault>(bounds.X.Min + rx * bounds.X.Length()); 69 p.Pos[1] = static_cast<vtkm::FloatDefault>(bounds.Y.Min + ry * bounds.Y.Length()); 70 p.Pos[2] = static_cast<vtkm::FloatDefault>(bounds.Z.Min + rz * bounds.Z.Length()); 71 p.ID = i; 72 seeds.push_back(p); 73 } 74 auto seedArray = vtkm::cont::make_ArrayHandle(seeds, vtkm::CopyFlag::Off); 75 76 //compute streamlines 77 vtkm::filter::Streamline streamline; 78 79 streamline.SetStepSize(stepSize); 80 streamline.SetNumberOfSteps(numSteps); 81 streamline.SetSeeds(seedArray); 82 83 streamline.SetActiveField(varName); 84 auto output = streamline.Execute(ds); 85 86 vtkm::io::VTKDataSetWriter wrt(outputFile); 87 wrt.WriteDataSet(output); 88 89 return 0; 90 } 91