1 /****************************************************************************
2 * MeshLab o o *
3 * A versatile mesh processing toolbox o o *
4 * _ O _ *
5 * Copyright(C) 2005 \/)\/ *
6 * Visual Computing Lab /\/| *
7 * ISTI - Italian National Research Council | *
8 * \ *
9 * All rights reserved. *
10 * *
11 * This program is free software; you can redistribute it and/or modify *
12 * it under the terms of the GNU General Public License as published by *
13 * the Free Software Foundation; either version 2 of the License, or *
14 * (at your option) any later version. *
15 * *
16 * This program is distributed in the hope that it will be useful, *
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
19 * GNU General Public License (http://www.gnu.org/licenses/gpl.txt) *
20 * for more details. *
21 * *
22 ****************************************************************************/
23
24 #include "utilsWatermark.h"
25 #include "qstring.h"
26 #include "filter_watermark.h"
27 #include "ERF.h"
28
29 using namespace std;
30
31 // Seed generator for the psedo random sequence in the watermark table
32 // The seed is obtained by multiplying the ascii codes of each char of the chosen string
ComputeSeed(QString string_code)33 unsigned int utilsWatermark::ComputeSeed( QString string_code )
34 {
35 unsigned int code = (unsigned int)(string_code[0].digitValue());
36
37 for (int i=1; i<string_code.length(); i++)
38 {
39 code *= (unsigned int)(string_code[i].digitValue());
40 }
41
42 return code;
43 }
44
CreateWatermarkTable(unsigned int seed)45 vector< vector<double> > utilsWatermark::CreateWatermarkTable(unsigned int seed)
46 {
47
48
49 double cellTablesize = 1.0;
50 static vcg::math::MarsenneTwisterRNG rnd;
51 rnd.initialize(seed);
52 double p = 0; // I valori generati di 'p' vanno da -RANGE_MARCHIO a RANGE_MARCHIO
53
54 int nsTheta = 360.0/cellTablesize;
55 int nsPhi = 180.0/cellTablesize;
56
57 vector<vector<double> > watermarkTable;
58 watermarkTable.resize(nsPhi);
59 for (int k = 0; k < nsPhi; ++k)
60 watermarkTable[k].resize(nsTheta);
61
62
63 int i,j;
64 for (i=0; i<nsPhi; i++)
65 {
66 for (j=0; j<nsTheta; j++)
67 {
68 p = 2.0*(rnd.generate01())*RANGE_MARCHIO - RANGE_MARCHIO;
69 watermarkTable[i][j] = p;
70 }
71 }
72
73 return watermarkTable;
74 }
75
76
77
78 // Spheric to Cartesian Coordinates System conversion
sph2cartesian(double R,double theta,double phi,float & x,float & y,float & z)79 void utilsWatermark::sph2cartesian(double R, double theta, double phi, float &x, float &y, float &z)
80 {
81 // Conversione Gradi --> Radianti
82 phi *= DEG_TO_RAD;
83 theta *= DEG_TO_RAD;
84
85 x = float(R*cos(theta)*sin(phi));
86 z = float(R*sin(theta)*sin(phi));
87 y = float(R*cos(phi));
88
89 }
90
91 // Cartesian to Spheric Coordinates System conversion
cartesian2sph(double x,double y,double z,double & R,double & theta,double & phi)92 void utilsWatermark::cartesian2sph(double x, double y, double z, double &R, double &theta, double &phi)
93 {
94 R = sqrt( x*x + y*y + z*z );
95 theta = atan2( z,x );
96 phi = acos( y / R );
97
98 // Conversione Radianti --> Gradi
99 phi *= RAD_TO_DEG;
100 theta *= RAD_TO_DEG;
101 }
102
103 // rounding
round(double x)104 double utilsWatermark::round(double x)
105 {
106 double app = floor(x);
107
108 if (x - app >= 0.5)
109 return (app + 1.0);
110 else
111 return app;
112 }
113
114 // evaluate threshold
thresholdRo(double muRhoH0,double varRhoH0,double muRhoH1,double varRhoH1,double Pf,double & Pm)115 double utilsWatermark::thresholdRo(double muRhoH0, double varRhoH0,
116 double muRhoH1, double varRhoH1,
117 double Pf, double &Pm)
118 {
119
120 Erf erf;
121 double Y;
122 double Z;
123
124 // Calcolo di T(ro)
125 /////////////////////////////////////////////////
126
127
128 Y = erf.inverfc( 2.0*Pf );
129
130 double Trho = sqrt(2.0*varRhoH0)*Y + muRhoH0;
131
132 // Calcolo di Pm
133 ///////////////////////////////////////////////////////
134
135 //Z = gsl_sf_erfc( sqrt( ((muRhoH1 - Trho)*(muRhoH1 - Trho))/(2.0*varRhoH1)) );
136 //Z = erfc( sqrt( ((muRhoH1 - Trho)*(muRhoH1 - Trho))/(2.0*varRhoH1)) );
137 Z = erf.erfc( sqrt( ((muRhoH1 - Trho)*(muRhoH1 - Trho))/(2.0*varRhoH1)) );
138
139 Pm = 0.5*Z;
140
141 return Trho;
142 }
143