← Back to branch summary

~ubuntu-branches/ubuntu/precise/octave-nurbs/precise

« back to all changes in this revision

Viewing changes to src/bspeval.cc

  • Committer: Bazaar Package Importer
  • Author(s): Thomas Weber
  • Date: 2011-04-20 21:53:46 UTC
  • mfrom: (3.1.2 sid)
  • Revision ID: james.westby@ubuntu.com-20110420215346-mtrs63hlmsh5yuaz
Tags: 1.3.3-1
* New upstream release
* Bump standards version to 3.9.1, no changes needed

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/bspeval.cc
1
1
/* Copyright (C) 2009 Carlo de Falco
2
2
  
3
 
   This program is free software; you can redistribute it and/or modify
 
3
   This program is free software: you can redistribute it and/or modify
4
4
   it under the terms of the GNU General Public License as published by
5
 
   the Free Software Foundation; either version 2 of the License, or
 
5
   the Free Software Foundation, either version 2 of the License, or
6
6
   (at your option) any later version.
7
 
 
 
7
 
8
8
   This program is distributed in the hope that it will be useful,
9
9
   but WITHOUT ANY WARRANTY; without even the implied warranty of
10
10
   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
11
11
   GNU General Public License for more details.
12
12
 
13
13
   You should have received a copy of the GNU General Public License
14
 
   along with this program; if not, write to the Free Software
15
 
   Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
 
14
   along with this program.  If not, see <http://www.gnu.org/licenses/>.
16
15
*/
17
16
 
18
17
#include <octave/oct.h>
19
18
#include "low_level_functions.h"
 
19
#include <omp.h>
20
20
 
21
21
static bool bspeval_bad_arguments(const octave_value_list& args);
22
22
 
40
40
       p - Evaluated points, matrix of size (dim,nu)\n\
41
41
")
42
42
{
43
 
  
44
 
  
45
 
 
46
 
  int       d = args(0).int_value();
47
 
  const Matrix    c = args(1).matrix_value();
48
 
  const RowVector k = args(2).row_vector_value();
49
 
  const NDArray   u = args(3).array_value();
50
 
  
51
 
  octave_idx_type nu = u.length();
52
 
  octave_idx_type mc = c.rows(),
53
 
    nc = c.cols();
54
 
 
55
 
  Matrix p(mc, nu, 0.0);
56
 
  RowVector N(d+1,0.0);
57
43
 
58
44
  octave_value_list retval;
59
 
  if (!error_state)
60
 
    {
61
 
      if (nc + d == k.length() - 1) 
62
 
        {        
63
 
          int s, tmp1;
64
 
          double tmp2;
65
 
          
66
 
          for (octave_idx_type col(0); col<nu; col++)
67
 
            {   
68
 
              s = findspan(nc-1, d, u(col), k);
69
 
              basisfun(s, u(col), d, k, N);    
70
 
              tmp1 = s - d;                
71
 
              for (octave_idx_type row(0); row<mc; row++)
72
 
                {
73
 
                  double tmp2 = 0.0;
74
 
                  for ( octave_idx_type i(0); i<=d; i++)                   
75
 
                    tmp2 +=  N(i)*c(row,tmp1+i);          
76
 
                  p(row,col) = tmp2;
77
 
                }             
78
 
            }   
79
 
        } 
80
 
      else 
81
 
        {
82
 
          error("inconsistent bspline data, d + columns(c) != length(k) - 1.");
83
 
        }
84
 
    }
85
 
  retval(0) = octave_value(p);
 
45
  if (!bspeval_bad_arguments (args))
 
46
    {      
 
47
      int             d = args(0).int_value();
 
48
      const Matrix    c = args(1).matrix_value();
 
49
      const RowVector k = args(2).row_vector_value();
 
50
      const NDArray   u = args(3).array_value();
 
51
      
 
52
      octave_idx_type nu = u.length();
 
53
      octave_idx_type mc = c.rows(),
 
54
        nc = c.cols();
 
55
      
 
56
      Matrix p(mc, nu, 0.0);
 
57
      
 
58
      if (!error_state)
 
59
        {
 
60
          if (nc + d == k.length() - 1) 
 
61
            {    
 
62
#pragma omp parallel default (none) shared (d, c, k, u, nu, mc, nc, p)
 
63
              {
 
64
                RowVector N(d+1,0.0);
 
65
                int s, tmp1;
 
66
                double tmp2;
 
67
#pragma omp for 
 
68
                for (octave_idx_type col=0; col<nu; col++)
 
69
                  {     
 
70
                    //printf ("thread %d, col %d\n", omp_get_thread_num (), col);
 
71
                    s = findspan (nc-1, d, u(col), k);
 
72
                    basisfun (s, u(col), d, k, N);    
 
73
                    tmp1 = s - d;                
 
74
                    for (octave_idx_type row(0); row<mc; row++)
 
75
                      {
 
76
                        tmp2 = 0.0;
 
77
                        for ( octave_idx_type i(0); i<=d; i++)                   
 
78
                          tmp2 +=  N(i)*c(row,tmp1+i);    
 
79
                        p(row,col) = tmp2;
 
80
                      }             
 
81
                  }   
 
82
              }// end omp 
 
83
            }
 
84
          else 
 
85
            {
 
86
              error("inconsistent bspline data, d + columns(c) != length(k) - 1.");
 
87
            }
 
88
          retval(0) = octave_value(p);
 
89
        }
 
90
    }      
86
91
  return retval;
87
92
} 
88
93
 
89
 
static bool bspeval_bad_arguments(const octave_value_list& args) 
 
94
static bool bspeval_bad_arguments (const octave_value_list& args) 
90
95
{ 
91
96
  if (args.length() != 4)
92
97
    {
93
 
      error("wrong number of input arguments.");
 
98
      error("bspeval: wrong number of input arguments.");
94
99
      return true;
95
100
    }
96
101
  if (!args(0).is_real_scalar()) 
97
102
    { 
98
 
      error("degree should be a scalar."); 
 
103
      error("bspeval: degree should be a scalar."); 
99
104
      return true; 
100
105
    } 
101
106
  if (!args(1).is_real_matrix()) 
102
107
    { 
103
 
      error("the control net should be a matrix of doubles."); 
 
108
      error("bspeval: the control net should be a matrix of doubles."); 
104
109
      return true; 
105
110
    } 
106
111
  if (!args(2).is_real_matrix()) 
107
112
    { 
108
 
      error("the knot vector should be a real vector."); 
 
113
      error("bspeval: the knot vector should be a real vector."); 
109
114
      return true; 
110
115
    } 
111
116
  if (!args(3).is_real_type()) 
112
117
    { 
113
 
      error("the set of parametric points should be an array of doubles."); 
 
118
      error("bspeval: the set of parametric points should be an array of doubles."); 
114
119
      return true; 
115
120
    } 
116
121
  return false;