diff --git a/examples/images/example2DBijectiveRotations.cpp b/examples/images/example2DBijectiveRotations.cpp
new file mode 100644
index 0000000000..5a77fab0dc
--- /dev/null
+++ b/examples/images/example2DBijectiveRotations.cpp
@@ -0,0 +1,147 @@
+/**
+* This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as
+ * published by the Free Software Foundation, either version 3 of the
+ * License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see .
+ *
+ **/
+
+#pragma once
+
+/**
+* @file CBDR.h
+ * @author S. Breuils, J.O. Lachaud, D. Coeurjolly ; stephane.breuils@univ-smb.fr
+ * @ingroup Examples
+ * @date 2024/08
+ *
+ * This file is part of the DGtal library.
+ */
+#include
+#include "DGtal/images/ImageContainerBySTLVector.h"
+#include "DGtal/helpers/StdDefs.h"
+#include "DGtal/base/Common.h"
+#include "DGtal/io/Color.h"
+#include "DGtal/io/readers/PPMReader.h"
+#include "DGtal/io/writers/GenericWriter.h"
+#include "DGtal/images/RigidTransformation2D.h"
+
+
+//////////////////////////////////////////////////////////////////////////////
+// Inclusions
+#include "DGtal/images/bijectiverotations/QSH.h"
+#include "DGtal/images/bijectiverotations/CDLR.h"
+#include "DGtal/images/bijectiverotations/RBC.h"
+#include "DGtal/images/bijectiverotations/OTC.h"
+#include "DGtal/images/bijectiverotations/CBDR.h"
+#include "DGtal/images/bijectiverotations/Rotationtables.h"
+
+
+using namespace std;
+using namespace DGtal;
+using namespace functors;
+using namespace Z2i;
+
+std::vector supportedBijectiveRotation = {
+ "OTC", "CBDR", "CDLR", "QSH" , "RBC"
+};
+
+// Generic function to handle rotation based on type
+template
+TImage
+performRotation( const TImage& img, TBijectiveRotation& obj )
+{
+ const auto& sbr = supportedBijectiveRotation;
+ std::string structName = obj.tostring();
+ TImage output( img.domain() );
+ if ( std::find( sbr.begin(), sbr.end(), structName ) != sbr.end() )
+ {
+ trace.beginBlock ( obj.tostring() );
+ output = obj.rotateImage( img );
+ trace.endBlock ();
+ }
+ else
+ throw std::runtime_error("Unsupported bijective rotation : " + structName);
+ return output;
+}
+
+void usage( const std::string& cmd )
+{
+ std::cout << "Usage: " << cmd << " [] [black|*white*] [detect]" << "\n"
+ << "\t Compute the rotated by an angle (in degrees) \n"
+ << "\t - in { OTC, CBDR, CDLR, QSH, *RBC* }\n"
+ << std::endl;
+}
+
+
+int main( int argc, char* argv[] )
+{
+ if ( argc < 3 ) { usage( argv[ 0 ] ); return 1; }
+ std::string fname = argv[ 1 ];
+ double degree = std::atof( argv[ 2 ] );
+ double angle = degree * M_PI / 180.0;
+ std::string method = ( argc > 3 ) ? argv[ 3 ] : "RBC";
+
+ typedef ImageContainerBySTLVector< Domain, Color > Image;
+ typedef ForwardRigidTransformation2D < Space > ForwardTrans;
+ typedef DomainRigidTransformation2D < Domain, ForwardTrans > MyDomainTransformer;
+ typedef MyDomainTransformer::Bounds Bounds;
+
+ Image image = DGtal::PPMReader::importPPM ( fname );
+ Image output = image;
+ int W=image.domain().myUpperBound[0]+1;
+ int H=image.domain().myUpperBound[1]+1;
+ trace.info() << "Image has size " << W << "x" << H << std::endl;
+ Point c(W/2, H/2);
+
+ if ( method == "QSH" )
+ {
+ DGtal::QSH rotQSH(angle,c);
+ output = performRotation> (image,rotQSH);
+ }
+ else if ( method == "CDLR" )
+ {
+ auto linf = std::make_shared>>();
+ DGtal::CDLR > rotDSL(angle, c, linf);
+ output = performRotation> (image,rotDSL);
+ }
+ else if ( method == "RBC" )
+ {
+ DGtal::RBC_vec rot_rbcvec(2*max(W,H));
+ rot_rbcvec.setAngle() = angle;
+ rot_rbcvec.center() = c;
+ DGtal::RBC rot_RBC(rot_rbcvec,angle,c);
+ output = performRotation> (image,rot_RBC);
+ }
+ else if ( method == "CBDR" )
+ {
+ auto linfCBDR = std::make_shared>>();
+ const int n = 4;
+ const int kmax=30;
+ DGtal::CBDR rot_CBDR(angle,c,n,kmax,linfCBDR,true);
+ output = performRotation> (image,rot_CBDR);
+ }
+
+ else if ( method == "OTC" )
+ {
+ int rwidth = 2;
+ std::vector< std::vector< int > > tableOTC = DGtal::functions::loadOTCTable("../tables/",rwidth);
+ DGtal::OTC rot_OTC( tableOTC, rwidth, c, W, H );
+ int alpha = int( round( angle * 180.0 / M_PI ) );
+ rot_OTC.set_angle( alpha ); // in degree
+ output = performRotation> (image,rot_OTC);
+ }
+
+
+ std::string out_fname = "rotated-image-" + method + "-" + std::to_string( int(degree) );
+ out_fname += ".ppm";
+ output >> out_fname;
+ return 0;
+}
\ No newline at end of file
diff --git a/src/DGtal/images/bijectiveRotations/CBDR.h b/src/DGtal/images/bijectiveRotations/CBDR.h
new file mode 100644
index 0000000000..06c8dc7795
--- /dev/null
+++ b/src/DGtal/images/bijectiveRotations/CBDR.h
@@ -0,0 +1,261 @@
+/**
+* This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as
+ * published by the Free Software Foundation, either version 3 of the
+ * License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see .
+ *
+ **/
+
+#pragma once
+
+/**
+* @file CBDR.h
+ * @author S. Breuils, J.O. Lachaud, D. Coeurjolly
+ *
+ * @date 2024/08
+ *
+ * This file is part of the DGtal library.
+ */
+
+#if defined(CBDR_RECURSES)
+#error Recursive header files inclusion detected in CBDR.h
+#else // defined(CBDR_RECURSES)
+/** Prevents recursive inclusion of headers. */
+#define CBDR_RECURSES
+
+#if !defined CBDR_h
+/** Prevents repeated inclusion of headers. */
+#define CBDR_h
+
+//////////////////////////////////////////////////////////////////////////////
+// Inclusions
+#include "GAVector.h"
+#include "CBDRSolver.h"
+#include "CBDRFastSolver.h"
+
+#include "Rotationtables.h"
+#include "NBijectiveReflectionGenerator.h"
+
+
+namespace DGtal {
+ /**
+ * Description of template struct CBDR
+ * \brief CBDR : Composition of Bijective Digitized Reflections,
+ * @tparam TSpace a 2 dimensional space.
+ * @tparam TInputValue type of the input point e.g., TSpace::RealPoint.
+ * @tparam TOutputValue type of the output point e.g., TSpace::Point
+ */
+ template
+struct CBDR {
+ ///Checking concepts
+ BOOST_CONCEPT_ASSERT(( concepts::CSpace ));
+ BOOST_STATIC_ASSERT(( TSpace::dimension == 2 ));
+ BOOST_STATIC_ASSERT(( TOutputValue::dimension == 2 ));
+ BOOST_STATIC_ASSERT(( TInputValue::dimension == 2 ));
+
+
+ typedef Reflection DigitizedReflection;
+ typedef std::vector>,GAVector>> BijectiveSearchTree;
+
+
+ /**
+ * CBDR Constructor.
+ * @param theta the angle given in radians.
+ * @param center the center of rotation.
+ * @param nbreflect the number of composition of bjijective reflections applied (2,4)
+ * @param km conditions the number of bijective reflection normal vectors
+ * @param policy either Linf, L2, Lcontinuity, see Policy
+ * @param precompute use precomputed table of the sorted bijective composition of bijective digitized reflections
+ * @param fast use the table that stores for each angle the composition that minimised the Linf metric distorsion
+ */
+ CBDR(const double theta,const typename TSpace::Point center,const size_t nbreflect,const size_t km, std::shared_ptr,CBDR_naiverotation>> policy,
+ const bool precompute=true, const bool fast =true ):nbijectiveGen(km),my_domain(typename TSpace::Point(0,0),typename TSpace::Point(100,100)),my_angle(theta),my_center(center),nbReflections(nbreflect),kmax(km),my_policy(policy),
+ usePrecomputedTable(precompute),useFastTable(fast),N(100),my_cbdr(initCBDRVec()) {
+ }
+
+
+
+
+ // whenever fast is set, load the table that stores for each angle the optimised set of bijective reflections
+ std::shared_ptr> initCBDRVec() {
+ if(useFastTable) {
+ std::cout << "use fast table"<> cbdrFastSolver(fastCBDRTable,my_angle,my_center,kmax);
+ return std::make_shared>(cbdrFastSolver.solve());
+ }else{
+ initcbdr_loadBijectiveReflectionSearchTree(my_domain);
+ /// solver part
+ CBDRSolver> cbdrSolver(my_angle, my_center,kmax,N);
+ return std::make_shared>(cbdrSolver.solve(my_domain,nbijectiveGen, vecBijectiveSearchTree,*my_policy));
+ }
+ }
+
+ void set_angle(const double newAngle) {
+ my_angle = newAngle;
+ if(fastTableFound) {
+ // does not use the policy since the table is already stored.
+ CBDRFastSolver> cbdrFastSolver(fastCBDRTable,my_angle,my_center,kmax);
+ my_cbdr=std::make_shared>(cbdrFastSolver.solve());
+ }else {
+ CBDRSolver> cbdrSolver(my_angle, my_center,kmax,N);
+ my_cbdr=std::make_shared>(cbdrSolver.solve(my_domain,nbijectiveGen, vecBijectiveSearchTree,*my_policy));
+ }
+ }
+
+ void setPolicy(const std::shared_ptr,CBDR_naiverotation>>& newPolicy) {
+ // need to search again for the precomputed tables
+ my_policy = newPolicy;
+ my_cbdr = initCBDRVec();
+ }
+
+
+ TOutputValue operator()( const TInputValue & aInput ) const
+ {
+ return my_cbdr->operator()(aInput-my_center)+my_center;
+ }
+
+ template
+ TImage rotateImage(TImage img) const {
+ return my_cbdr->rotateImage(img);
+ }
+
+ std::string tostring() const {
+ return {"CBDR"};
+ }
+
+ /// @return the centre of rotation
+ inline TOutputValue center() const{return my_center;}
+
+ public:
+ NBijectiveGenerator nbijectiveGen;
+ std::vector vecBijectiveSearchTree;
+ HyperRectDomain my_domain;
+
+ size_t kmax;
+ size_t nbReflections;
+ bool usePrecomputedTable;
+ bool useFastTable;
+ int N;/// number of sample rotation angle
+ double my_angle;
+ TOutputValue my_center;
+ std::shared_ptr,CBDR_naiverotation>> my_policy;
+ std::vector> fastCBDRTable;
+ std::shared_ptr> my_cbdr;
+
+
+ private:
+ bool fastTableFound;
+
+ /// load fast optimised table in case the table is found
+ std::vector> loadFastOptimisedTable(const std::string& fastTableName);
+ void initcbdr_loadBijectiveReflectionSearchTree(const HyperRectDomain& my_domain);
+ std::vector> initFastPrecomputationTable(
+ const HyperRectDomain& points,
+ NBijectiveGenerator& nbijectiveVectors,
+ std::vector& vecBijectiveSearchTree);
+ };
+
+
+ template
+ std::vector> CBDR::loadFastOptimisedTable(const std::string &fastTableName) {
+ auto vecTable = DGtal::functions::loadFastCBDRTable>(fastTableName);
+ std::vector> fastCBDRtab(vecTable.size());
+ std::transform(vecTable.begin(),vecTable.end(), fastCBDRtab.begin(), [](std::tuple>>,double,double>& x) {
+ return DGtal::CBDR_naiverotation>(std::get<0>(x));
+ });
+ return fastCBDRtab;
+ }
+
+ template
+ void CBDR::initcbdr_loadBijectiveReflectionSearchTree(const HyperRectDomain& my_domain) {
+
+ for(size_t nbReflec = 2 ; nbReflec <= nbReflections ; nbReflec+=2) {
+ BijectiveSearchTree vecBijNormals;
+ std::string tableName("CBDRTable_"+std::to_string(nbReflec)+"_"+std::to_string(kmax)+".txt");
+ // check that the file exists when usePrecomputedTable is set to true
+ if(usePrecomputedTable ) {
+ struct stat buffer;
+ bool exists = (stat (tableName.c_str(), &buffer) == 0);
+ if(!exists) {
+ usePrecomputedTable=false;
+ }
+ }
+
+ if(usePrecomputedTable){
+ std::cout << "loading precomputed table ..."<(tableName,nbReflec,kmax);
+ } else {
+ std::cout << "does not use precomputed table"<
+ std::vector> CBDR::
+ initFastPrecomputationTable(const HyperRectDomain &points, NBijectiveGenerator &nbijectiveVectors,
+ std::vector &vecBijectiveSearchTree) {
+ std::vector> fastCBDRtab;
+ auto Linf = std::make_shared,DGtal::CBDR_naiverotation>>();
+ auto LContinuity = std::make_shared,DGtal::CBDR_naiverotation>>();
+
+ std::string tablefilename ="CBDROptimisedTable_"+std::to_string(nbReflections)+"_"+std::to_string(kmax)+".txt";
+ std::ofstream file(tablefilename);
+ for(int alpha = 0 ; alpha < 91; ++alpha) {
+ double currentAngle = (alpha * M_PI) / 180.0;
+ CBDRSolver> cbdrsolv(currentAngle, my_center,kmax,N);
+ std::vector> bestReflections =cbdrsolv.solve(points,nbijectiveVectors,
+ vecBijectiveSearchTree,*my_policy);
+
+ CBDR_naiverotation currentRotation(bestReflections);
+ fastCBDRtab.push_back(currentRotation);
+
+ // display vector of reflections
+ for(size_t i = 0 ; ievaluate(points,currentRotation,currentAngle,my_center);
+ file << ";"<evaluate(points,currentRotation,currentAngle,my_center);
+ file << std::endl;
+
+ }
+ file.close();
+ std::cout << "table written !"<.
+ *
+ **/
+
+#pragma once
+
+/**
+* @file CBDRFastSolver.h
+ * @author S. Breuils, J.O. Lachaud, D. Coeurjolly
+ *
+ * @date 2024/08
+ *
+ * This file is part of the DGtal library.
+ */
+
+#if defined(CBDRFASTSOLVER_RECURSES)
+#error Recursive header files inclusion detected in CBDRFastSolver.h
+#else // defined(CBDRFASTSOLVER_RECURSES)
+/** Prevents recursive inclusion of headers. */
+#define CBDRFASTSOLVER_RECURSES
+
+#if !defined CBDRFASTSOLVER_h
+/** Prevents repeated inclusion of headers. */
+#define CBDRFASTSOLVER_h
+
+//////////////////////////////////////////////////////////////////////////////
+// Inclusions
+#include "NBijectiveReflectionGenerator.h"
+namespace DGtal {
+ template
+ class CBDRFastSolver {
+ typedef CBDR_naiverotation BijectiveReflections;
+ typedef Reflection DigitizedReflection;
+ public:
+ CBDRFastSolver(const std::vector& vecOptimisedReflections,const double rotAngle,
+ const typename TDomain::Point center,
+ const int km):my_vecOptimisedReflections(vecOptimisedReflections),my_angle(rotAngle),my_center(center),kmax(km){}
+
+ BijectiveReflections solve() {
+ double angleInDegrees = my_angle * 180.0 / M_PI;
+ angleInDegrees = fmod(angleInDegrees, 360.0);
+ if (angleInDegrees < 0) {
+ angleInDegrees += 360.0;
+ }
+ return my_vecOptimisedReflections[static_cast(std::round(angleInDegrees))];
+ }
+ protected:
+ std::vector my_vecOptimisedReflections;
+ double my_angle;
+ typename TDomain::Point my_center;
+ int kmax;
+ };
+}
+
+
+
+#endif //CBDRFASTSOLVER
+#undef CBDRFASTSOLVER_RECURSES
+#endif // else defined(CBDRFASTSOLVER_RECURSES)
\ No newline at end of file
diff --git a/src/DGtal/images/bijectiveRotations/CBDRSolver.h b/src/DGtal/images/bijectiveRotations/CBDRSolver.h
new file mode 100644
index 0000000000..6b2be68bdc
--- /dev/null
+++ b/src/DGtal/images/bijectiveRotations/CBDRSolver.h
@@ -0,0 +1,178 @@
+/**
+* This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as
+ * published by the Free Software Foundation, either version 3 of the
+ * License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see .
+ *
+ **/
+
+#pragma once
+
+/**
+* @file CBDRSolver.h
+ * @author S. Breuils, J.O. Lachaud, D. Coeurjolly
+ *
+ * @date 2024/08
+ *
+ * This file is part of the DGtal library.
+ */
+
+#if defined(CBDRSOLVER_RECURSES)
+#error Recursive header files inclusion detected in CBDRSolver.h
+#else // defined(CBDRSOLVER_RECURSES)
+/** Prevents recursive inclusion of headers. */
+#define CBDRSOLVER_RECURSES
+
+#if !defined CBDRSOLVER_h
+/** Prevents repeated inclusion of headers. */
+#define CBDRSOLVER_h
+
+//////////////////////////////////////////////////////////////////////////////
+// Inclusions
+#include
+#include "Policy.h"
+#include "NBijectiveReflectionGenerator.h"
+namespace DGtal{
+ /**
+ * Description of template struct 'CBDR Solver'
+ * \brief CBDR solver, use a policy to choose the composition of digitized reflections that minimises an error
+ * @tparam TSpace a 2 dimensional space.
+ * @tparam TDomain a 2 dimensional domain.
+ */
+ template
+ class CBDRSolver_GAvec{
+ public:
+ typedef CBDR_naiverotation BijectiveReflections;
+ typedef functors::ForwardRigidTransformation2D RealRotation;
+ typedef std::vector>,GAVector>> bijectiveReflect;
+ typedef ErrorVectorField ErrorRealVectors;
+ typedef Reflection DigitizedReflection;
+ typedef std::vector> VectorField;
+
+
+ CBDRSolver_GAvec(const size_t km, const double rotAngle, const typename TDomain::Point center):nBijectiveGenerator(km),my_angle(rotAngle),my_center(center){}
+
+ /// find the composition of bijective reflections that minimize either
+ /// Linf, L2 or Lcontinuity of error vector field given by the policy
+ /// @return {vector that minimizes an error given by policy, the error}
+ std::pair>,double> outputCompositionReflection(const TDomain& set2d,
+ typename bijectiveReflect::iterator& lowerAngle,
+ typename bijectiveReflect::iterator& upperAngle,
+ const Policy& policy){
+
+ // first error
+ auto it = lowerAngle;
+
+ typename bijectiveReflect::iterator itMinError;
+
+ std::vector bijectiveReflections;
+ std::vector> firstReflectionsIndex=(*lowerAngle).first;
+
+ for(size_t i = 0 ;i reflections(bijectiveReflections);
+
+ // finally get the average error
+ double minError = policy.evaluate(set2d, reflections,my_angle, my_center); //outputError(firsterrors,policy); // should now be policy
+
+ itMinError = lowerAngle;
+
+ for(it = lowerAngle+1 ; it != upperAngle ; ++it ){
+ std::vector currentbijectiveReflections;
+ std::vector> currentReflectionsIndex=(*it).first;
+ for(size_t i = 0 ;i currentreflections(currentbijectiveReflections);
+
+ // // finally get the error
+ double error = policy.evaluate(set2d, currentreflections,my_angle, my_center);
+ if(fabs(error) nBijectiveGenerator;
+ double my_angle;
+ typename TDomain::Point my_center;
+ };
+
+
+ /**
+ * Description of template struct 'CBDR Solver'
+ * \brief Aim: CBDR solver, use a policy to choose the composition of digitized reflections that minimises an error
+ * @tparam TSpace a 2 dimensional space.
+ * @tparam TDomain a 2 dimensional domain.
+ */
+ template
+ struct CBDRSolver {
+ typedef CBDR_naiverotation BijectiveReflections;
+ typedef functors::ForwardRigidTransformation2D RealRotation;
+ typedef std::vector>,GAVector>> bijectiveReflect;
+ typedef ErrorVectorField ErrorRealVectors;
+ typedef Reflection DigitizedReflection;
+ typedef std::vector> VectorField;
+ typedef std::vector>,GAVector>> BijectiveSearchTree;
+
+ CBDRSolver(const double rotAngle, const typename TDomain::Point center,const int km, const int NSamples):my_angle(rotAngle),my_center(center),kmax(km),N(NSamples){}
+
+ void setPolicy(const Policy customPolicy){my_policy(customPolicy);}
+
+ std::vector> solve(const TDomain& points,NBijectiveGenerator& nbijectiveVectors, std::vector& vecBijectiveSearchTree,const Policy& policy){
+ /// assume the sorted table was already generated either from a file or through one of the function of NbictiveGenerator
+ std::vector> bestParam;
+ double errorMin = points.myUpperBound[0];
+ for(size_t nReflection = 0 ; nReflection < vecBijectiveSearchTree.size();++nReflection){
+ typename std::vector>,GAVector>>::iterator lowerAngle;
+ typename std::vector>,GAVector>>::iterator upperAngle;
+ int numberOfCompositions=N;
+
+ nbijectiveVectors.getKNearestBijectiveComposition(
+ vecBijectiveSearchTree[nReflection],
+ lowerAngle,
+ upperAngle,
+ numberOfCompositions,
+ my_angle);
+
+
+ CBDRSolver_GAvec rotationSolver(kmax,my_angle,my_center);
+ std::pair>,double> bestParam_Error =rotationSolver.outputCompositionReflection(points,lowerAngle,upperAngle,policy);
+ if(bestParam_Error.second < errorMin) {
+ bestParam = bestParam_Error.first;
+ errorMin = bestParam_Error.second;
+ }
+ }
+ std::vector bestGAVectors;
+ for(GAVector indexB:bestParam){
+ bestGAVectors.push_back(DigitizedReflection(indexB));
+ }
+ return bestGAVectors;
+ }
+
+
+ protected:
+ int kmax;
+ int N;/// number of sample rotation angle
+ double my_angle;
+ typename TDomain::Point my_center;
+
+ };
+}
+
+
+#endif //CBDRSOLVER
+#undef CBDRSOLVER_RECURSES
+#endif // else defined(CBDRSOLVER_RECURSES)
+
diff --git a/src/DGtal/images/bijectiveRotations/CBDR_naiverotation.h b/src/DGtal/images/bijectiveRotations/CBDR_naiverotation.h
new file mode 100644
index 0000000000..ee9b01eb65
--- /dev/null
+++ b/src/DGtal/images/bijectiveRotations/CBDR_naiverotation.h
@@ -0,0 +1,104 @@
+/**
+* This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as
+ * published by the Free Software Foundation, either version 3 of the
+ * License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see .
+ *
+ **/
+
+#pragma once
+
+/**
+* @file CBDR_naiverotation.h
+ * @author S. Breuils, J.O. Lachaud, D. Coeurjolly
+ *
+ * @date 2024/08
+ *
+ * This file is part of the DGtal library.
+ */
+
+#if defined(CBDRNAIVEROTATION_RECURSES)
+#error Recursive header files inclusion detected in CBDR_naiverotation.h
+#else // defined(CBDRNAIVEROTATION_RECURSES)
+/** Prevents recursive inclusion of headers. */
+#define CBDRNAIVEROTATION_RECURSES
+
+#if !defined CBDRNAIVEROTATION_h
+/** Prevents repeated inclusion of headers. */
+#define CBDRNAIVEROTATION_h
+
+//////////////////////////////////////////////////////////////////////////////
+// Inclusions
+#include "DigitizedReflection.h"
+#include "GAVector.h"
+namespace DGtal {
+ /// vec since the parameters are the vectors of digitized reflections
+ template
+struct CBDR_naiverotation{
+ typedef Reflection DigitizedReflection;
+
+ std::vector bijectiveNormalVectors;
+
+ CBDR_naiverotation(const std::vector& bijectiveReflections={} ):bijectiveNormalVectors(bijectiveReflections){ }
+
+ explicit CBDR_naiverotation(const std::vector>& bijectiveGAvec ){
+ bijectiveNormalVectors.resize(bijectiveGAvec.size());
+ std::transform(bijectiveGAvec.begin(), bijectiveGAvec.end(), bijectiveNormalVectors.begin(),
+ [](const GAVector& p) { return DigitizedReflection(p); });
+ }
+
+ /**
+ * Operator
+ * @return apply the composition of bijective digitized reflection to the point.
+ */
+ TOutputValue operator()( const TInputValue & aInput ) const
+ {
+ if(bijectiveNormalVectors.size()<=0){
+ return (TOutputValue)aInput;
+ }
+ int i =0;
+ DigitizedReflection firstReflection(bijectiveNormalVectors[i]);
+ TOutputValue resultat = firstReflection(aInput);
+
+ i+=1;
+ for(; i < bijectiveNormalVectors.size(); ++i)
+ {
+ DigitizedReflection currentReflection(bijectiveNormalVectors[i]);
+ resultat = currentReflection(resultat);
+ }
+ return resultat;
+ }
+
+ template
+ TImage rotateImage(TImage img) const {
+ typedef typename TImage::Domain TDomain;
+ typedef DGtal::functors::DomainRigidTransformation2D < typename TImage::Domain, DGtal::CBDR_naiverotation> MyDomainTransformer;
+ typedef std::pair < typename TSpace::Point, typename TSpace::Point > Bounds;
+
+ MyDomainTransformer domainTransformer ( *this );
+ Bounds bounds = domainTransformer ( img.domain() );
+ TDomain transformedDomain ( bounds.first, bounds.second );
+ TImage rotatedImage ( transformedDomain );
+
+ for (typename TDomain::ConstIterator it = img.domain().begin(); it != img.domain().end(); ++it )
+ {
+ rotatedImage.setValue((*this)(*it),img(*it));
+ }
+ return rotatedImage;
+ }
+
+ };
+}
+
+
+#endif //CBDRNAIVEROTATION
+#undef CBDRNAIVEROTATION_RECURSES
+#endif // else defined(CBDRNAIVEROTATION_RECURSES)
diff --git a/src/DGtal/images/bijectiveRotations/CDLR.h b/src/DGtal/images/bijectiveRotations/CDLR.h
new file mode 100644
index 0000000000..c9950c4112
--- /dev/null
+++ b/src/DGtal/images/bijectiveRotations/CDLR.h
@@ -0,0 +1,203 @@
+/**
+* This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as
+ * published by the Free Software Foundation, either version 3 of the
+ * License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see .
+ *
+ **/
+
+#pragma once
+
+/**
+* @file CDLR.h
+ * @author S. Breuils, J.O. Lachaud, D. Coeurjolly
+ *
+ * @date 2024/07/9
+ *
+ * This file is part of the DGtal library.
+ */
+
+#if defined(CDLR_RECURSES)
+#error Recursive header files inclusion detected in CDLR.h
+#else // defined(CDLR_RECURSES)
+/** Prevents recursive inclusion of headers. */
+#define CDLR_RECURSES
+
+#if !defined CDLR_h
+/** Prevents repeated inclusion of headers. */
+#define CDLR_h
+
+#include "DGtal/base/Common.h"
+#include
+#include "CDLR_naiverotation.h"
+#include "DigitizedReflection.h"
+#include "Policy.h"
+namespace DGtal {
+
+ /**
+ * Description of template struct CDLR
+ * \brief CDLR : Rotation with Discrete Line Reflections,
+ * @tparam TSpace a 2 dimensional space.
+ * @tparam TInputValue type of the input point e.g., TSpace::RealPoint.
+ * @tparam TOutputValue type of the output point e.g., TSpace::Point
+ */
+ template < typename TSpace, typename TInputValue = typename TSpace::Point, typename TOutputValue = typename TSpace::Point>
+ class CDLR {
+ public:
+ /**
+ * CDLR Constructor.
+ * @param ang the angle given in radians.
+ * @param ptCenter the center of rotation.
+ * @param policy either Linf, L2, Lcontinuity, see Policy
+ */
+ CDLR( double ang, TOutputValue ptCenter, std::shared_ptr,CDLR_naiverotation>> policy );
+
+ /// @return the angle of rotation.
+ inline double angle() const{return my_angle;};
+ /// @return set the angle of rotation and call the composition of reflections solver.
+ inline void set_angle(const double new_angle) {
+ my_angle=new_angle;
+ my_naive_rdlr_rotation.set_angle(new_angle);
+ dslSolver(new_angle,my_center);
+ }
+
+ /// initialisation function to find the composition of Discrete Line Reflection that minimises the sum of
+ /// linf and lcontinuity
+ void dslSolver(double ang, TOutputValue ptCenter);
+
+
+
+ /// @return the centre of rotation
+ inline TOutputValue center() const{return my_center;};
+ /// @return a reference to the centre of rotation
+ inline TOutputValue& center(){return my_center;};
+
+ /// @param p a lattice point
+ /// @return the rotation of the point \a p according to the current
+ /// angle and center.
+ TOutputValue rotate( const TInputValue& p ) const;
+ TOutputValue operator()( const TInputValue& p ) const;
+
+ template
+ Img rotateImage(Img img) const;
+
+ std::string tostring() const {
+ return {"CDLR"};
+ }
+
+ protected:
+ /// The angle of rotation
+ double my_angle;
+ /// The center of rotation
+ TOutputValue my_center;
+
+
+
+ CDLR_naiverotation my_naive_rdlr_rotation;
+ std::shared_ptr,CDLR_naiverotation>> my_policy;
+
+ private:
+ /// DSL specific variables, see Andres paper
+ double a;
+ double b;
+ double my_minAngle;
+ };
+
+ template
+ CDLR::CDLR( const double ang, const TOutputValue ptCenter, std::shared_ptr,CDLR_naiverotation>> policy ):my_angle(ang),my_center(ptCenter),my_naive_rdlr_rotation(ang,ptCenter,0.),my_policy(policy){
+ dslSolver(ang,ptCenter);
+ }
+
+ template
+ void CDLR::dslSolver(double ang, TOutputValue ptCenter) {
+
+ a = std::sin(ang/2.0);
+ b = std::cos(ang/2.0);
+
+ TInputValue origin(0,0);
+ DGtal::functors::ForwardRigidTransformation2D::RealPoint,DGtal::functors::Identity> euclideanRot(my_center,my_angle,origin);
+
+ double errorMinAlpha = 200.*200;
+ int N = 200;
+
+ // create the domain
+ TOutputValue A(0,0);
+ TOutputValue B(1000,1000);
+ HyperRectDomain my_domain(A,B);
+
+
+ /// look for the pair of reflections that minimizes the Linf and Lcontinuity norm
+ std::vector mix_errors;
+ std::vector angles;
+ for(int k=0;k linf_per_image;
+ my_naive_rdlr_rotation.set_startingAngle(alphaCourant);
+ double error = my_policy->evaluate(my_domain,my_naive_rdlr_rotation,my_angle,my_center);
+ mix_errors.push_back(error);
+ }
+ for(int idxImages = 0 ; idxImages < mix_errors.size() ; idxImages++) {
+ if(mix_errors[idxImages] < errorMinAlpha) {
+ errorMinAlpha = mix_errors[idxImages];
+ my_minAngle = angles[idxImages];
+ }
+ }
+ }
+
+
+ template
+ TOutputValue CDLR::rotate( const TInputValue& p ) const{
+ return this->operator()(p);
+ }
+
+
+
+
+
+
+
+ template
+ TOutputValue CDLR::operator()( const TInputValue& p ) const{
+ return my_naive_rdlr_rotation.reflect(my_minAngle+(my_angle),my_center,my_naive_rdlr_rotation.reflect(my_minAngle,my_center,p-my_center))+my_center;
+ }
+
+
+ template
+ template
+ TImage CDLR::rotateImage( const TImage img ) const{
+ typedef typename TImage::Domain TDomain;
+ typedef DGtal::functors::DomainRigidTransformation2D < typename TImage::Domain, DGtal::CDLR> MyDomainTransformer;
+ typedef std::pair < typename TSpace::Point, typename TSpace::Point > Bounds;
+
+ typename TSpace::Point bottomLeft(-2,-2);
+ typename TSpace::Point topRight(2,2);
+ MyDomainTransformer domainTransformer ( *this );
+ Bounds bounds = domainTransformer ( img.domain() );
+ TDomain transformedDomain ( bounds.first+bottomLeft, bounds.second+topRight );
+ TImage rotatedImage ( transformedDomain );
+
+ for (typename TDomain::ConstIterator it = img.domain().begin(); it != img.domain().end(); ++it )
+ {
+ rotatedImage.setValue((*this)(*it),img(*it));
+ }
+ return rotatedImage;
+ }
+
+
+}
+
+
+#endif //CDLR
+#undef CDLR_RECURSES
+#endif // else defined(CDLR_RECURSES)
diff --git a/src/DGtal/images/bijectiveRotations/CDLR_naiverotation.h b/src/DGtal/images/bijectiveRotations/CDLR_naiverotation.h
new file mode 100644
index 0000000000..8aa66bc7f6
--- /dev/null
+++ b/src/DGtal/images/bijectiveRotations/CDLR_naiverotation.h
@@ -0,0 +1,139 @@
+
+/**
+* This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as
+ * published by the Free Software Foundation, either version 3 of the
+ * License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see .
+ *
+ **/
+
+#pragma once
+
+/**
+* @file CDLR_naiverotation.h
+ * @author S. Breuils, J.O. Lachaud, D. Coeurjolly
+ *
+ * @date 2024/07
+ *
+ * This file is part of the DGtal library.
+ */
+
+#if defined(CDLR_NAIVEROTATION_RECURSES)
+#error Recursive header files inclusion detected in RDSL.h
+#else // defined(CDLR_NAIVEROTATION_RECURSES)
+/** Prevents recursive inclusion of headers. */
+#define CDLR_NAIVEROTATION_RECURSES
+
+#if !defined CDLR_NAIVEROTATION_h
+/** Prevents repeated inclusion of headers. */
+#define CDLR_NAIVEROTATION_h
+
+#include "DGtal/base/Common.h"
+#include
+#include "PointUtils.h"
+#include "DigitizedReflection.h"
+
+namespace DGtal {
+
+
+ template < typename TSpace, typename TInputValue = typename TSpace::Point, typename TOutputValue = typename TSpace::Point>
+ struct CDLR_naiverotation {
+
+ inline int X(int y, const double a, const double b, const int k) const {
+ return std::ceil((2.0*k-1)/2.0 - (a/b)*(y)); // a = sin(theta) ; b=cos(theta)
+ }
+
+ TOutputValue reflect(const double angle, TOutputValue center, const TInputValue& p ) const{
+ double a = std::sin(angle/2.0);
+ double b = std::cos(angle/2.0);
+ TOutputValue pcentered = p;
+
+ int k = std::floor(pcentered[0] + (a/b)*pcentered[1] + 0.5);
+
+ TOutputValue X1 = TOutputValue(X(std::ceil(a*b*k),a,b,k),std::ceil(a*b*k));
+ TOutputValue X2 = TOutputValue(X(std::floor(a*b*k),a,b,k),std::floor(a*b*k));
+
+ const double line2 = a*X1[0]-b*X1[1];
+ if(line2(-b/2)){
+ return TOutputValue(X(2*X1[1]-pcentered[1],a,b,k),2*X1[1]-pcentered[1]);
+ }else{
+ const double line3 = a*X2[0]-b*X2[1];
+ if(line3(-b/2)){
+ return TOutputValue(X(2*X2[1]-pcentered[1],a,b,k),2*X2[1]-pcentered[1]);
+ }
+ else{
+ return TOutputValue(X(X1[1]+X2[1]-pcentered[1],a,b,k),X1[1]+X2[1]-pcentered[1]);
+ }
+ }
+
+
+ }
+
+
+
+
+ CDLR_naiverotation( double ang=0., TOutputValue ptCenter=TOutputValue(0,0), double starting_angle =0. ):my_angle(ang),my_center(ptCenter),my_startingAngle(starting_angle) {
+ }
+
+ /// @return the angle of rotation.
+ inline double angle() const{return my_angle;};
+
+ /// @return the starting angle.
+ inline double startingAngle() const{return my_startingAngle;};
+
+ ///set the angle of rotation and call the composition of reflections solver.
+ inline void set_angle(const double new_angle){my_angle=new_angle;}
+
+ /// set the starting angle not the rotation angle
+ inline void set_startingAngle(const double new_startingAngle){my_startingAngle=new_startingAngle;}
+
+
+
+
+ /// @return the centre of rotation
+ inline TOutputValue center() const{return my_center;};
+ /// @return a reference to the centre of rotation
+ inline TOutputValue& center(){return my_center;};
+
+ /// @param p a lattice point
+ /// @return the rotation of the point \a p according to the current
+ /// angle and center.
+ TOutputValue rotate( const TInputValue& p ) const {
+ return reflect(my_startingAngle + my_angle, my_center,
+ reflect(
+ my_startingAngle, my_center, p));
+ }
+ TOutputValue operator()( const TInputValue& p ) const {
+ return rotate(p);
+ }
+
+
+
+
+ protected:
+ /// The angle of rotation.
+ double my_angle;
+ double my_startingAngle;
+
+ /// The center of rotation.
+ TOutputValue my_center;
+
+
+ };
+
+}
+
+
+
+#endif //CDLR_NAIVEROTATION
+#undef CDLR_NAIVEROTATION_RECURSES
+#endif // else defined(CDLR_NAIVEROTATION_RECURSES)
+
diff --git a/src/DGtal/images/bijectiveRotations/DigitizedReflection.h b/src/DGtal/images/bijectiveRotations/DigitizedReflection.h
new file mode 100644
index 0000000000..af121fdc96
--- /dev/null
+++ b/src/DGtal/images/bijectiveRotations/DigitizedReflection.h
@@ -0,0 +1,70 @@
+/**
+* This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as
+ * published by the Free Software Foundation, either version 3 of the
+ * License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see .
+ *
+ **/
+
+#pragma once
+
+/**
+ * @file DigitizedReflection.h
+ * @author S. Breuils, J.O. Lachaud, D. Coeurjolly
+ *
+ * @date 2024/07/9
+ *
+ * This file is part of the DGtal library.
+ */
+
+#if defined(DigitizedReflection_RECURSES)
+#error Recursive header files inclusion detected in DigitizedReflection.h
+#else // defined(DigitizedReflection_RECURSES)
+/** Prevents recursive inclusion of headers. */
+#define DigitizedReflection_RECURSES
+
+#if !defined DigitizedReflection_h
+/** Prevents repeated inclusion of headers. */
+#define DigitizedReflection_h
+
+//////////////////////////////////////////////////////////////////////////////
+// Inclusions
+#include "GAVector.h"
+
+namespace DGtal {
+ template
+ struct Reflection
+ {
+ GAVector normalVector;
+
+ explicit Reflection ( const GAVector & m=GAVector())
+ :normalVector(m){}
+
+ /**
+ * Operator
+ * @return the reflected and digitized point.
+ */
+ inline
+ TOutputValue operator()( const TInputValue & aInput ) const
+ {
+ DGtal::functors::VectorRounding < TInputValue, TOutputValue > roundingOpe;
+ Z2i::RealPoint m_r = Z2i::RealPoint(normalVector.my_gavec[0],normalVector.my_gavec[1]);
+ Z2i::RealPoint x_r = Z2i::RealPoint(aInput[0],aInput[1]);
+ Z2i::RealPoint p=x_r - 2.0*((x_r[0]*m_r[0] + x_r[1]*m_r[1])/(m_r[0]*m_r[0] + m_r[1]*m_r[1]))*m_r;
+ return roundingOpe(p);
+ }
+
+ };
+}
+#endif //DigitizedReflection
+
+#undef DigitizedReflection_RECURSES
+#endif // else defined(DigitizedReflection_RECURSES)
\ No newline at end of file
diff --git a/src/DGtal/images/bijectiveRotations/ErrorBijectiveRotation.h b/src/DGtal/images/bijectiveRotations/ErrorBijectiveRotation.h
new file mode 100644
index 0000000000..828479a826
--- /dev/null
+++ b/src/DGtal/images/bijectiveRotations/ErrorBijectiveRotation.h
@@ -0,0 +1,102 @@
+/**
+* This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as
+ * published by the Free Software Foundation, either version 3 of the
+ * License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see .
+ *
+ **/
+
+#pragma once
+
+/**
+* @file ErrorBijectiveRotation.h
+ * @author S. Breuils, J.O. Lachaud, D. Coeurjolly
+ *
+ * @date 2024/07/9
+ *
+ * This file is part of the DGtal library.
+ */
+
+#if defined(ErrorBijectiveRotation_RECURSES)
+#error Recursive header files inclusion detected in ErrorBijectiveRotation.h
+#else // defined(ErrorBijectiveRotation_RECURSES)
+/** Prevents recursive inclusion of headers. */
+#define ErrorBijectiveRotation_RECURSES
+
+#if !defined ErrorBijectiveRotation_h
+/** Prevents repeated inclusion of headers. */
+#define ErrorBijectiveRotation_h
+
+//////////////////////////////////////////////////////////////////////////////
+// Inclusions
+#include
+
+namespace DGtal{
+ template
+ struct ErrorVectorField{
+
+ typedef functors::ForwardRigidTransformation2D RealRotation;
+ typedef std::vector> VectorField;
+
+ protected:
+ TBijectiveReflections normalVectors;
+ typename TSpace::Point my_center;
+ RealRotation targetRotation;
+ RealRotation originCenteredRotation;
+
+ public:
+ ErrorVectorField(const TBijectiveReflections& reflections,const double theta,const typename TSpace::Point center):normalVectors(reflections),targetRotation(center,theta,{0,0}),my_center(center),originCenteredRotation({0,0},theta,{0,0})
+ {}
+
+ /// compute vector of errors for each pixel between either
+ /// - each digitized reflections and the real reflection
+ /// - or the the composition of digitized reflections and the final rotation
+ VectorField getOutputVectorFieldFromContour(const TDomain& set2dContour, bool continuityVecField =false ){
+ VectorField outVecField;
+
+ for (typename TDomain::ConstIterator it = set2dContour.begin(); it != set2dContour.end(); ++it ) {
+ std::vector pixelError;
+ typename TSpace::Point p = *it;
+
+ typename TSpace::Point preflections = normalVectors(p-my_center)+my_center;
+ typename TSpace::RealPoint protation = targetRotation(p);
+
+ // compute the error field
+ typename TSpace::RealPoint error = protation - preflections ;
+ pixelError.push_back(error);
+
+ // compute the eucliean rotation of the neighbors of p
+ if(continuityVecField) {
+ // for the 8-Neighbor, compute the rotation
+ for(int veci = -1 ; veci <2 ; ++veci) {
+ for(int vecj = -1 ; vecj<2 ; ++vecj) {
+ if(veci!=0 || vecj!=0){
+ typename TSpace::RealPoint vecij_rot = originCenteredRotation({static_cast(veci),static_cast(vecj)});
+ typename TSpace::RealPoint neigh_rot = protation+vecij_rot;
+ pixelError.push_back(neigh_rot-preflections);
+ }
+
+ }
+ }
+
+ }
+
+
+ outVecField.push_back(pixelError);
+ }
+ return outVecField;
+ }
+ };
+}
+
+#endif //ErrorBijectiveRotation
+#undef ErrorBijectiveRotation_RECURSES
+#endif // else defined(ErrorBijectiveRotation_RECURSES)
diff --git a/src/DGtal/images/bijectiveRotations/ErrorVectorField.h b/src/DGtal/images/bijectiveRotations/ErrorVectorField.h
new file mode 100644
index 0000000000..7ea4228f40
--- /dev/null
+++ b/src/DGtal/images/bijectiveRotations/ErrorVectorField.h
@@ -0,0 +1,110 @@
+/**
+* This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as
+ * published by the Free Software Foundation, either version 3 of the
+ * License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see .
+ *
+ **/
+
+#pragma once
+
+/**
+* @file ErrorVectorField.h
+ * @author S. Breuils, J.O. Lachaud, D. Coeurjolly
+ *
+ * @date 2024/08
+ *
+ * This file is part of the DGtal library.
+ */
+
+#if defined(ERRORVECTORFIELD_RECURSES)
+#error Recursive header files inclusion detected in ErrorVectorField.h
+#else // defined(ERRORVECTORFIELD_RECURSES)
+/** Prevents recursive inclusion of headers. */
+#define ERRORVECTORFIELD_RECURSES
+
+#if !defined ERRORVECTORFIELD_h
+/** Prevents repeated inclusion of headers. */
+#define ERRORVECTORFIELD_h
+
+//////////////////////////////////////////////////////////////////////////////
+// Inclusions
+#include
+
+namespace DGtal{
+ template
+ struct ErrorVectorField{
+
+ typedef functors::ForwardRigidTransformation2D RealRotation;
+ typedef std::vector> VectorField;
+
+ protected:
+ TBijectiveReflections normalVectors;
+ typename TSpace::Point my_center;
+ RealRotation targetRotation;
+ RealRotation originCenteredRotation;
+
+ public:
+ ErrorVectorField(const TBijectiveReflections& reflections,const double theta,const typename TSpace::Point center):normalVectors(reflections),targetRotation(center,theta,{0,0}),my_center(center),originCenteredRotation({0,0},theta,{0,0})
+ {}
+
+ /// compute vector of errors for each pixel between either
+ /// - each digitized reflections and the real reflection
+ /// - or the the composition of digitized reflections and the final rotation
+ VectorField getOutputVectorFieldFromContour(const TDomain& set2dContour, bool continuityVecField =false ){
+ VectorField outVecField;
+ // std::cout << "DSL output vector error : normalVectors "<
+#include
+
+namespace DGtal {
+ template
+ struct GAVector{
+ TInputValue my_gavec;
+
+ /// @brief zero multivector
+ explicit GAVector( const TInputValue& pt= TInputValue(0,0))
+ : my_gavec( pt ) {}
+
+ /// @brief scalar part of the a mv, i.e. compute the dot product part of the geometric product
+ typename TInputValue::Component dot( GAVector other ) const
+ {
+ return my_gavec[0] * other.my_gavec[0] + my_gavec[1] * other.my_gavec[1];
+ }
+
+ /// @brief bivector part of the a mv, i.e. compute the outer product part of the geometric product
+ typename TInputValue::Component bivectorPart( const GAVector& other ) const
+ {
+ return my_gavec[0] * other.my_gavec[1] - my_gavec[1] * other.my_gavec[0];
+ }
+
+ /// @brief geometric product with a scalar
+ GAVector operator*( typename TInputValue::Component f ) const
+ {
+ return GAVector( TInputValue(my_gavec[0] * f, my_gavec[1] * f) );
+ }
+
+ /// @brief geometric product between two vectors
+ GAVector operator*(const GAVector& v )// geometric product
+ {
+ return GAVector( TInputValue(this->dot(v), this->bivectorPart(v)) );
+ }
+ bool operator<( const GAVector& other ) const
+ {
+ return bivectorPart( other ) > 0;
+ }
+
+ double angleToXAxis() const
+ {
+ return atan2( my_gavec[1], my_gavec[0] );
+ }
+
+ bool operator==( const GAVector& other ) const
+ {
+ return bivectorPart( other ) == 0;
+ }
+
+
+ };
+}
+
+#endif //GAVector
+
+#undef GAVector_RECURSES
+#endif // else defined(GAVector_RECURSES)
diff --git a/src/DGtal/images/bijectiveRotations/NBijectiveReflectionGenerator.h b/src/DGtal/images/bijectiveRotations/NBijectiveReflectionGenerator.h
new file mode 100644
index 0000000000..3eed98e44e
--- /dev/null
+++ b/src/DGtal/images/bijectiveRotations/NBijectiveReflectionGenerator.h
@@ -0,0 +1,290 @@
+/**
+* This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as
+ * published by the Free Software Foundation, either version 3 of the
+ * License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see .
+ *
+ **/
+
+#pragma once
+
+/**
+* @file NBijectiveReflectionGenerator.h
+ * @author S. Breuils, J.O. Lachaud, D. Coeurjolly
+ *
+ * @date 2024/08
+ *
+ * This file is part of the DGtal library.
+ */
+
+#if defined(NBIJECTIVEREFLECTIONGENERATOR_RECURSES)
+#error Recursive header files inclusion detected in NBijectiveReflectionGenerator.h
+#else // defined(CBDRFASTSOLVER_RECURSES)
+/** Prevents recursive inclusion of headers. */
+#define NBIJECTIVEREFLECTIONGENERATOR_RECURSES
+
+#if !defined NBIJECTIVEREFLECTIONGENERATOR_h
+/** Prevents repeated inclusion of headers. */
+#define NBIJECTIVEREFLECTIONGENERATOR_h
+
+//////////////////////////////////////////////////////////////////////////////
+// Inclusions
+#include
+#include
+#include
+#include
+#include
+#include
+#include
+#include
+#include
+#include
+#include "GAVector.h"
+#include "DigitizedReflection.h"
+#include "CBDR_naiverotation.h"
+
+namespace DGtal{
+ template
+ struct NBijectiveGenerator{
+
+ explicit NBijectiveGenerator(const size_t km):kmax(km){
+ for(int k=0; k({x_kpp_k,y_kpp_k}));
+ BijectiveVectors.push_back(GAVector({x_k_kpp,y_k_kpp}));
+ BijectiveVectors.push_back(GAVector({x_1_2kpp,y_1_2kpp}));
+ BijectiveVectors.push_back(GAVector({x_2kpp_1,y_2kpp_1}));
+ BijectiveVectors.push_back(GAVector({x_kpp_k_m,y_kpp_k_m}));
+ BijectiveVectors.push_back(GAVector({x_k_kpp_m,y_k_kpp_m}));
+ BijectiveVectors.push_back(GAVector({x_1_2kpp_m,y_1_2kpp_m}));
+ BijectiveVectors.push_back(GAVector({x_2kpp_1_m,y_2kpp_1_m}));
+ }
+ }
+
+
+
+ /// @brief compose bijective digitized reflections
+ /// @param Avector the composition of bijective reflection
+ /// @return vector of indices of bijective digitized reflections
+ std::vector,GAVector>> composeBijectiveReflections(
+ std::vector,GAVector>>& Avector ){
+ std::vector,GAVector>> result;
+ result.reserve(Avector.size()*BijectiveVectors.size());
+
+ for(std::pair,GAVector> normalVectorsAndResult : Avector){
+ std::pair,GAVector> resultatCourant = normalVectorsAndResult;
+ for(size_t indexB = 0 ; indexB outVec = normalVectorsAndResult.second*BijectiveVectors[indexB];
+ if(outVec.my_gavec[0] <0){
+ // resultatCourant.first.insert(resultatCourant.first.begin(),indexB);
+ // resultatCourant.second = outVec*-1;
+ }else{
+ resultatCourant.first.push_back(indexB);
+ resultatCourant.second = outVec;
+ result.push_back(resultatCourant);
+ resultatCourant=normalVectorsAndResult;
+ }
+
+ }
+ }
+ // unique
+ std::sort(result.begin(), result.end(), [](const std::pair,GAVector>& b1, const std::pair,GAVector>& b2) {
+ return b1.second < b2.second;});
+
+
+ return result;
+ }
+
+ /// @brief predicate to compare two composition of bijective reflections
+ /// @param b1 first pair of composition of bijective reflections
+ /// @param b2 second pair of composition of bijective reflections
+ /// @return composition leads to same error?
+ bool sameBijectiveComposition(const std::pair,GAVector>& b1, const std::pair,GAVector>& b2) const{
+ // both reflections thanks to reflection composer
+ if(fabs(b1.second.angleToXAxis()-b1.second.angleToXAxis())>1e-4){
+ return false;
+ }
+ std::vector> normals1 = getGAVectorFromIndex(b1.first);
+ std::vector> normals2 = getGAVectorFromIndex(b2.first);
+ CBDR_naiverotation reflectionsCompose1(normals1);
+ CBDR_naiverotation reflectionsCompose2(normals2);
+
+ typename TSpace::Point zeroPt(0,0);
+
+ for(int i =0 ; i<100;++i){
+ for(int j =0 ; j<100;++j){
+ double x=i-50;
+ double y=j-50;
+ typename TSpace::RealPoint pointCourant = {x,y};
+ typename TSpace::Point out1 = reflectionsCompose1(pointCourant);
+ typename TSpace::Point out2 = reflectionsCompose2(pointCourant);
+ if((out1-out2)!= zeroPt ){
+ return false;
+ }
+ }
+ }
+ return true;
+ }
+
+
+ std::vector>,GAVector>> vecBijNormals_index_2_GAVector(const std::vector,GAVector>>& vecBijNormalsIndex) {
+ std::vector>,GAVector>> vecBijNormalsIndex_GAVector;
+ for(std::pair,GAVector> pairIndexGAVector : vecBijNormalsIndex) {
+ std::vector> resGavec;
+
+ for(int index : pairIndexGAVector.first) {
+ resGavec.push_back(BijectiveVectors[index]);
+ }
+ vecBijNormalsIndex_GAVector.push_back({resGavec,pairIndexGAVector.second});
+ }
+ return vecBijNormalsIndex_GAVector;
+
+ }
+
+ /// @brief generate all n =2,4 reflection composition
+ /// @param n number of reflection to compose
+ std::vector,GAVector>> n_bijectiveReflections_get_NormalVectorsAngles(size_t n)
+ {
+ std::vector,GAVector>> vecBijNormals;
+ for(int i =0 ; i< BijectiveVectors.size() ; ++i){
+ vecBijNormals.push_back({{i},BijectiveVectors[i]});
+ }
+ if(n>3) {
+ --n;
+ }
+
+ for ( size_t j = 1; j < n; ++j )
+ {
+ vecBijNormals = composeBijectiveReflections( vecBijNormals);
+
+ // remove duplicates with predicate defined below
+ auto last = std::unique(vecBijNormals.begin(),vecBijNormals.end(),[this](const std::pair,GAVector>& b1, const std::pair,GAVector>& b2) {
+ return sameBijectiveComposition(b1, b2);});
+ vecBijNormals.erase(last,vecBijNormals.end());
+ }
+ if(n==3){
+ /// interpret as the composition of 4 vectors : 3 bijective and 1 trivial bijective reflection
+ for(auto& pairIndicesGAVec : vecBijNormals)
+ pairIndicesGAVec.first.insert(pairIndicesGAVec.first.begin(),0);
+ }
+
+ std::sort(vecBijNormals.begin(), vecBijNormals.end(), [](const std::pair,GAVector>& b1, const std::pair,GAVector>& b2) {
+ return b1.second.angleToXAxis() < b2.second.angleToXAxis();});
+
+
+ return vecBijNormals;
+ }
+
+ void displayNormalVectorsAndAngles(const std::vector,GAVector>>& vecBijNormals){
+ for(size_t i =0;i finalNormalVec = vecBijNormals[i].second;
+
+ if(2.0*finalNormalVec.angleToXAxis() < M_PI_4/8){
+ std::cout << "m=[";
+ for(int indexBijective : currentIndicesNormalVector){
+ std::cout <<"("<>,GAVector