Il est possible d'écrire un programme dédié à la construction d'un arbre en utilisant la bibliothèque Triskele. Ce programme prend 3 arguments :

• le nom du fichier contenant l'image en entrée qui construira un arbre MAX
• le nom du fichier qui contiendra l'image produite
• une liste de valeurs (séparée par des virgules) pour réaliser les élagages

Il faut commencer par les déclarations d'usage.

La première inclusion sera utilisée pour convertir la liste de paramètre en valeur de seuil pour élagueur l'arbre.

Les autres ne sont pas indispensables, car elles seront dans tous les cas incluses par Triskele.

// is_any_of
#include <boost/algorithm/string.hpp>
 
// pixels type
#include "gdal.h"
// log functions
#include "obelixDebug.hpp"
// Point and Size (2D or 3D)
#include "obelixGeo.hpp"
// read/write images
#include "IImage.hpp"
// nodata pixels definition in images
#include "Border.hpp"

De même, on peut déclarer les classes Triskele qui pourront être utilisées.

#include "ArrayTree/triskeleArrayTreeBase.hpp"
#include "ArrayTree/ArrayTreeBuilder.hpp"
#include "Attributes/WeightAttributes.hpp"
#include "Attributes/AreaAttributes.hpp"
#include "AttributeProfiles.hpp"

En fonction de ce qui a été déclaré, il peut être utile de simplifier les appels en omettant le contexte d'utilisation (par exemple en évitant de préfixer avec std::).

using namespace std;
using namespace boost;
using namespace obelix;
using namespace obelix::triskele;

Dans notre exemple, nous nous limiterons à des images 2D et la construction d'un arbre MAX.

static const int GeoDimT (2);
static const TreeType treeType = MAX;

Nous déclarons la fonction filtre qui réalisera l'élagage et l'écriture du résultat. Elle sera donc dépendante de la nature des pixels lus InPixelT, de l'image en entrée et de la liste des seuils de coupure.

template<typename PixelT> void filter (const IImage<GeoDimT> &inputImage, IImage<GeoDimT> &outputImage, const vector <DimImgPack> &thresholds);

Le point d'entré du programme, commencera par analyser les arguments

• argv[0] : nom du programme
• argv[1] : nom du fichier contenant l'image en entrée
• argv[2] : nom du fichier contenant l'image produite
• argv[3] : liste des seuils séparé par des virgules

int
main (int argc, char** argv, char** envp) {
  // get parameters
  if (argc != 4) {
    cout << "Usage: " << argv[0] << " input-filename outpout-filename thresholds,..." << endl;
    exit (1);
  }
  // set thresholds vetcor
  vector <string> values;
  split (values, argv[3], is_any_of (","));
  vector <DimImgPack> thresholds (values.size ());
  transform (values.begin (), values.end (), thresholds.begin (), [] (const string &val) { return stoi (val); });

Ensuite nous lisons les métadonnées de l'image pour déterminer la nature des pixels et créer le fichier qui contiendra l'image produite.

  IImage<GeoDimT> inputImage (argv[1]);
  // get metadata
  inputImage.readImage ();
  // get image size (width x height)
  const Size<GeoDimT> size (inputImage.getSize ());
 
  IImage<GeoDimT> outputImage (argv[2]);
  // set metadata
  outputImage.createImage (size, inputImage.getDataType (), thresholds.size (),
			   inputImage, NullPoint2D);

Enfin nous réalisons le traitement en invoquant la fonction filtre avec le bon type.

  switch (inputImage.getDataType ()) {
  case GDT_Byte:
    filter<uint8_t> (inputImage, outputImage, thresholds); break;
  case GDT_UInt16:
    filter<uint16_t> (inputImage, outputImage, thresholds); break;
  case GDT_Int16:
    filter<int16_t> (inputImage, outputImage, thresholds); break;
  case GDT_UInt32:
    filter<uint32_t> (inputImage, outputImage, thresholds); break;
  case GDT_Int32:
    filter<int32_t> (inputImage, outputImage, thresholds); break;
  case GDT_Float32:
    filter<float> (inputImage, outputImage, thresholds); break;
  case GDT_Float64:
    filter<double> (inputImage, outputImage, thresholds); break;
 
  default :
    cerr << "unknown type!" << endl; break;
    return 1;
  }
  return 0;
}

Voici la fonction principale du programme.

template<typename InPixelT>
void
filter (const IImage<GeoDimT> &inputImage, IImage<GeoDimT> &outputImage, const vector <DimImgPack> &thresholds) {
 
  // switch on debug mode
  // Log::debug = true;

Nous lisons les pixels de l'image.

  // empty raster
  Raster<InPixelT, GeoDimT> raster;
  // true image size (width x height)
  const Size<GeoDimT> size (inputImage.getSize ());
  // read first band (0) in 2D mode from origine [0,0] to end [width, height].
  inputImage.readBand (raster, 0, NullPoint2D, size);

Nous déclarons l'objet arbre qui contiendra le chaînage des feuilles vers la racine.

  // no border (i.e. all pixels are take in account)
  const Border<GeoDimT> border (size, false);
  // neighbors take in account (default connectivity C4)
  const GraphWalker<GeoDimT> graphWalker (border);
  // tree builder base on raster, connectivity and type of tree
  ArrayTreeBuilder<PixelT, PixelT, GeoDimT> atb (raster, graphWalker, treeType);
  // get the number of core
  const DimCore coreCount (boost::thread::hardware_concurrency ());
  // declare empty tree and set the number of thread for all algorithms
  Tree<GeoDimT> tree (coreCount);
  // reserve nodes weight attributes
  WeightAttributes<PixelT, GeoDimT> weightAttributes (tree, getDecrFromTreetype (treeType));

Nous déclenchons la construction de l'arbre.

• ATB contient les éléments de l'images (matrice de pixels, connectivité, type d'arbre)
• tree coquille vide qui va contenir l'arbre
• weightAttributes est le premier tableau d'attributs associés aux nœuds (niveaux de gris)

  atb.buildTree (tree, weightAttributes);

Nous créons un second tableau d'attributs associés aux nœuds (taille en pixels).

  // create area attribute
  const AreaAttributes<GeoDimT> areaAttributes (tree);

récupération des cardinalités.

  // output channels count
  const DimChannel outputChannelCount (thresholds.size ());
  // number of pixels (width x height)
  const DimImg pixelsCount (size.getPixelsCount ());

Nous produisons des bandes de sorti (élagage de l'arbre)

  // write node grayscale as attribute profiles
  AttributeProfiles<PixelT, GeoDimT> attributeProfiles (tree);
  // copy weightAttributes to attributeProfiles
  atb.setAttributProfiles (attributeProfiles);

Nous réalisons toutes les coupures pour tous les seuils définis en arguments du programme.

  // result allocation memory
  vector <vector <OutPixelT> > allBands (outputChannelCount, vector<OutPixelT> (pixelsCount, 0));
  // do all cuts
  areaAttributes.cut (allBands, attributeProfiles, thresholds);

Nous recopions chaque matrice de pixels dans dans le fichier de sorti.

  // for all cuts
  for (DimChannel channel (0); channel < outputChannelCount; ++channel)
    // write output raster
    outputImage.writeBand (allBands[channel].data (), channel);
}

Si l'on fait l'hypthèse d'être dans le répertoire build avec la production d'une bibliothèse libtriskeledebug.so a partir de l'option Debug de cmake, voici a quoi peut ressembler la commende de compilation :

g++ -DNO_OTB -I../triskele/include -I/usr/include/hdf5/serial -I/usr/include/gdal -pthread -g -DENABLE_LOG -std=gnu++11 ../triskele/test/TestMaxTree.cpp -o testMaxTree  -L . -ltriskeledebug -lhdf5_serial -lhdf5_serial_cpp -lstdc++ -lpthread -lboost_system -lboost_chrono -lboost_thread -lboost_program_options -lboost_date_time -lboost_serialization -lboost_filesystem -lboost_unit_test_framework -lgdal -ltbb -ltinyxml

Sinon, il est possible d'enrichir le fichier CMakeLists.txt de la manière suivante :

  message ("sample max tree")
  add_executable (testMaxTree "test/TestMaxTree.cpp")
  target_compile_definitions (testMaxTree PUBLIC -DNO_OTB)
  target_link_libraries (testMaxTree ${allLibs})

Voici un exemple d'utilisation

./testMaxTree strasbourg.tif output.tif 10,100,10000