[prev in list] [next in list] [prev in thread] [next in thread] 

List:       kde-commits
Subject:    playground/graphics/tamponi/painterlyframework/utilities/src
From:       Emanuele Tamponi <emanuele () valinor ! it>
Date:       2007-12-01 15:48:08
Message-ID: 1196524088.508853.13311.nullmailer () svn ! kde ! org
[Download RAW message or body]

SVN commit 743647 by tamponi:

Add a bit of apidox and some code in the ProfileGenerator.


 M  +13 -1     function.h  
 M  +69 -12    ipg/profile_generator.cpp  


--- trunk/playground/graphics/tamponi/painterlyframework/utilities/src/function.h \
#743646:743647 @@ -16,7 +16,19 @@
 #define PTS 10000
 
 /**
-
+ This class provides a wrapper around mathematical functions.
+ You can use it to define and then do mathematical operations on mathematical \
functions. + Consider for example the functions:
+   y = x^2 + 2
+   y = x^6 + 2x^4 + 5x
+ They are polynomials and can be defined as:
+   QVector\<double\> fc(3, 0); fc[0] = 2; fc[2] = 1;
+   QVector\<double\> gc(7, 0); gc[1] = 5; gc[4] = 2; gc[6] = 1;
+   Function f = Polynomial(fc);
+   Function g = Polynomial(gc);
+ You can now use f and g as they were mathematical functions. You can sum, divide, \
multiply or subtract them, + integrate them, do polynomial interpolation on them. The \
integration feature is provided for the interval + [-1, 1] only on generic functions. \
                You can convert them to Polynomial using the toPolynomial method.
  */
 class Function {
 
--- trunk/playground/graphics/tamponi/painterlyframework/utilities/src/ipg/profile_generator.cpp \
#743646:743647 @@ -13,17 +13,47 @@
 #include <gsl/gsl_matrix.h>
 #include <gsl/gsl_vector.h>
 
+double XYZtoRGB[9] = {
+    3.2406, -1.5372, -0.4986,
+   -0.9689,  1.8758,  0.0415,
+    0.0557, -0.2040,  1.0570
+};
+
+double bradfordAdaptation[9] = {
+    0.8951,  0.2664, -0.1614,
+   -0.7502,  1.7135,  0.0367,
+    0.0389, -0.0685,  1.0296
+};
+
+double bradfordAdaptationi[9] = {
+    0.986993, -0.147054, 0.159963,
+    0.432305,  0.518360, 0.049291,
+   -0.008529,  0.040043, 0.968487
+};
+
+double D65White[3] = {
+    0.950468, 1.0, 1.088970
+};
+
 ProfileGenerator::ProfileGenerator(unsigned int N, const Function &H)
     : m_N(N), m_H(H)
 {
     // 1- Load CMFs
     // 2- Generate GQ for each CMF*m_hFunct
-    // 3- Build the H matrix
-    // 4- Build Chromaticity Adaptation matrix
-    // 5- Build R as M^-1 * A^-1 * H
+    // 3- Build the R matrix (from Reflectance to XYZ-H)
+    // 4- Build Chromatic Adaptation (from H to D65) matrix C
+    // 5- Build m_T as M * C * R where M is the transformation matrix from XYZ-D65 \
to RGB  // 6- Order wavelenghts
     // 7- Save matrix and wavelenght order ( in save() )
 
+    m_T = gsl_matrix_calloc(3, 3*N);
+
+    gsl_matrix *R = gsl_matrix_calloc(3, 3*N);
+    gsl_matrix *C = gsl_matrix_calloc(3, 3);
+    gsl_matrix_view M = gsl_matrix_view_array(XYZtoRGB, 3, 3);
+
+    gsl_matrix *tmp = gsl_matrix_calloc(3, 3);
+
     // 1
     Function x("resources/CMFs/2deg_x.txt");
     Function y("resources/CMFs/2deg_y.txt");
@@ -35,21 +65,48 @@
     GaussianQuadrature wz(N, z*H);
 
     // 3
-    m_T = gsl_matrix_alloc(3, 3*N);
-    gsl_matrix_set_zero(m_T);
-
     double k = 0;
-    for (int j = 0; j < N; j++)
+    for (int j = 0; j < N; j++) {
+        gsl_matrix_set(R, 0, 0*N+j, wx.weights()[j]);
+        gsl_matrix_set(R, 1, 1*N+j, wy.weights()[j]);
+        gsl_matrix_set(R, 2, 2*N+j, wz.weights()[j]);
         k += wy.weights()[j];
-    for (int j = 0; j < N; j++) {
-        gsl_matrix_set(m_T, 0, 0*N, wx.weights()[j]);
-        gsl_matrix_set(m_T, 1, 1*N, wy.weights()[j]);
-        gsl_matrix_set(m_T, 2, 2*N, wz.weights()[j]);
     }
-    gsl_matrix_scale(m_T, 1.0/k);
+    gsl_matrix_scale(R, 1.0/k);
 
     // 4
+    gsl_vector *HW = gsl_vector_calloc(3);
+    gsl_vector *rhoS = gsl_vector_calloc(3);
+    gsl_vector *rhoD = gsl_vector_calloc(3);
+    gsl_matrix *rho = gsl_matrix_calloc(3, 3);
 
+    gsl_matrix_view MA  = gsl_matrix_view_array(bradfordAdaptation,  3, 3);
+    gsl_matrix_view MAi = gsl_matrix_view_array(bradfordAdaptationi, 3, 3);
+
+    gsl_vector_view D65W = gsl_vector_view_array(D65White, 3);
+    for (int i = 0; i < 3; i++) {
+        double *x = gsl_vector_ptr(HW, i);
+        for (int j = 0; j < N; j++)
+            *x += gsl_matrix_get(R, i, i*N+j);
+    }
+
+    gsl_blas_dgemv(CblasNoTrans, 1.0, MA, HW,           0.0, rhoS);
+    gsl_blas_dgemv(CblasNoTrans, 1.0, MA, &D65W.vector, 0.0, rhoD);
+
+    for (int i = 0; i < 3; i++)
+        gsl_matrix_set(rho, i, i, gsl_vector_get(rhoD, i) / gsl_vector_get(rhoD, \
i)); +
+    gsl_blas_dgemm(CBlasNoTrans, CBlasNoTrans, 1.0, MAi, rho, 0.0, tmp);
+    gsl_blas_dgemm(CBlasNoTrans, CBlasNoTrans, 1.0, MA,  rho, 0.0, C);
+
+    gsl_matrix_free(rho);
+    gsl_vector_free(rhoD);
+    gsl_vector_free(rhoS);
+    gsl_vector_free(HW);
+
+    // 5
+    gsl_blas_dgemm(CBlasNoTrans, CBlasNoTrans, 1.0, M,   C, 0.0, tmp);
+    gsl_blas_dgemm(CBlasNoTrans, CBlasNoTrans, 1.0, tmp, R, 0.0, m_T);
 }
 
 ProfileGenerator::~ProfileGenerator()


[prev in list] [next in list] [prev in thread] [next in thread]