begin process at 2012 05 27 15:24:05
  Trouver un code source :
 
dans
 
Accueil > 

Code

 > 

Maths & Algorithmes

 > DAUBCHIES D4 WAVELET : TRANSFORM ET INVERSE TRANSFORM

DAUBCHIES D4 WAVELET : TRANSFORM ET INVERSE TRANSFORM


 Information sur la source

Note :
Aucune note
Catégorie :Maths & Algorithmes Niveau :Initié Date de création :24/09/2004 Date de mise à jour :26/09/2004 12:27:58 Vu :4 748

Auteur : hitcher

Ecrire un message privé
Commentaire sur cette source (5)
Ajouter un commentaire et/ou une note

 Description

Cette class permet d'effectuer la transformation (et l'inverse) en utilisant l'ondelette de daubechies a 4 moments.


Remarque: pour traiter l'effet bord du a ce type d'ondelette, j'ai utiliser la technique de "wrap around" tel que T+1=T-1 et donc T+2=T-2, au lieu du T+1=0 et T+2=1

Source

  • #ifndef DAUB4_H
  • #define DAUB4_H
  • #include <cmath>
  • static const double c0=0.6830127;
  • static const double c1=1.1830127;
  • static const double c2=0.3169873;
  • static const double c3=-0.1830127;
  • static const double b0=c3;
  • static const double b1=-c2;
  • static const double b2=c1;
  • static const double b3=-c0;
  • static const double sqr= pow(2,-0.5);
  • class Daub
  • {
  • private:
  • int T;
  • public:
  • Daub(int n, double data[]):T(n){}
  • static void InvTrans(double* data, int n);
  • static void Transform(double* data, int n);
  • ~Daub(){}
  • };
  • void Daub::InvTrans(double* data, int n)
  • {
  • int i,j;
  • double* dt=new double[2*n];
  • for(i=0; i<2*n; i++){dt[i]=0;}
  • for(i=0, j=0; i<2*n; i+=2, j+=1){dt[i]=data[j];}
  • for(i=0, j=0; i<n; i+=2, j+=1)
  • {
  • if(i<2)
  • {
  • data[0]=sqr*(c0*dt[0]+c2*dt[1]+b0*dt[n]+b2*dt[n+1]);
  • data[1]=sqr*(c1*dt[0]+c3*dt[1]+b1*dt[n]+b3*dt[n+1]);
  • }
  • else
  • {
  • data[i] =sqr*(c0*dt[i]+c2*dt[i-2]+b0*dt[n+i]+b2*dt[n+i-2]);
  • data[i+1]=sqr*(c1*dt[i]+c3*dt[i-2]+b1*dt[n+i]+b3*dt[n+i-2]);
  • }
  • }
  • delete [] dt;
  • }
  • void Daub::Transform(double* data, int n)
  • {
  • int i,g, wrp1, wrp2;
  • double* dt= new double[n];
  • if(n!=2){wrp1=n-2;wrp2=n-3;}
  • else{wrp1=0;wrp2=1;}
  • for(i=0, g=0 ; i<n ;i+=2, g++)
  • {
  • if(g==n/2-1)
  • {
  • dt[g]=sqr*(c0*data[i]+c1*data[i+1]+c2*data[wrp1]+c3*data[wrp2]);
  • dt[g+n/2]=sqr*(b0*data[i]+b1*data[i+1]+b2*data[wrp1]+b3*data[wrp2]);
  • }
  • else
  • {
  • dt[g]=sqr*(c0*data[i]+c1*data[i+1]+c2*data[i+2]+c3*data[i+3]);
  • dt[g+n/2]=sqr*(b0*data[i]+b1*data[i+1]+b2*data[i+2]+b3*data[i+3]);
  • }
  • }
  • for(i=0; i<n ;i++){data[i]=dt[i];}
  • delete [] dt;
  • }
  • #endif

#ifndef DAUB4_H
#define DAUB4_H

#include <cmath>

	static const double c0=0.6830127;
	static const double c1=1.1830127;
	static const double c2=0.3169873;
	static const double c3=-0.1830127;
	static const double b0=c3;
	static const double b1=-c2;
	static const double b2=c1;
	static const double b3=-c0;
	static const double sqr= pow(2,-0.5);
 
class Daub
{
	private:
	int T;

	public:
	
	Daub(int n, double data[]):T(n){}
	static void InvTrans(double* data, int n);
	static void Transform(double* data, int n);
	~Daub(){}
};

void Daub::InvTrans(double* data, int n)
{

	int i,j;
	double* dt=new double[2*n];

	for(i=0; i<2*n; i++){dt[i]=0;}

	for(i=0, j=0; i<2*n; i+=2, j+=1){dt[i]=data[j];}

	for(i=0, j=0; i<n; i+=2, j+=1)
	{
		if(i<2)
		{
		data[0]=sqr*(c0*dt[0]+c2*dt[1]+b0*dt[n]+b2*dt[n+1]);
		data[1]=sqr*(c1*dt[0]+c3*dt[1]+b1*dt[n]+b3*dt[n+1]);
		}
		else
		{
		data[i]  =sqr*(c0*dt[i]+c2*dt[i-2]+b0*dt[n+i]+b2*dt[n+i-2]);
		data[i+1]=sqr*(c1*dt[i]+c3*dt[i-2]+b1*dt[n+i]+b3*dt[n+i-2]);
		}
	}
delete [] dt;
}

void Daub::Transform(double* data, int n)
{
	int i,g, wrp1, wrp2;
	double* dt= new double[n];

	if(n!=2){wrp1=n-2;wrp2=n-3;}
	else{wrp1=0;wrp2=1;}

	for(i=0, g=0 ; i<n ;i+=2, g++)
	{
		if(g==n/2-1)
		{
	dt[g]=sqr*(c0*data[i]+c1*data[i+1]+c2*data[wrp1]+c3*data[wrp2]);
	dt[g+n/2]=sqr*(b0*data[i]+b1*data[i+1]+b2*data[wrp1]+b3*data[wrp2]);
		}
		else
		{
	dt[g]=sqr*(c0*data[i]+c1*data[i+1]+c2*data[i+2]+c3*data[i+3]);
	dt[g+n/2]=sqr*(b0*data[i]+b1*data[i+1]+b2*data[i+2]+b3*data[i+3]);
		}
	}
	for(i=0; i<n ;i++){data[i]=dt[i];}
	delete [] dt;
}

#endif

 Conclusion

#include <iostream.h>
#include "Daub4.h"
//example de débruitage en ondelette

int main()
{
const int N=32;
int m,n;
//série bruité
double tabint[N] = {0.61739 ,2.96222,
3.26528 ,2.27141,3.82775,4.8849 2,7.23847,7.38413,
11.1395 ,8.95568,10.6950,11.399 8,12.0592,13.6894,
14.6041 ,15.9120,17.1826,17.125 8,18.5656,18.7819,
20.3416 ,20.4525,21.5375,23.845 5,26.9481,26.8730,
26.1749 ,27.1218,28.9285,31.304 1,31.1834,31.5983,
};

//calcule du nombre d'itération nécessaire pour effectuer la transformation à toutes les échelles

m=static_cast<int>(log(N)/log(2)-1);

//on effectue la transformation compléte

for(int i=0;i<=m;i++)
{
n=static_cast<int>(N/pow(2,i));
Daub::Transform(tabint, n);
}

//on annule le premier et certain des derniers de transformation
tabint[1]=0;
for(i=15;i<=31;i++){ta bint[i]=0;}

//on reconstruit le signal avec inverse transform
for(i=m;i>=0;i--)
{
n=static_cast<int> (N/pow(2,i));
Daub::InvTrans(tabint, n);
}

//affichage des résultats
for(i=0 ; i<N;i++){cout<< tabint[i]<<endl;}
return 0;
}


 Historique

26 septembre 2004 12:27:59 :
j'ai ajouté un exemple et implémenté les modifications sugérées qui ne faisait "tiquer" le debbuger.

 Sources du même auteur

OPTIMISATION PAR GOLDEN SECTION SEARCH À UNE DIMENSION
Source avec Zip ROMBERG INTEGRATION
Source avec Zip ESTIMATION LOI STABLE (MC KULLOCH) POUR UN DLL POUR EXCEL

 Sources de la même categorie

Source avec Zip UN EXAMPLE D'APPLICATION EN CUDA DE L'ALGORITHME DE SCAN POU... par oguzaras
Source avec Zip Source avec une capture CHIFFREMENT DE VIGENERE par lajouad
Source avec Zip Source avec une capture ANALYSE SYNTAXIQUE par lajouad
Source avec Zip Source avec une capture STRUCTURE D'UNE MATRICE PAR LES LISTE LINÉAIRE (NON CONTUGUS... par benzarabel
Source avec Zip Source avec une capture DESSINER UNE ARBRE BINAIRE( MODE CONSOLE): par benzarabel

Commentaires et avis

Commentaire de coucou747 le 25/09/2004 18:09:55 administrateur CS

je ne comrpends pas grand chose, a ton explication, tu pourrais détailler stp ?

Commentaire de plus_plus_fab le 26/09/2004 01:40:14

salut,

quelques commentaires en vrac :

tu fais du C++ :
-#include <cmath>, met le en dessous des #ifndef , #define.  Et utilises l'espace de nomage std.
- tu devrais mettre tes variables static const dans la classe (private).
- Le constructeur est faux
    Daub(int n, double data[]) : T(n) {}
- indata est inutile, et T doit etre private.
- Invtrans et transform n'ont pas a etre static.
- pour l'efficacité, il faut à tout prix éviter les temporaires
   dt[g]=sqr*(c0*data[i]+c1*data[i+1]+c2*data[wrp1]+c3*data[wrp2]) en crée      un bon paquet par exemple !
- et pour finir, utlilise double*& ou lieu de double* dans la déclaration de transform, ça t'évitera de copier à la fin.

et pour aller + loin, utilise la généricité au lieu d'imposer un type double.

bonne continuation ...

Commentaire de hitcher le 26/09/2004 12:12:52

Pour explication cet opérateurs mathématiques transforme un signal de longueur 2^i , avec i appartenant à Nen 2 série de longeurs 2^(i-1)
L'une contenant les coefficients d'echelle et l'autre les coefficients de transformation.
Cet opérateur est très utilisé un débruitage, reconnaissance de forme et filtrage,....
pour plus d'info sur les ondelettes: http://www.bearcave.com/misl/misl_tech/wavelets/


Pour ce qui est du code en lui même, effecitvement je suis totalement autodidacte et je n'arrive pas à mettre en oeuvre certaine des suggestions. Par contre la non utilisation des templates est voulu car l'output ne peut en aucun cas être int.

Commentaire de plus_plus_fab le 26/09/2004 22:44:15

avec un exemple, c'est plus clair ! et ça laisse apparaitre la fragilité de la conception ...
En gros, tu donnes cette classe a quelqu'un qui connait quand meme le principe, il est incapable de faire quoi que ce soit ! (mon cas)

mon idee serait de faire une conception totalement différente qui ressemblerait à ceci :

template <class Alg>
class Signal
{
private:
      typedef typename Alg::value_type value_type;
      typedef int integer_type;
      typedef typename Alg::seuillage_type seuillage_type;
      typedef void (typename Alg::*ptr_to_fseuil)(seuillage_type);

public:
     Signal(integer_type size, const value_type* signal_input);
     ~Signal();
     void transform(const Alg&);
     void restore(ptr_to_fseuil);
     const T* coef_begin()const;
     const T* coef_end()const;
     const T* sig_res_begin()const;
     const T* sig_res_end()const;
     const T* sig_init_begin()const;
     const T* sig_init_end()const;

private:
      T* coefs;
      T* sig_res;
      const T* sig_init;
      
};

tu n'aurais plus qu'a ecrire la class Alg, et ça fonctionnerai que ce soit de la transformée de fourrier, de l'ondelettes de Haar, ou à 4 moments, ...    
cette class "Alg" (Daub),  devra fournir les types et méthodes suivantes :
types : value_type, seuillage_type.
methodes :
void transform(const value_type* input, value_type* output)const;
void restore(const value_type* input, value_type* output, ptr_to_fseuil)const;
Le constructeur sera vide, il n'y aura pas de membres pour cette classe.

pour moi, restore applique une fonction de seuillage (présente dans la classe Daub (Alg)) sur les coefficients calculées (input), puis restaure le signal à partir des coefficients seuillés (output).

Tu te sens de le faire ? Je peux t'aider en tout cas, ça m'intéresse, contacte moi en privé.

Quand ce sera fonctionnel, il faudra optimiser les calculs à coup d'élimination de temporaires et de valarray<>, on en est pas la ...

@+
      


Commentaire de plus_plus_fab le 26/09/2004 22:55:17

Si on veut separer le seuillage de la "restoration", on peut faire comme ça : (je pense que c'est mieux)

void restore(const value_type* input, value_type* output)const;

void seuillage(value_type* input_output, seuillage_type seuil)const;
j'ignore pour l'instant quel parametre invoquer pour effecteur un seuillage, donc c'est un peu au pif ça ...

 Ajouter un commentaire




Nos sponsors


Sondage...

Comparez les prix

CalendriCode

Mai 2012
LMMJVSD
 123456
78910111213
14151617181920
21222324252627
28293031   

Consulter la suite du CalendriCode

A découvrir



 
Développement réalisé par Nicolas SOREL (Nix) avec l'aide de : Cyril DURAND et Emmanuel (EBArtSoft), Merci à Vincent pour ses précieux conseils.
CodeS-SourceS.com© Toute reproduction même partielle est interdite sauf accord écrit du Webmaster
CodeS-SourceS.com© est une marque déposée tous droits réservés

Google Coop CodeS-SourceS Google Coop CodeS-SourceS
Temps d'éxécution de la page : 10,078 sec (3)

Nous contacter | Annoncer sur CodeS-SourceS | Mentions légales