rt300@34
|
1 //
|
rt300@34
|
2 // hilbert.cpp
|
rt300@34
|
3 // sonicZoom
|
rt300@34
|
4 //
|
rt300@34
|
5 // Created by Robert Tubb on 04/04/2013.
|
rt300@34
|
6 //
|
rt300@34
|
7 //
|
rt300@34
|
8
|
rt300@34
|
9 #include "hilbert.h"
|
rt300@34
|
10
|
rt300@43
|
11 // UTILS
|
rt300@34
|
12 //--------------------------------------------------------------------
|
rt300@34
|
13 unsigned int rotl(unsigned int value, int shift, int digits) {
|
rt300@34
|
14 unsigned int ones = (1 << digits) - 1;
|
rt300@34
|
15 if(shift >= digits)
|
rt300@34
|
16 shift -= digits;
|
rt300@34
|
17
|
rt300@34
|
18 unsigned int out = (value << shift) | (value >> (digits - shift));
|
rt300@34
|
19 return out & ones;
|
rt300@34
|
20 }
|
rt300@34
|
21 //--------------------------------------------------------------------
|
rt300@34
|
22 unsigned int rotr(unsigned int value, int shift, int digits) {
|
rt300@34
|
23 int ones = (1 << digits) - 1;
|
rt300@34
|
24 // can only handle -digits to 2*digits-1
|
rt300@34
|
25 if(shift >= digits) shift = shift - digits;
|
rt300@34
|
26 if(shift < 0) shift = digits + shift;
|
rt300@34
|
27 unsigned int out = (value >> shift) | (value << (digits - shift));
|
rt300@34
|
28 return out & ones;
|
rt300@34
|
29 }
|
rt300@34
|
30 //--------------------------------------------------------------------
|
rt300@34
|
31 unsigned int bin2dec(vector<bool> binary){
|
rt300@34
|
32 // takes in a vector of booleans. we assume that the leftmost bool (highest index) is the MSB
|
rt300@34
|
33 unsigned int dec = 0;
|
rt300@34
|
34 int i,pwr;
|
rt300@34
|
35 int N = binary.size();
|
rt300@34
|
36
|
rt300@34
|
37 for (i = 0; i<N; i++){
|
rt300@34
|
38 pwr = N-i-1;
|
rt300@34
|
39 dec += binary[i]*(1 << pwr);
|
rt300@34
|
40
|
rt300@34
|
41 }
|
rt300@34
|
42
|
rt300@34
|
43 return dec;
|
rt300@34
|
44 }
|
rt300@34
|
45 //--------------------------------------------------------------------
|
rt300@34
|
46 vector<bool> dec2bin(unsigned int number, unsigned int size){
|
rt300@34
|
47 vector<bool> bin;
|
rt300@34
|
48 int i, pwr;
|
rt300@34
|
49 for(i=0;i<size;i++){
|
rt300@34
|
50 bin.push_back(0);
|
rt300@34
|
51 }
|
rt300@34
|
52 if(number >= (1 << size)){
|
rt300@34
|
53 cout << " bad size for Hilbert::dec2bin" << "\n";
|
rt300@34
|
54 return bin;
|
rt300@34
|
55 }
|
rt300@34
|
56 for(i=0;i<size;i++){
|
rt300@34
|
57 pwr = size-i;
|
rt300@34
|
58 bin[i] = (number/(1 << pwr)) % 2;
|
rt300@34
|
59 bin.push_back(0);
|
rt300@34
|
60 }
|
rt300@35
|
61 return bin;
|
rt300@34
|
62
|
rt300@34
|
63
|
rt300@34
|
64 }
|
rt300@35
|
65 //--------------------------------------------------------------------
|
rt300@34
|
66 void printBool(vector <bool> vcode){
|
rt300@34
|
67 vector<bool>::iterator bit;
|
rt300@34
|
68
|
rt300@34
|
69 for(bit=vcode.begin(); bit!=vcode.end() ; bit++){
|
rt300@34
|
70
|
rt300@34
|
71 cout << *bit;
|
rt300@34
|
72 }
|
rt300@34
|
73 cout << "\n";
|
rt300@34
|
74
|
rt300@34
|
75 }
|
rt300@34
|
76 //--------------------------------------------------------------------
|
rt300@43
|
77 //====================================================================
|
rt300@43
|
78 //--------------------------------------------------------------------
|
rt300@43
|
79
|
rt300@43
|
80 Hilbert::Hilbert(){
|
rt300@43
|
81 N = 7;
|
rt300@43
|
82 P = 5;
|
rt300@43
|
83 codeLength = (int)pow(2.0,P);
|
rt300@43
|
84
|
rt300@43
|
85
|
rt300@43
|
86 makeCode();
|
rt300@43
|
87 makeRotationRules();
|
rt300@43
|
88
|
rt300@43
|
89 }
|
rt300@43
|
90
|
rt300@43
|
91 //--------------------------------------------------------------------
|
rt300@43
|
92
|
rt300@43
|
93 Hilbert::Hilbert(int aN, int aP){
|
rt300@34
|
94 N = aN;
|
rt300@34
|
95 P = aP;
|
rt300@34
|
96 codeLength = (int)pow(2.0,P);
|
rt300@43
|
97
|
rt300@43
|
98
|
rt300@34
|
99 makeCode();
|
rt300@34
|
100 makeRotationRules();
|
rt300@34
|
101
|
rt300@43
|
102 }
|
rt300@43
|
103 //--------------------------------------------------------------------
|
rt300@43
|
104
|
rt300@43
|
105 void Hilbert::changeCurve(int aN, int aP){
|
rt300@43
|
106 N = aN;
|
rt300@43
|
107 P = aP;
|
rt300@43
|
108 codeLength = (int)pow(2.0,P);
|
rt300@34
|
109
|
rt300@34
|
110
|
rt300@43
|
111 makeCode();
|
rt300@43
|
112 makeRotationRules();
|
rt300@34
|
113
|
rt300@34
|
114 }
|
rt300@34
|
115 //--------------------------------------------------------------
|
rt300@34
|
116 void Hilbert::makeRotationRules(){
|
rt300@43
|
117 theRotations.clear();
|
rt300@43
|
118 theEntryVertices.clear();
|
rt300@34
|
119
|
rt300@43
|
120 const int *rp;
|
rt300@43
|
121 const int *ep;
|
rt300@43
|
122 if(P == 5){
|
rt300@43
|
123 rp = rotations5;
|
rt300@43
|
124 ep = entryVertices5;
|
rt300@43
|
125 }else if(P == 4){
|
rt300@43
|
126 rp = rotations4;
|
rt300@43
|
127 ep = entryVertices4;
|
rt300@43
|
128 }else if(P == 3){
|
rt300@43
|
129 rp = rotations3;
|
rt300@43
|
130 ep = entryVertices3;
|
rt300@43
|
131 }else if(P == 2){
|
rt300@43
|
132 rp = rotations2;
|
rt300@43
|
133 ep = entryVertices2;
|
rt300@43
|
134 }else if(P == 1){ // will this even work?
|
rt300@43
|
135 rp = rotations1;
|
rt300@43
|
136 ep = entryVertices1;
|
rt300@43
|
137 }else{
|
rt300@43
|
138 cout << "ERROR: makeRotationRules: bad numb of MIDI params";
|
rt300@43
|
139 return;
|
rt300@43
|
140 }
|
rt300@43
|
141
|
rt300@43
|
142 for(int i=0; i < codeLength; i++){ // fill vector
|
rt300@43
|
143
|
rt300@43
|
144 theEntryVertices.push_back(ep[i]);
|
rt300@43
|
145 theRotations.push_back(rp[i]);
|
rt300@43
|
146
|
rt300@43
|
147 }
|
rt300@34
|
148
|
rt300@34
|
149 }
|
rt300@34
|
150 //--------------------------------------------------------------
|
rt300@34
|
151 void Hilbert::makeCode(){
|
rt300@34
|
152
|
rt300@34
|
153 ////////////////////////////////~~~,,,,,,,,.........
|
rt300@34
|
154 ////////// gray code generation
|
rt300@34
|
155 //////
|
rt300@34
|
156 ///
|
rt300@43
|
157 // choose the gray code sequence and rotation based on the numnber of dimensions P
|
rt300@34
|
158
|
rt300@34
|
159 // ALL SEQUENCES SHOULD END WITH 1 0 0 ... i.e. MSB changes
|
rt300@34
|
160
|
rt300@34
|
161 // max MRL 5 digit
|
rt300@43
|
162
|
rt300@43
|
163 theGrayCode.clear();
|
rt300@43
|
164 theGrayCodeD.clear();
|
rt300@43
|
165
|
rt300@43
|
166 int trans5[] = {2,3,4,0,2,1,4,3,2,0,4,3,2,1,4,0,2,3,4,0,2,1,4,3,2,0,4,3,2,1,4,0};
|
rt300@34
|
167
|
rt300@34
|
168 // other ones
|
rt300@34
|
169 int trans4[] = {1,2,0,3,0,2,1,2,3,0,3,2,1,3,1,0};
|
rt300@34
|
170 int trans3[] = {1,2,1,0,1,2,1,0};
|
rt300@34
|
171 int trans2[] = {1,0,1,0};
|
rt300@43
|
172 int trans1[] = {0,0};
|
rt300@34
|
173 // balanced 5
|
rt300@34
|
174 int transB5[] = {1,2,3,4,5,1,5,2,3,5,2,4,2,3,4,1,4,3,2,3,1,5,3,4,1,5,2,5,3,4,1,3};
|
rt300@34
|
175
|
rt300@43
|
176 int *tp;
|
rt300@43
|
177
|
rt300@34
|
178
|
rt300@34
|
179 if(P == 5){
|
rt300@43
|
180 tp = trans5;
|
rt300@43
|
181 }else if(P == 4){
|
rt300@43
|
182 tp = trans4;
|
rt300@43
|
183 }else if(P == 3){
|
rt300@43
|
184 tp = trans3;
|
rt300@43
|
185 }else if(P == 2){
|
rt300@43
|
186 tp = trans2;
|
rt300@43
|
187 }else if(P == 1){ // will this even work?
|
rt300@43
|
188 tp = trans1;
|
rt300@34
|
189 }
|
rt300@43
|
190
|
rt300@34
|
191
|
rt300@34
|
192 int code[codeLength][P]; // start with normal array
|
rt300@34
|
193 for(int j=0; j<P; j++){
|
rt300@34
|
194 code[0][j] = 0;
|
rt300@34
|
195 }
|
rt300@34
|
196
|
rt300@34
|
197 for(int i=0; i < codeLength-1; i++){ // don't need last 3
|
rt300@34
|
198 for(int j=0; j<P; j++){
|
rt300@43
|
199 if (j == abs(tp[i])){
|
rt300@34
|
200 code[i+1][j] = !code[i][j];
|
rt300@34
|
201 }else{
|
rt300@34
|
202 code[i+1][j] = code[i][j];
|
rt300@34
|
203 }
|
rt300@34
|
204 }
|
rt300@34
|
205
|
rt300@34
|
206 }
|
rt300@34
|
207
|
rt300@34
|
208 for(int i=0; i < codeLength; i++){ // fill vector
|
rt300@34
|
209
|
rt300@34
|
210 theGrayCode.push_back(vector<bool>());
|
rt300@34
|
211
|
rt300@34
|
212 for(int j=0; j<P; j++){
|
rt300@34
|
213 theGrayCode[i].push_back(code[i][j]);
|
rt300@34
|
214
|
rt300@34
|
215 }
|
rt300@34
|
216
|
rt300@34
|
217 }
|
rt300@34
|
218
|
rt300@34
|
219 for(int i=0;i<codeLength;i++){
|
rt300@34
|
220 theGrayCodeD.push_back(bin2dec(theGrayCode[i]));
|
rt300@34
|
221 //cout << bin2dec(theGrayCode[i]) << "\n";
|
rt300@34
|
222 //cout << theGrayCodeD[i] << "\n";
|
rt300@34
|
223 }
|
rt300@34
|
224
|
rt300@34
|
225 }
|
rt300@34
|
226
|
rt300@34
|
227 //--------------------------------------------------------------------
|
rt300@34
|
228 //--------------------------------------------------------------------
|
rt300@34
|
229 vector<int> Hilbert::calculateParamsFromIndex(unsigned long long index){
|
rt300@34
|
230 // set %set start and direction of biggest cube
|
rt300@34
|
231 unsigned int i;
|
rt300@34
|
232 unsigned int mask = (1 << P) - 1;
|
rt300@34
|
233 unsigned int entryPoint = 0;
|
rt300@34
|
234 int direction = P - 1, directionInv = P-1;
|
rt300@34
|
235 unsigned int vertex=0, subindex=0, newe=0, newd=0, entryPointInv = 0;
|
rt300@34
|
236
|
rt300@34
|
237 unsigned int pbin[N];
|
rt300@34
|
238
|
rt300@34
|
239 vector<int> params;
|
rt300@34
|
240 for(i=0; i<P;i++){
|
rt300@34
|
241 params.push_back(0);
|
rt300@34
|
242 }
|
rt300@34
|
243 for(i=0;i<N;i++){
|
rt300@34
|
244 pbin[i] = 0;
|
rt300@34
|
245
|
rt300@34
|
246 }
|
rt300@34
|
247
|
rt300@34
|
248 // for loop thru bit levels
|
rt300@36
|
249
|
rt300@34
|
250 for(int blev = N-1; blev >= 0; blev--){
|
rt300@34
|
251 // get next highest bits of index
|
rt300@34
|
252
|
rt300@34
|
253 subindex = index >> blev*P;
|
rt300@34
|
254 subindex = mask & subindex;
|
rt300@34
|
255
|
rt300@35
|
256 //cout << "subindex: " << subindex << "\n";
|
rt300@34
|
257 // which vertex corresponds to this index?
|
rt300@34
|
258 if(subindex < theGrayCode.size()) vertex = theGrayCodeD[subindex];
|
rt300@34
|
259
|
rt300@35
|
260 //cout << "rotated: " << vertex << "\n";
|
rt300@34
|
261
|
rt300@34
|
262
|
rt300@34
|
263 //% inverse rotate this T_e,d
|
rt300@34
|
264 entryPointInv = rotr(entryPoint, direction+1,P);
|
rt300@34
|
265 directionInv = (P - direction - 2);
|
rt300@34
|
266 vertex = rotate(vertex,entryPointInv,directionInv);
|
rt300@34
|
267
|
rt300@35
|
268 //cout << "vertex: " << vertex << "\n";
|
rt300@34
|
269
|
rt300@34
|
270 //% build up bits of params
|
rt300@34
|
271
|
rt300@34
|
272 for(int j=0;j<P;j++){
|
rt300@34
|
273
|
rt300@34
|
274 params[j] += ((vertex >> P-j-1) & 1)*(1<<blev);
|
rt300@34
|
275 }
|
rt300@43
|
276 if(subindex > theEntryVertices.size()){
|
rt300@43
|
277 cout << "ERROR: subindex too large for theEntryVertices\n";
|
rt300@34
|
278 }
|
rt300@43
|
279 newe = theEntryVertices[subindex];
|
rt300@43
|
280 newd = theRotations[subindex];
|
rt300@34
|
281 // next rotations...
|
rt300@34
|
282 int lse2 = rotr(newe, P-direction-1,P);
|
rt300@34
|
283 entryPoint = entryPoint ^ lse2;
|
rt300@34
|
284 direction = (direction + newd + 1) % P;
|
rt300@35
|
285 //cout << "------------\n";
|
rt300@34
|
286 }
|
rt300@34
|
287
|
rt300@34
|
288 return params;
|
rt300@34
|
289 }
|
rt300@34
|
290 //--------------------------------------------------------------------
|
rt300@34
|
291 int Hilbert::rotate(int vertex, int entryPoint, int direction){
|
rt300@34
|
292 vertex = vertex ^ entryPoint;
|
rt300@34
|
293 vertex = rotr(vertex,direction+1,P);
|
rt300@34
|
294 return vertex;
|
rt300@34
|
295 }
|
rt300@34
|
296 //--------------------------------------------------------------------
|
rt300@34
|
297
|
rt300@34
|
298
|
rt300@34
|
299 //--------------------------------------------------------------------
|
rt300@35
|
300 unsigned long long Hilbert::calculateIndexFromParams(vector<int> params){
|
rt300@34
|
301 // set %set start and direction of biggest cube
|
rt300@34
|
302
|
rt300@35
|
303 unsigned long long h = 0;
|
rt300@34
|
304
|
rt300@34
|
305 unsigned int i;
|
rt300@34
|
306 unsigned int entryPoint = 0;
|
rt300@34
|
307 int direction = P - 1;
|
rt300@35
|
308 unsigned int vertex=0, newe=0, newd=0;
|
rt300@35
|
309 unsigned long long subindex; // needs to be long cos we lewt shift by alot
|
rt300@34
|
310 vector<int> subindices;
|
rt300@34
|
311
|
rt300@34
|
312
|
rt300@34
|
313 // for loop thru bit levels
|
rt300@36
|
314
|
rt300@34
|
315 for(int blev = N-1; blev >= 0; blev--){
|
rt300@34
|
316 //% get next highest bit of param
|
rt300@34
|
317
|
rt300@34
|
318 vector<bool> vertexb;
|
rt300@34
|
319 for(i=0;i<P;i++){
|
rt300@34
|
320 vertexb.push_back((params[i]/(1 << blev)) % 2);
|
rt300@34
|
321
|
rt300@34
|
322 }
|
rt300@35
|
323 //printBool(vertexb);
|
rt300@34
|
324 vertex = bin2dec(vertexb);
|
rt300@35
|
325 //cout << "INV vertex: " << vertex << "\n";
|
rt300@34
|
326 // rotate it
|
rt300@34
|
327 vertex = rotate(vertex,entryPoint,direction);
|
rt300@34
|
328
|
rt300@34
|
329 // % get new index of thisun
|
rt300@35
|
330 //cout << "INV rotated: " << vertex << "\n";
|
rt300@34
|
331
|
rt300@34
|
332
|
rt300@34
|
333 for(i=0;i<codeLength;i++){
|
rt300@34
|
334 if(vertex == theGrayCodeD[i]){
|
rt300@34
|
335 subindex = i;
|
rt300@34
|
336 }
|
rt300@34
|
337 }
|
rt300@35
|
338 //cout << "INV subindex: " << subindex << "\n";
|
rt300@34
|
339 // work out next rotations
|
rt300@43
|
340 if(subindex > theEntryVertices.size()){
|
rt300@43
|
341 cout << "ERROR: subindex too large for theEntryVertices\n";
|
rt300@43
|
342 }
|
rt300@43
|
343 newe = theEntryVertices[subindex];
|
rt300@43
|
344 newd = theRotations[subindex];
|
rt300@34
|
345
|
rt300@34
|
346 // next rotations...
|
rt300@34
|
347 int lse2 = rotr(newe, P-direction-1,P);
|
rt300@34
|
348 entryPoint = entryPoint ^ lse2;
|
rt300@34
|
349 direction = (direction + newd + 1) % P;
|
rt300@34
|
350
|
rt300@34
|
351 // now build up h from subindices
|
rt300@34
|
352 h += subindex << (blev*P);
|
rt300@35
|
353 //cout << "------------\n";
|
rt300@34
|
354 }
|
rt300@34
|
355
|
rt300@34
|
356 return h;
|
rt300@34
|
357 }
|
rt300@34
|
358 //--------------------------------------------------------------------
|