source: trunk/src/Locator.cpp@ 2730

Last change on this file since 2730 was 2730, checked in by Takeshi Nakazato, 12 years ago

New Development: No

JIRA Issue: Yes CAS-4770

Ready for Test: Yes

Interface Changes: No

What Interface Changed: Please list interface changes

Test Programs: List test programs

Put in Release Notes: Yes/No

Module(s): Module Names change impacts.

Description: Describe your changes here...

Rewrite implementations for locator and interpolator.
Documentation (doxygen format) is added to header files.


File size: 3.0 KB
Line 
1//
2// C++ Implementation: Locator
3//
4// Description:
5//
6//
7// Author: Takeshi Nakazato <takeshi.nakazato@nao.ac.jp>, (C) 2012
8//
9// Copyright: See COPYING file that comes with this distribution
10//
11//
12#include <assert.h>
13
14#include "Locator.h"
15
16namespace asap {
17Locator::Locator()
18 : x_(0),
19 n_(0),
20 ascending_(true),
21 copy_(false)
22{
23}
24
25Locator::Locator(double *v, unsigned int n, bool copystorage)
26 : x_(0),
27 n_(0),
28 ascending_(true),
29 copy_(false)
30{
31 set(v, n, copystorage);
32}
33
34Locator::~Locator()
35{
36 if (copy_ && x_)
37 delete[] x_;
38}
39
40void Locator::set(double *v, unsigned int n, bool copystorage)
41{
42 if (copy_) {
43 if (!copystorage || n > n_) {
44 delete[] x_;
45 x_ = 0;
46 }
47 }
48 copy_ = copystorage;
49 n_ = n;
50 if (copy_) {
51 if (!x_)
52 x_ = new double[n_];
53 for (unsigned int i = 0; i < n_; i++)
54 x_[i] = v[i];
55 }
56 else {
57 x_ = v;
58 }
59 ascending_ = (x_[0] <= x_[n_-1]);
60}
61
62unsigned int Locator::bisection(double x, unsigned int left, unsigned int right)
63{
64 unsigned int jl = left;
65 unsigned int ju = right;
66
67 if (ascending_) {
68 // ascending order
69 if (x <= x_[0])
70 return 0;
71 else if (x > x_[n_-1])
72 return n_;
73
74 unsigned int jm;
75 while (ju - jl > 1) {
76 jm = (ju + jl) / 2;
77 if (x > x_[jm])
78 jl = jm;
79 else
80 ju = jm;
81 }
82 }
83 else {
84 // descending order
85 if (x >= x_[0])
86 return 0;
87 else if (x < x_[n_-1])
88 return n_;
89
90 unsigned int jm;
91 while (ju - jl > 1) {
92 jm = (ju + jl) / 2;
93 if (x < x_[jm])
94 jl = jm;
95 else
96 ju = jm;
97 }
98 }
99
100 return ju;
101}
102
103void Locator::hunt(double x, unsigned int &left, unsigned int &right)
104{
105 unsigned int inc = 1;
106 if (ascending_) {
107 // ascending order
108 if (x >= x_[left]) {
109 // forward hunt
110 if (left >= n_ - 1) {
111 right = n_;
112 return;
113 }
114 right = left + inc;
115 while (x >= x_[right]) {
116 left = right;
117 inc *= 2;
118 right = left + inc;
119 if (right > n_ - 1) {
120 right = n_;
121 break;
122 }
123 }
124 }
125 else {
126 // backward hunt
127 if (left == 0) {
128 right = 0;
129 return;
130 }
131 right = left;
132 left -= inc;
133 while (x < x_[left]) {
134 right = left;
135 inc *= 2;
136 if (inc >= right) {
137 left = 0;
138 break;
139 }
140 left = right - inc;
141 }
142 }
143 }
144 else {
145 // descending order
146 if (x < x_[left]) {
147 // forward hunt
148 if (left >= n_ - 1) {
149 right = n_;
150 return;
151 }
152 right = left + inc;
153 while (x < x_[right]) {
154 left = right;
155 inc *= 2;
156 right = left + inc;
157 if (right > n_ - 1) {
158 right = n_;
159 break;
160 }
161 }
162 }
163 else {
164 // backward hunt
165 if (left == 0)
166 return;
167 right = left;
168 left -= inc;
169 while (x >= x_[left]) {
170 right = left;
171 inc *= 2;
172 if (inc >= right) {
173 left = 0;
174 break;
175 }
176 left = right - inc;
177 }
178 }
179 }
180}
181}
Note: See TracBrowser for help on using the repository browser.