EMAN2
reconstructor_tools.cpp
Go to the documentation of this file.
1/*
2 * Author: David Woolford, 07/25/2007 (woolford@bcm.edu)
3 * Copyright (c) 2000-2007 Baylor College of Medicine
4 *
5 * This software is issued under a joint BSD/GNU license. You may use the
6 * source code in this file under either license. However, note that the
7 * complete EMAN2 and SPARX software packages have some GPL dependencies,
8 * so you are responsible for compliance with the licenses of these packages
9 * if you opt to use BSD licensing. The warranty disclaimer below holds
10 * in either instance.
11 *
12 * This complete copyright notice must be included in any revised version of the
13 * source code. Additional authorship citations may be added, but existing
14 * author citations must be preserved.
15 *
16 * This program is free software; you can redistribute it and/or modify
17 * it under the terms of the GNU General Public License as published by
18 * the Free Software Foundation; either version 2 of the License, or
19 * (at your option) any later version.
20 *
21 * This program is distributed in the hope that it will be useful,
22 * but WITHOUT ANY WARRANTY; without even the implied warranty of
23 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24 * GNU General Public License for more details.
25 *
26 * You should have received a copy of the GNU General Public License
27 * along with this program; if not, write to the Free Software
28 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
29 *
30 * */
31
32//#define DEBUG_POINT 1
33
34#include <cstring>
35#include <math.h>
36#include <gsl/gsl_sf_bessel.h>
37#include "reconstructor_tools.h"
38
39
40using namespace EMAN;
41
42const string FourierInserter3DMode1::NAME = "nearest_neighbor";
43const string FourierInserter3DMode2l::NAME = "trilinear";
44const string FourierInserter3DMode2::NAME = "gauss_2";
45const string FourierInserter3DMode3::NAME = "gauss_3";
46//const string FourierInserter3DMode4::NAME = "gauss_4";
47const string FourierInserter3DMode5::NAME = "gauss_5";
48const string FourierInserter3DMode6::NAME = "gauss_var";
49const string FourierInserter3DMode7::NAME = "gridding_5";
50const string FourierInserter3DMode11::NAME = "gridding_7";
51const string FourierInserter3DMode8::NAME = "experimental";
52const string FourierInserter3DMode9::NAME = "kaiser_bessel";
53const string FourierInserter3DMode10::NAME = "kaiser_bessel_derived";
54
55template <> Factory < FourierPixelInserter3D >::Factory()
56{
57 force_add<FourierInserter3DMode1>();
58 force_add<FourierInserter3DMode2l>();
59 force_add<FourierInserter3DMode2>();
60 force_add<FourierInserter3DMode3>();
61// force_add<FourierInserter3DMode4>();
62 force_add<FourierInserter3DMode5>();
63 force_add<FourierInserter3DMode6>();
64 force_add<FourierInserter3DMode7>();
65 force_add<FourierInserter3DMode8>();
66// force_add(&FourierInserter3DMode8::NEW);
67 force_add<FourierInserter3DMode9>();
68 force_add<FourierInserter3DMode10>();
69 force_add<FourierInserter3DMode11>();
70}
71
72
73
74void FourierPixelInserter3D::init()
75{
76 if ( params.has_key("data") )
77 {
78 data = params["data"];
79 if ( data == 0 )
80 throw NotExistingObjectException("data", "error the data pointer was 0 in FourierPixelInserter3D::init");
81 }
82 else throw NotExistingObjectException("data", "the data pointer was not defined in FourierPixelInserter3D::init");
83
84 if ( params.has_key("norm"))
85 {
86 norm = params["norm"];
87 if ( norm == 0 )
88 throw NotExistingObjectException("norm", "error the norm pointer was 0 in FourierPixelInserter3D::init");
89 }
90 else throw NotExistingObjectException("norm", "the norm pointer was not defined in FourierPixelInserter3D::init");
91
92 nx=data->get_xsize();
93 ny=data->get_ysize();
94 nz=data->get_zsize();
95 nxyz=(size_t)nx*ny*nz;
96 nx2=nx/2-1;
97 ny2=ny/2;
98 nz2=nz/2;
99
100 if (data->has_attr("subvolume_x0") && data->has_attr("subvolume_full_nx")) {
101 subx0=data->get_attr("subvolume_x0");
102 suby0=data->get_attr("subvolume_y0");
103 subz0=data->get_attr("subvolume_z0");
104 fullnx=data->get_attr("subvolume_full_nx");
105 fullny=data->get_attr("subvolume_full_ny");
106 fullnz=data->get_attr("subvolume_full_nz");
107 }
108 else {
109 subx0=suby0=subz0=-1;
110 }
111}
112
113bool FourierInserter3DMode1::insert_pixel(const float& xx, const float& yy, const float& zz, const std::complex<float> dt, const float& weight)
114{
115 int x0 = (int) floor(xx + 0.5f);
116 int y0 = (int) floor(yy + 0.5f);
117 int z0 = (int) floor(zz + 0.5f);
118
119 size_t off;
120 if (subx0<0) off=data->add_complex_at(x0,y0,z0,dt*weight);
121 else off=data->add_complex_at(x0,y0,z0,subx0,suby0,subz0,fullnx,fullny,fullnz,dt*weight);
122 if (static_cast<int>(off)!=nxyz) norm[off/2]+=weight;
123 else return false;
124
125 return true;
126}
127
128bool FourierInserter3DMode2::insert_pixel(const float& xx, const float& yy, const float& zz, const std::complex<float> dt,const float& weight)
129{
130 int x0 = (int) floor(xx);
131 int y0 = (int) floor(yy);
132 int z0 = (int) floor(zz);
133
134 // note that subnx differs in the inserters. In the reconstructors it subx0 is 0 for the full volume. Here it is -1
135 if (subx0<0) { // normal full reconstruction
136 if (x0<-nx2-1 || y0<-ny2-1 || z0<-nz2-1 || x0>nx2 || y0>ny2 || z0>nz2 ) return false;
137
138 int x1=x0+1;
139 int y1=y0+1;
140 int z1=z0+1;
141// if (x0<-nx2) x0=-nx2;
142// if (x1>nx2) x1=nx2;
143// if (y0<-ny2) y0=-ny2;
144// if (y1>ny2) y1=ny2;
145// if (z0<-nz2) z0=-nz2;
146// if (z1>nz2) z1=nz2;
147
148// float h=2.0/((1.0+pow(Util::hypot3sq(xx,yy,zz),.5))*EMConsts::I2G);
149 float h=1.0f/EMConsts::I2G;
150 //size_t idx;
151 float r, gg;
152// int pc=0;
153 for (int k = z0 ; k <= z1; k++) {
154 for (int j = y0 ; j <= y1; j++) {
155 for (int i = x0; i <= x1; i ++) {
156 r = Util::hypot3sq((float) i - xx, j - yy, k - zz);
157 gg = Util::fast_exp(-r *h)*weight;
158
159 size_t off;
160 off=data->add_complex_at_fast(i,j,k,dt*gg);
161 if (off!=nxyz) norm[off/2]+=gg;
162// if (off!=nxyz) norm[off/2]+=weight; // experiment 6/1/20, use true Gaussian kernel, not just weight
163 }
164 }
165 }
166 return true;
167 }
168 else { // for subvolumes, not optimized yet
169 //size_t idx;
170 float r, gg;
171 int pc=0;
172 for (int k = z0 ; k <= z0 + 1; k++) {
173 for (int j = y0 ; j <= y0 + 1; j++) {
174 for (int i = x0; i <= x0 + 1; i ++) {
175 r = Util::hypot3sq((float) i - xx, j - yy, k - zz);
176 gg = Util::fast_exp(-r / EMConsts::I2G)*weight;
177
178 size_t off;
179 if (subx0<0) off=data->add_complex_at(i,j,k,dt*gg);
180 else off=data->add_complex_at(i,j,k,subx0,suby0,subz0,fullnx,fullny,fullnz,dt*gg);
181 if (static_cast<int>(off)!=nxyz) { norm[off/2]+=gg; pc+=1; }
182 }
183 }
184 }
185
186 if (pc>0) return true;
187 return false;
188 }
189}
190
191bool FourierInserter3DMode2l::insert_pixel(const float& xx, const float& yy, const float& zz, const std::complex<float> dt,const float& weight)
192{
193 int x0 = (int) floor(xx);
194 int y0 = (int) floor(yy);
195 int z0 = (int) floor(zz);
196
197 // note that subnx differs in the inserters. In the reconstructors it subx0 is 0 for the full volume. Here it is -1
198 if (subx0<0) { // normal full reconstruction
199 if (x0<-nx2-1 || y0<-ny2-1 || z0<-nz2-1 || x0>nx2 || y0>ny2 || z0>nz2 ) return false;
200
201 int x1=x0+1;
202 int y1=y0+1;
203 int z1=z0+1;
204// if (x0<-nx2) x0=-nx2;
205// if (x1>nx2) x1=nx2;
206// if (y0<-ny2) y0=-ny2;
207// if (y1>ny2) y1=ny2;
208// if (z0<-nz2) z0=-nz2;
209// if (z1>nz2) z1=nz2;
210
211// float h=2.0/((1.0+pow(Util::hypot3sq(xx,yy,zz),.5))*EMConsts::I2G);
212 float h=1.0f/EMConsts::I2G;
213 //size_t idx;
214 float r, gg;
215// int pc=0;
216 for (int k = z0 ; k <= z1; k++) {
217 for (int j = y0 ; j <= y1; j++) {
218 for (int i = x0; i <= x1; i ++) {
219// r = Util::hypot3sq((float) i - xx, j - yy, k - zz);
220// gg = Util::fast_exp(-r *h)*weight;
221 gg=(1.0-fabs(i-xx))*(1.0-fabs(j-yy))*(1.0-fabs(k-zz))*weight;
222
223 size_t off;
224 off=data->add_complex_at_fast(i,j,k,dt*gg);
225 if (off!=nxyz) norm[off/2]+=gg;
226// if (off!=nxyz) norm[off/2]+=weight; // experiment 6/1/20, use true Gaussian kernel, not just weight
227#ifdef DEBUG_POINT
228 if (k==23 && j==0 && i==113) {
229 std::complex<float> pv = data->get_complex_at(113,0,23);
230 float pnv = norm[off/2];
231 printf("insert: %1.4g\t%1.4g\t%1.4g\t%1.4g\t%1.4g\t%1.4g\n",pv.real()/pnv,pv.imag()/pnv,pnv,dt.real(),dt.imag(),gg);
232 }
233#endif
234
235 }
236 }
237 }
238 return true;
239 }
240 else { // for subvolumes, not optimized yet
241 //size_t idx;
242 float r, gg;
243 int pc=0;
244 for (int k = z0 ; k <= z0 + 1; k++) {
245 for (int j = y0 ; j <= y0 + 1; j++) {
246 for (int i = x0; i <= x0 + 1; i ++) {
247// r = Util::hypot3sq((float) i - xx, j - yy, k - zz);
248// gg = Util::fast_exp(-r / EMConsts::I2G)*weight;
249 gg=(1.0-fabs(i-xx))*(1.0-fabs(j-yy))*(1.0-fabs(k-zz))*weight;
250
251 size_t off;
252 if (subx0<0) off=data->add_complex_at(i,j,k,dt*gg);
253 else off=data->add_complex_at(i,j,k,subx0,suby0,subz0,fullnx,fullny,fullnz,dt*gg);
254 if (static_cast<int>(off)!=nxyz) { norm[off/2]+=gg; pc+=1; }
255 }
256 }
257 }
258
259 if (pc>0) return true;
260 return false;
261 }
262}
263
264
265bool FourierInserter3DMode3::insert_pixel(const float& xx, const float& yy, const float& zz, const std::complex<float> dt,const float& weight)
266{
267 int x0 = (int) floor(xx-.5);
268 int y0 = (int) floor(yy-.5);
269 int z0 = (int) floor(zz-.5);
270
271 if (subx0<0) { // normal full reconstruction
272 if (x0<-nx2-2 || y0<-ny2-2 || z0<-nz2-2 || x0>nx2+1 || y0>ny2+1 || z0>nz2+1 ) return false;
273
274 // no error checking on add_complex_fast, so we need to be careful here
275 int x1=x0+2;
276 int y1=y0+2;
277 int z1=z0+2;
278 if (x0<-nx2) x0=-nx2;
279 if (x1>nx2) x1=nx2;
280 if (y0<-ny2) y0=-ny2;
281 if (y1>ny2) y1=ny2;
282 if (z0<-nz2) z0=-nz2;
283 if (z1>nz2) z1=nz2;
284
285// float h=2.0/((1.0+pow(Util::hypot3sq(xx,yy,zz),.5))*EMConsts::I2G);
286// float h=2.0/EMConsts::I3G;
287 float h=32.0f/((8.0f+Util::hypot3(xx,yy,zz))*EMConsts::I3G);
288// float w=weight;
289 float w=weight/(1.0f+6.0f*Util::fast_exp(-h)+12*Util::fast_exp(-h*2.0f)+8*Util::fast_exp(-h*3.0f)); // approx normalization so higer radii aren't upweighted relative to lower due to wider Gaussian
290 //size_t idx;
291 float r, gg;
292// int pc=0;
293 for (int k = z0 ; k <= z1; k++) {
294 for (int j = y0 ; j <= y1; j++) {
295 for (int i = x0; i <= x1; i ++) {
296 r = Util::hypot3sq((float) i - xx, (float)j - yy, (float)k - zz);
297// gg=weight;
298 gg = Util::fast_exp(-r *h)*w;
299// gg = Util::fast_exp(-r / EMConsts::I2G)*weight;
300// gg = sqrt(Util::fast_exp(-r / EMConsts::I2G))*weight;
301
302 size_t off;
303 off=data->add_complex_at_fast(i,j,k,dt*gg);
304 norm[off/2]+=gg;
305 }
306 }
307 }
308 return true;
309 }
310 printf("region writing not supported in mode 3\n");
311 return false;
312}
313
314
315bool FourierInserter3DMode5::insert_pixel(const float& xx, const float& yy, const float& zz, const std::complex<float> dt,const float& weight)
316{
317 int x0 = (int) floor(xx-2.5);
318 int y0 = (int) floor(yy-2.5);
319 int z0 = (int) floor(zz-2.5);
320
321 if (subx0<0) { // normal full reconstruction
322 if (x0<-nx2-4 || y0<-ny2-4 || z0<-nz2-4 || x0>nx2+3 || y0>ny2+3 || z0>nz2+3 ) return false;
323
324 // no error checking on add_complex_fast, so we need to be careful here
325 int x1=x0+5;
326 int y1=y0+5;
327 int z1=z0+5;
328 if (x0<-nx2) x0=-nx2;
329 if (x1>nx2) x1=nx2;
330 if (y0<-ny2) y0=-ny2;
331 if (y1>ny2) y1=ny2;
332 if (z0<-nz2) z0=-nz2;
333 if (z1>nz2) z1=nz2;
334
335// float h=2.0/((1.0+pow(Util::hypot3sq(xx,yy,zz),.5))*EMConsts::I2G);
336 float h=1.0f/EMConsts::I5G;
337 float w=weight;
338
339 // Not sure exactly what this was doing? Using wider Gaussian at high radius?
340// float h=32.0f/((8.0f+Util::hypot3(xx,yy,zz))*EMConsts::I3G);
341// float w=weight/(1.0f+6.0f*Util::fast_exp(-h)+12*Util::fast_exp(-h*2.0f)+8*Util::fast_exp(-h*3.0f)+
342// 6.0f*Util::fast_exp(-h*4.0f)+24.0f*Util::fast_exp(-h*5.0f)+24.0f*Util::fast_exp(-h*6.0f)+12.0f*Util::fast_exp(-h*8.0f)+
343// 24.0f*Util::fast_exp(-h*9.0f)+8.0f*Util::fast_exp(-h*12.0f)); // approx normalization so higer radii aren't upweighted relative to lower due to wider Gaussian
344 //size_t idx;
345 float r, gg;
346// int pc=0;
347 for (int k = z0 ; k <= z1; k++) {
348 for (int j = y0 ; j <= y1; j++) {
349 for (int i = x0; i <= x1; i ++) {
350 r = Util::hypot3sq((float) i - xx, j - yy, k - zz);
351// gg=weight;
352 gg = Util::fast_exp(-r *h);
353// gg = Util::fast_exp(-r / EMConsts::I2G)*weight;
354// gg = sqrt(Util::fast_exp(-r / EMConsts::I2G))*weight;
355
356 size_t off;
357 off=data->add_complex_at_fast(i,j,k,dt*gg*w);
358 norm[off/2]+=gg*w; // This would use a Gaussian WEIGHT with square kernel
359// norm[off/2]+=w; // This would use a Gaussian KERNEL rather than WEIGHT
360
361#ifdef RECONDEBUG
362 std::complex<double> v1=dt*gg*w,v2=gg*w;
363
364 if (k<5 && j<5&& i<5&& k>=0 && j>=0 && i>=0) {
365 int idx=i*2+j*10+k*50;
366 ddata[idx]+=v1.real();
367 ddata[idx+1]+=v1.imag();
368 dnorm[idx]+=v2.real();
369 dnorm[idx+1]+=v2.imag();
370 }
371#endif
372 }
373 }
374 }
375 return true;
376 }
377 printf("region writing not supported in mode 5\n");
378 return false;
379}
380
381
382bool FourierInserter3DMode6::insert_pixel(const float& xx, const float& yy, const float& zz, const std::complex<float> dt,const float& weight)
383{
384 int x0 = (int) floor(xx-2.5);
385 int y0 = (int) floor(yy-2.5);
386 int z0 = (int) floor(zz-2.5);
387
388 if (subx0<0) { // normal full reconstruction
389 if (x0<-nx2-4 || y0<-ny2-4 || z0<-nz2-4 || x0>nx2+3 || y0>ny2+3 || z0>nz2+3 ) return false;
390
391 // no error checking on add_complex_fast, so we need to be careful here
392 int x1=x0+5;
393 int y1=y0+5;
394 int z1=z0+5;
395 if (x0<-nx2) x0=-nx2;
396 if (x1>nx2) x1=nx2;
397 if (y0<-ny2) y0=-ny2;
398 if (y1>ny2) y1=ny2;
399 if (z0<-nz2) z0=-nz2;
400 if (z1>nz2) z1=nz2;
401
402// float h=2.0/((1.0+pow(Util::hypot3sq(xx,yy,zz),.5))*EMConsts::I2G);
403// float h=1.0f/EMConsts::I5G;
404// float h=1.0f/(Util::hypot3sq(xx/nx2,yy/ny2,zz/nz2)*EMConsts::I5G*2.0+.1); // gaussian kernel is used as a weight not a kernel. We increase radius of integration with resolution. Changed away from this on 11/10/17
405 float h=1.0f/((Util::hypot3sq(xx,yy,zz)/4000)+.15); // gaussian kernel is used as a weight not a kernel. We increase radius of integration with resolution. Changed away from this on 11/10/17
406 h=h<0.1?0.1:h; // new formula has h from 0.2 to 5
407// if (yy==0 &&zz==0) printf("%0.0f\t%0.0f\t%0.0f\t%g\n",xx,yy,zz,h);
408// printf("%1.0f\t%1.0f\t%1.0f\t%1.4f\n",xx,yy,zz,h);
409 float w=weight;
410
411// float h=32.0f/((8.0f+Util::hypot3(xx,yy,zz))*EMConsts::I3G);
412// float w=weight/(1.0f+6.0f*Util::fast_exp(-h)+12*Util::fast_exp(-h*2.0f)+8*Util::fast_exp(-h*3.0f)+
413// 6.0f*Util::fast_exp(-h*4.0f)+24.0f*Util::fast_exp(-h*5.0f)+24.0f*Util::fast_exp(-h*6.0f)+12.0f*Util::fast_exp(-h*8.0f)+
414// 24.0f*Util::fast_exp(-h*9.0f)+8.0f*Util::fast_exp(-h*12.0f)); // approx normalization so higer radii aren't upweighted relative to lower due to wider Gaussian
415 //size_t idx;
416 float r, gg;
417// int pc=0;
418 for (int k = z0 ; k <= z1; k++) {
419 for (int j = y0 ; j <= y1; j++) {
420 for (int i = x0; i <= x1; i ++) {
421 r = Util::hypot3sq((float) i - xx, j - yy, k - zz);
422// gg=weight;
423 gg = Util::fast_exp(-r *h);
424// if (gg<.00001) continue; // skip tiny weights for speed
425// gg = Util::fast_exp(-r / EMConsts::I2G)*weight;
426// gg = sqrt(Util::fast_exp(-r / EMConsts::I2G))*weight;
427
428 size_t off;
429 off=data->add_complex_at_fast(i,j,k,dt*gg*w);
430// if (i==16 && j==7 && k==4) printf("%g\t%g\t%g\t%g\t%g\n",dt.real(),dt.imag(),gg*w,gg,w);
431 norm[off/2]+=gg*w; // This would use a Gaussian WEIGHT with square kernel
432// norm[off/2]+=w; // This would use a Gaussian KERNEL rather than WEIGHT
433
434#ifdef RECONDEBUG
435 std::complex<double> v1=dt*gg*w,v2=gg*w;
436
437 if (k<5 && j<5&& i<5&& k>=0 && j>=0 && i>=0) {
438 int idx=i*2+j*10+k*50;
439 ddata[idx]+=v1.real();
440 ddata[idx+1]+=v1.imag();
441 dnorm[idx]+=v2.real();
442 dnorm[idx+1]+=v2.imag();
443 }
444#endif
445 }
446 }
447 }
448 return true;
449 }
450 printf("region writing not supported in mode 5\n");
451 return false;
452}
453
454// Kernel determined by examples/kernel_opt
455const float FourierInserter3DMode7::kernel[9][9][9] = {
4560.3293228,0.2999721,0.2230405,0.1263839,0.0412195,-0.0116440,-0.0287601,-0.0214417,0.0000000,
4570.2999721,0.2731577,0.2028874,0.1146381,0.0369511,-0.0111675,-0.0266030,-0.0197382,0.0000000,
4580.2230405,0.2028874,0.1501083,0.0839287,0.0258589,-0.0098285,-0.0208839,-0.0152399,0.0000000,
4590.1263839,0.1146381,0.0839287,0.0455768,0.0122079,-0.0078769,-0.0135036,-0.0094886,0.0000000,
4600.0412195,0.0369511,0.0258589,0.0122079,0.0007022,-0.0056645,-0.0066442,-0.0042390,0.0000000,
461-0.0116440,-0.0111675,-0.0098285,-0.0078769,-0.0056645,-0.0035595,-0.0018569,-0.0007100,0.0000000,
462-0.0287601,-0.0266030,-0.0208839,-0.0135036,-0.0066442,-0.0018569,0.0004320,0.0008101,0.0000000,
463-0.0214417,-0.0197382,-0.0152398,-0.0094886,-0.0042390,-0.0007100,0.0008101,0.0008532,0.0000000,
4640.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
4650.2999721,0.2731577,0.2028874,0.1146381,0.0369511,-0.0111675,-0.0266030,-0.0197382,0.0000000,
4660.2731577,0.2486666,0.1844966,0.1039437,0.0330965,-0.0106899,-0.0246017,-0.0181663,0.0000000,
4670.2028874,0.1844966,0.1363422,0.0759894,0.0230832,-0.0093549,-0.0192968,-0.0140161,0.0000000,
4680.1146381,0.1039437,0.0759894,0.0410967,0.0107719,-0.0074293,-0.0124545,-0.0087120,0.0000000,
4690.0369511,0.0330965,0.0230832,0.0107719,0.0004174,-0.0052781,-0.0061010,-0.0038740,0.0000000,
470-0.0111675,-0.0106899,-0.0093549,-0.0074294,-0.0052781,-0.0032679,-0.0016753,-0.0006270,0.0000000,
471-0.0266030,-0.0246017,-0.0192968,-0.0124545,-0.0061010,-0.0016753,0.0004302,0.0007650,0.0000000,
472-0.0197382,-0.0181663,-0.0140161,-0.0087120,-0.0038740,-0.0006270,0.0007650,0.0007953,0.0000000,
4730.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
4740.2230405,0.2028874,0.1501083,0.0839287,0.0258589,-0.0098285,-0.0208839,-0.0152399,0.0000000,
4750.2028874,0.1844966,0.1363422,0.0759894,0.0230832,-0.0093549,-0.0192968,-0.0140161,0.0000000,
4760.1501083,0.1363422,0.1003225,0.0552536,0.0158836,-0.0080490,-0.0150930,-0.0107868,0.0000000,
4770.0839287,0.0759894,0.0552535,0.0294215,0.0070640,-0.0062153,-0.0096798,-0.0066651,0.0000000,
4780.0258589,0.0230832,0.0158836,0.0070640,-0.0002933,-0.0042462,-0.0046694,-0.0029154,0.0000000,
479-0.0098285,-0.0093549,-0.0080490,-0.0062153,-0.0042462,-0.0024995,-0.0012019,-0.0004126,0.0000000,
480-0.0208839,-0.0192968,-0.0150930,-0.0096798,-0.0046694,-0.0012019,0.0004197,0.0006424,0.0000000,
481-0.0152398,-0.0140161,-0.0107868,-0.0066651,-0.0029154,-0.0004126,0.0006424,0.0006405,0.0000000,
4820.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
4830.1263839,0.1146381,0.0839287,0.0455768,0.0122079,-0.0078769,-0.0135036,-0.0094886,0.0000000,
4840.1146381,0.1039437,0.0759894,0.0410967,0.0107719,-0.0074294,-0.0124545,-0.0087120,0.0000000,
4850.0839287,0.0759894,0.0552535,0.0294215,0.0070640,-0.0062153,-0.0096798,-0.0066651,0.0000000,
4860.0455767,0.0410967,0.0294215,0.0149533,0.0025724,-0.0045680,-0.0061197,-0.0040604,0.0000000,
4870.0122078,0.0107719,0.0070640,0.0025724,-0.0010797,-0.0028931,-0.0028473,-0.0017048,0.0000000,
488-0.0078770,-0.0074294,-0.0062153,-0.0045680,-0.0028931,-0.0015217,-0.0006148,-0.0001526,0.0000000,
489-0.0135036,-0.0124545,-0.0096798,-0.0061197,-0.0028473,-0.0006148,0.0003888,0.0004754,0.0000000,
490-0.0094886,-0.0087120,-0.0066651,-0.0040604,-0.0017048,-0.0001527,0.0004754,0.0004373,0.0000000,
4910.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
4920.0412195,0.0369511,0.0258589,0.0122079,0.0007022,-0.0056645,-0.0066442,-0.0042390,0.0000000,
4930.0369511,0.0330965,0.0230832,0.0107719,0.0004174,-0.0052781,-0.0061010,-0.0038740,0.0000000,
4940.0258589,0.0230832,0.0158836,0.0070640,-0.0002933,-0.0042462,-0.0046694,-0.0029154,0.0000000,
4950.0122078,0.0107719,0.0070640,0.0025724,-0.0010797,-0.0028931,-0.0028473,-0.0017048,0.0000000,
4960.0007022,0.0004174,-0.0002933,-0.0010797,-0.0015774,-0.0015967,-0.0011988,-0.0006268,0.0000000,
497-0.0056645,-0.0052781,-0.0042463,-0.0028931,-0.0015967,-0.0006377,-0.0001117,0.0000590,0.0000000,
498-0.0066442,-0.0061010,-0.0046694,-0.0028473,-0.0011988,-0.0001117,0.0003294,0.0003044,0.0000000,
499-0.0042390,-0.0038740,-0.0029154,-0.0017048,-0.0006268,0.0000590,0.0003044,0.0002424,0.0000000,
5000.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
501-0.0116440,-0.0111675,-0.0098285,-0.0078769,-0.0056645,-0.0035595,-0.0018569,-0.0007100,0.0000000,
502-0.0111675,-0.0106899,-0.0093549,-0.0074294,-0.0052781,-0.0032679,-0.0016753,-0.0006270,0.0000000,
503-0.0098285,-0.0093549,-0.0080490,-0.0062153,-0.0042463,-0.0024995,-0.0012019,-0.0004126,0.0000000,
504-0.0078770,-0.0074294,-0.0062153,-0.0045680,-0.0028931,-0.0015217,-0.0006148,-0.0001526,0.0000000,
505-0.0056645,-0.0052781,-0.0042463,-0.0028931,-0.0015967,-0.0006377,-0.0001117,0.0000590,0.0000000,
506-0.0035595,-0.0032679,-0.0024995,-0.0015217,-0.0006377,-0.0000555,0.0001796,0.0001647,0.0000000,
507-0.0018569,-0.0016753,-0.0012019,-0.0006148,-0.0001117,0.0001796,0.0002446,0.0001630,0.0000000,
508-0.0007100,-0.0006270,-0.0004126,-0.0001527,0.0000590,0.0001647,0.0001630,0.0000978,0.0000000,
5090.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
510-0.0287600,-0.0266030,-0.0208839,-0.0135036,-0.0066442,-0.0018569,0.0004320,0.0008101,0.0000000,
511-0.0266030,-0.0246017,-0.0192968,-0.0124545,-0.0061010,-0.0016753,0.0004302,0.0007650,0.0000000,
512-0.0208839,-0.0192968,-0.0150930,-0.0096798,-0.0046694,-0.0012019,0.0004197,0.0006424,0.0000000,
513-0.0135036,-0.0124545,-0.0096798,-0.0061197,-0.0028473,-0.0006148,0.0003888,0.0004754,0.0000000,
514-0.0066442,-0.0061010,-0.0046694,-0.0028473,-0.0011988,-0.0001117,0.0003294,0.0003044,0.0000000,
515-0.0018569,-0.0016753,-0.0012019,-0.0006148,-0.0001117,0.0001796,0.0002446,0.0001630,0.0000000,
5160.0004320,0.0004302,0.0004197,0.0003888,0.0003294,0.0002446,0.0001503,0.0000679,0.0000000,
5170.0008101,0.0007650,0.0006424,0.0004754,0.0003044,0.0001630,0.0000679,0.0000181,0.0000000,
5180.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
519-0.0214417,-0.0197382,-0.0152398,-0.0094886,-0.0042390,-0.0007100,0.0008101,0.0008532,0.0000000,
520-0.0197382,-0.0181663,-0.0140161,-0.0087120,-0.0038740,-0.0006270,0.0007650,0.0007953,0.0000000,
521-0.0152398,-0.0140161,-0.0107868,-0.0066651,-0.0029154,-0.0004126,0.0006424,0.0006405,0.0000000,
522-0.0094886,-0.0087120,-0.0066651,-0.0040604,-0.0017048,-0.0001527,0.0004754,0.0004373,0.0000000,
523-0.0042390,-0.0038740,-0.0029154,-0.0017048,-0.0006268,0.0000590,0.0003044,0.0002424,0.0000000,
524-0.0007100,-0.0006270,-0.0004126,-0.0001527,0.0000590,0.0001647,0.0001630,0.0000978,0.0000000,
5250.0008101,0.0007650,0.0006424,0.0004754,0.0003044,0.0001630,0.0000679,0.0000181,0.0000000,
5260.0008532,0.0007953,0.0006405,0.0004373,0.0002424,0.0000978,0.0000181,-0.0000082,0.0000000,
5270.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
5280.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
5290.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
5300.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
5310.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
5320.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
5330.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
5340.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
5350.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
5360.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000
537};
538
539bool FourierInserter3DMode7::insert_pixel(const float& xx, const float& yy, const float& zz, const std::complex<float> dt,const float& weight)
540{
541 int x0 = Util::fast_floor(xx-1.5);
542 int y0 = Util::fast_floor(yy-1.5);
543 int z0 = Util::fast_floor(zz-1.5);
544
545 if (subx0<0) { // normal full reconstruction
546 if (x0<-nx2-4 || y0<-ny2-4 || z0<-nz2-4 || x0>nx2+3 || y0>ny2+3 || z0>nz2+3 ) return false;
547
548 // no error checking on add_complex_fast, so we need to be careful here
549 int x1=x0+4;
550 int y1=y0+4;
551 int z1=z0+4;
552 if (x0<-nx2) x0=-nx2;
553 if (x1>nx2) x1=nx2;
554 if (y0<-ny2) y0=-ny2;
555 if (y1>ny2) y1=ny2;
556 if (z0<-nz2) z0=-nz2;
557 if (z1>nz2) z1=nz2;
558
559 float w=weight;
560
561 float r, gg;
562 for (int k = z0 ; k <= z1; k++) {
563 for (int j = y0 ; j <= y1; j++) {
564 for (int i = x0; i <= x1; i ++) {
565 gg=FourierInserter3DMode7::kernel[abs(Util::fast_floor((i-xx)*3.0f+0.5))][abs(Util::fast_floor((j-yy)*3.0f+0.5))][abs(Util::fast_floor((k-zz)*3.0f+0.5))];
566
567 size_t off;
568 off=data->add_complex_at_fast(i,j,k,dt*gg*w);
569 norm[off/2]+=w;
570// if (i==67&&j==19&&k==1) printf("%1.1f %1.1f %1.1f\t%d %d %d\t%d %d %d\t%f\t%f\t%f\t%f\t%f\t%f\n",xx,yy,zz,i,j,k,abs(Util::fast_floor((i-xx)*3.0f+0.5)),abs(Util::fast_floor((j-yy)*3.0f+0.5)),abs(Util::fast_floor((k-zz)*3.0f+0.5)),gg,w,dt.real(),dt.imag(),norm[off/2],data->get_value_at(off));
571
572 }
573 }
574 }
575 return true;
576 }
577 printf("region writing not supported in mode \n");
578 return false;
579}
580
581
582void FourierInserter3DMode8::init()
583{
584 FourierPixelInserter3D::init();
585// int P = (int)((1.0+0.25)*nx+1);
586// float r = (float)(nx+1)/(float)P;
587// mFreqCutoff = 2;
588// mDFreq = 0.2f;
589// if (W != 0) delete [] W;
590// W = Util::getBaldwinGridWeights(mFreqCutoff, (float)P, r,mDFreq,0.5f,0.2f);
591
592}
593
594bool FourierInserter3DMode8::insert_pixel(const float& xx, const float& yy, const float& zz, const std::complex<float> dt,const float& weight)
595{
596 static FILE *out400=NULL,*out862=NULL,*out962,*out872,*out4093=NULL;
597
598 int x0 = (int) floor(xx);
599 int y0 = (int) floor(yy);
600 int z0 = (int) floor(zz);
601
602 if (out400==NULL) {
603 out400=fopen("pxl4_0_0.txt","w");
604 out862=fopen("pxl8_6_2.txt","w");
605 out962=fopen("pxl9_6_2.txt","w");
606 out872=fopen("pxl8_7_2.txt","w");
607 out4093=fopen("pxl40_9_3.txt","w");
608 }
609
610 // note that subnx differs in the inserters. In the reconstructors it subx0 is 0 for the full volume. Here it is -1
611 if (subx0<0) { // normal full reconstruction
612 if (x0<-nx2-1 || y0<-ny2-1 || z0<-nz2-1 || x0>nx2 || y0>ny2 || z0>nz2 ) return false;
613
614 int x1=x0+1;
615 int y1=y0+1;
616 int z1=z0+1;
617// if (x0<-nx2) x0=-nx2;
618// if (x1>nx2) x1=nx2;
619// if (y0<-ny2) y0=-ny2;
620// if (y1>ny2) y1=ny2;
621// if (z0<-nz2) z0=-nz2;
622// if (z1>nz2) z1=nz2;
623
624// float h=2.0/((1.0+pow(Util::hypot3sq(xx,yy,zz),.5))*EMConsts::I2G);
625 float h=1.0f/EMConsts::I2G;
626 //size_t idx;
627 float r, gg;
628// int pc=0;
629 for (int k = z0 ; k <= z1; k++) {
630 for (int j = y0 ; j <= y1; j++) {
631 for (int i = x0; i <= x1; i ++) {
632 r = Util::hypot3sq((float) i - xx, j - yy, k - zz);
633// gg=weight;
634 gg = Util::fast_exp(-r *h)*weight;
635// gg = Util::fast_exp(-r / EMConsts::I2G)*weight;
636// gg = sqrt(Util::fast_exp(-r / EMConsts::I2G))*weight;
637
638 size_t off;
639 off=data->add_complex_at_fast(i,j,k,dt*gg);
640// off=data->add_complex_at(i,j,k,dt*gg);
641 if (off!=nxyz) norm[off/2]+=gg;
642
643 if (i==4&&j==0&&k==0) { fprintf(out400,"%1.4f\t%1.4f\t%1.4f\t%1.4f\t%1.4f\n",dt.real(),dt.imag(),gg,std::abs(dt),std::arg(dt)); fflush(out400); }
644 if (i==8&&j==6&&k==2) { fprintf(out862,"%1.4f\t%1.4f\t%1.4f\t%1.4f\t%1.4f\n",dt.real(),dt.imag(),gg,std::abs(dt),std::arg(dt)); fflush(out862); }
645 if (i==9&&j==6&&k==2) { fprintf(out962,"%1.4f\t%1.4f\t%1.4f\t%1.4f\t%1.4f\n",dt.real(),dt.imag(),gg,std::abs(dt),std::arg(dt)); fflush(out962); }
646 if (i==8&&j==7&&k==2) { fprintf(out872,"%1.4f\t%1.4f\t%1.4f\t%1.4f\t%1.4f\n",dt.real(),dt.imag(),gg,std::abs(dt),std::arg(dt)); fflush(out872); }
647 if (i==40&&j==9&&k==3) { fprintf(out4093,"%1.4f\t%1.4f\t%1.4f\t%1.4f\t%1.4f\n",dt.real(),dt.imag(),gg,std::abs(dt),std::arg(dt)); fflush(out4093); }
648 }
649 }
650 }
651 return true;
652 }
653 else { // for subvolumes, not optimized yet
654 //size_t idx;
655 float r, gg;
656 int pc=0;
657 for (int k = z0 ; k <= z0 + 1; k++) {
658 for (int j = y0 ; j <= y0 + 1; j++) {
659 for (int i = x0; i <= x0 + 1; i ++) {
660 r = Util::hypot3sq((float) i - xx, j - yy, k - zz);
661 gg = Util::fast_exp(-r / EMConsts::I2G)*weight;
662
663 size_t off;
664 if (subx0<0) off=data->add_complex_at(i,j,k,dt*gg);
665 else off=data->add_complex_at(i,j,k,subx0,suby0,subz0,fullnx,fullny,fullnz,dt*gg);
666 if (static_cast<int>(off)!=nxyz) { norm[off/2]+=gg; pc+=1; }
667 }
668 }
669 }
670
671 if (pc>0) return true;
672 return false;
673 }
674// int x0 = (int) floor(qx);
675// int y0 = (int) floor(qy);
676// int z0 = (int) floor(qz);
677
678// int sizeW = (int)(1+2*mFreqCutoff/mDFreq);
679// int sizeWmid = sizeW/2;
680
681// for (int z = z0-mFreqCutoff; z < z0+mFreqCutoff; ++z){
682// for (int y = y0-mFreqCutoff; y < y0+mFreqCutoff; ++y){
683// for (int x = x0-mFreqCutoff; x < x0+mFreqCutoff; ++x){
684// if ( x < 0 || x >= nx ) continue;
685// if ( y < 0 || y >= ny ) continue;
686// if ( z < 0 || z >= nz ) continue;
687// float dist = (float)((x-x0)*(x-x0)+(y-y0)*(y-y0)+(z-z0)*(z-z0));
688// dist = sqrtf(dist);
689// // We enforce a spherical kernel
690// if ( dist > mFreqCutoff ) continue;
691// int idx = (int)(sizeWmid + dist/mDFreq);
692// if (idx >= sizeW) throw;
693// float residual = dist/mDFreq - (int)(dist/mDFreq);
694// if ( fabs(residual) > 1) throw;
695//
696// float factor = W[idx]*(1.0f-residual)+W[idx+1]*residual*weight;
697//
698// size_t k = z*nx*ny + y*nx + 2*x;
699//
700// // float c = Util::agauss(1, x-x0,y-y0,z-z0, EMConsts::I2G);
701// rdata[k] += fq[0]*factor;
702// rdata[k+1] += fq[1]*factor;
703//
704//
705// norm[k/2] += weight;
706//
707// }
708// }
709// }
710
711 return true;
712}
713
714bool FourierInserter3DMode9::insert_pixel(const float& xx, const float& yy, const float& zz, const std::complex<float> dt,const float& weight)
715{
716 int N = 4; // kernel width
717
718 int x0 = (int) floor(xx-N/2);
719 int y0 = (int) floor(yy-N/2);
720 int z0 = (int) floor(zz-N/2);
721
722 if (subx0<0) { // normal full reconstruction
723 if (x0<-nx2-7 || y0<-ny2-7 || z0<-nz2-7 || x0>nx2+6 || y0>ny2+6 || z0>nz2+6 ) return false;
724
725 // no error checking on add_complex_fast, so we need to be careful here
726 int x1=x0+N;
727 int y1=y0+N;
728 int z1=z0+N;
729
730 if (x0<-nx2) x0=-nx2;
731 if (x1>nx2) x1=nx2;
732 if (y0<-ny2) y0=-ny2;
733 if (y1>ny2) y1=ny2;
734 if (z0<-nz2) z0=-nz2;
735 if (z1>nz2) z1=nz2;
736
737 float w=weight;
738 float a=15.0;
739 float r, kb;
740
741 for (int k = z0 ; k <= z1; k++) {
742 for (int j = y0 ; j <= y1; j++) {
743 for (int i = x0; i <= x1; i ++) {
744 r = Util::hypot3sq((float) i - xx, j - yy, k - zz);
745 kb = gsl_sf_bessel_i0_scaled(M_PI * a * sqrt(1.0f - Util::square((r/(nx2-1))-1))) /
746 gsl_sf_bessel_i0_scaled(M_PI * a);
747 size_t off;
748 off = data->add_complex_at_fast(i,j,k,dt*kb*w);
749 norm[off/2]+=w;
750 }
751 }
752 }
753 return true;
754 }
755 printf("region writing not supported in mode 9\n");
756 return false;
757}
758
759// imprecise KBD kernel/window
760bool FourierInserter3DMode10::insert_pixel(const float& xx, const float& yy, const float& zz, const std::complex<float> dt,const float& weight)
761{
762 const int N = 8; // kernel width
763
764 int x0 = (int) floor(xx-N/2);
765 int y0 = (int) floor(yy-N/2);
766 int z0 = (int) floor(zz-N/2);
767
768 if (subx0<0) { // normal full reconstruction
769 if (x0<-nx2-7 || y0<-ny2-7 || z0<-nz2-7 || x0>nx2+6 || y0>ny2+6 || z0>nz2+6 ) return false;
770
771 // no error checking on add_complex_fast, so we need to be careful here
772 int x1=x0+N;
773 int y1=y0+N;
774 int z1=z0+N;
775
776 if (x0<-nx2) x0=-nx2;
777 if (x1>nx2) x1=nx2;
778 if (y0<-ny2) y0=-ny2;
779 if (y1>ny2) y1=ny2;
780 if (z0<-nz2) z0=-nz2;
781 if (z1>nz2) z1=nz2;
782
783 float w=weight;
784 float ws [ N/2 + 1 ];
785 float alpha = 32.0;
786 float wm = 0.0;
787
788 // compute 1D weights... not exactly correct, but somewhat close.
789 for ( int p = 0; p <= N/2; p++) {
790 double tmp = gsl_sf_bessel_i0_scaled(M_PI * alpha * sqrt(1.0f - Util::square((((N/2)+p)/(N-1))-1))) / gsl_sf_bessel_i0_scaled(M_PI * alpha);
791 ws[p] = (float) tmp;
792 wm += (float) tmp;
793 }
794
795 float r, kb, dn;
796 for (int k = z0 ; k <= z1; k++) {
797 for (int j = y0 ; j <= y1; j++) {
798 for (int i = x0; i <= x1; i ++) {
799 r = Util::hypot3sq((float) i - xx, j - yy, k - zz);
800 kb = 0.0;
801 //quasi radial...true cumulative radial weights are much more time consuming to code.
802 for (int p = 0; p <= Util::get_min(r,(float) N/2); p++) {
803 kb += ws[p];
804 }
805 dn = sqrt(kb/wm);
806 size_t off = data->add_complex_at_fast(i,j,k,dt*dn*w);
807 norm[off/2]+=w;
808 }
809 }
810 }
811 return true;
812 }
813 printf("region writing not supported in mode 10\n");
814 return false;
815}
816
817const float FourierInserter3DMode11::kernel[12][12][12] = {
8180.3756136,0.3403310,0.2474435,0.1296879,0.0245121,-0.0418448,-0.0631182,-0.0511850,-0.0262455,-0.0054893,0.0036277,0.0000000,
8190.3403310,0.3082430,0.2237843,0.1167696,0.0212886,-0.0388006,-0.0578490,-0.0467286,-0.0238678,-0.0049178,0.0033657,0.0000000,
8200.2474435,0.2237843,0.1615623,0.0828749,0.0129467,-0.0306454,-0.0438698,-0.0349361,-0.0175894,-0.0034180,0.0026655,0.0000000,
8210.1296879,0.1167696,0.0828749,0.0402484,0.0028031,-0.0198868,-0.0258289,-0.0198082,-0.0095752,-0.0015321,0.0017471,0.0000000,
8220.0245121,0.0212886,0.0129467,0.0028031,-0.0054675,-0.0095100,-0.0091333,-0.0059726,-0.0023192,0.0001220,0.0008690,0.0000000,
823-0.0418448,-0.0388006,-0.0306454,-0.0198868,-0.0095100,-0.0018223,0.0022627,0.0032331,0.0023989,0.0011163,0.0002267,0.0000000,
824-0.0631182,-0.0578490,-0.0438698,-0.0258289,-0.0091333,0.0022627,0.0071368,0.0068533,0.0041017,0.0013575,-0.0001089,0.0000000,
825-0.0511850,-0.0467286,-0.0349361,-0.0198082,-0.0059726,0.0032331,0.0068533,0.0061573,0.0035159,0.0010550,-0.0001880,0.0000000,
826-0.0262455,-0.0238678,-0.0175894,-0.0095752,-0.0023192,0.0023989,0.0041017,0.0035159,0.0019496,0.0005548,-0.0001273,0.0000000,
827-0.0054893,-0.0049178,-0.0034180,-0.0015321,0.0001220,0.0011163,0.0013575,0.0010550,0.0005548,0.0001483,-0.0000403,0.0000000,
8280.0036277,0.0033657,0.0026655,0.0017471,0.0008690,0.0002267,-0.0001089,-0.0001880,-0.0001273,-0.0000403,0.0000092,0.0000000,
8290.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
8300.3403310,0.3082430,0.2237843,0.1167696,0.0212886,-0.0388006,-0.0578490,-0.0467286,-0.0238678,-0.0049178,0.0033657,0.0000000,
8310.3082430,0.2790696,0.2022997,0.1050786,0.0184299,-0.0359603,-0.0530014,-0.0426434,-0.0216941,-0.0043990,0.0031230,0.0000000,
8320.2237843,0.2022996,0.1458102,0.0744136,0.0110401,-0.0283535,-0.0401427,-0.0318353,-0.0159554,-0.0030383,0.0024746,0.0000000,
8330.1167696,0.1050786,0.0744136,0.0358774,0.0020773,-0.0183257,-0.0235546,-0.0179765,-0.0086344,-0.0013296,0.0016238,0.0000000,
8340.0212885,0.0184299,0.0110400,0.0020773,-0.0051866,-0.0086670,-0.0082158,-0.0053132,-0.0020141,0.0001649,0.0008100,0.0000000,
835-0.0388006,-0.0359603,-0.0283535,-0.0183257,-0.0086670,-0.0015308,0.0022350,0.0030945,0.0022783,0.0010565,0.0002142,0.0000000,
836-0.0578490,-0.0530014,-0.0401427,-0.0235546,-0.0082158,0.0022350,0.0066778,0.0063753,0.0038095,0.0012622,-0.0000978,0.0000000,
837-0.0467286,-0.0426434,-0.0318353,-0.0179765,-0.0053133,0.0030945,0.0063753,0.0057001,0.0032491,0.0009749,-0.0001724,0.0000000,
838-0.0238678,-0.0216940,-0.0159554,-0.0086344,-0.0020141,0.0022783,0.0038095,0.0032491,0.0017978,0.0005110,-0.0001172,0.0000000,
839-0.0049178,-0.0043990,-0.0030383,-0.0013296,0.0001649,0.0010565,0.0012622,0.0009749,0.0005110,0.0001362,-0.0000371,0.0000000,
8400.0033657,0.0031230,0.0024746,0.0016238,0.0008100,0.0002142,-0.0000978,-0.0001724,-0.0001172,-0.0000371,0.0000085,0.0000000,
8410.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
8420.2474435,0.2237843,0.1615623,0.0828749,0.0129467,-0.0306454,-0.0438698,-0.0349361,-0.0175894,-0.0034180,0.0026655,0.0000000,
8430.2237843,0.2022997,0.1458102,0.0744136,0.0110401,-0.0283535,-0.0401427,-0.0318353,-0.0159554,-0.0030383,0.0024746,0.0000000,
8440.1615622,0.1458101,0.1044314,0.0522464,0.0061329,-0.0222221,-0.0302624,-0.0236374,-0.0116456,-0.0020445,0.0019643,0.0000000,
8450.0828748,0.0744135,0.0522463,0.0244677,0.0002465,-0.0141589,-0.0175347,-0.0131427,-0.0061593,-0.0008027,0.0012941,0.0000000,
8460.0129467,0.0110400,0.0061328,0.0002465,-0.0044008,-0.0064289,-0.0058001,-0.0035858,-0.0012204,0.0002714,0.0006519,0.0000000,
847-0.0306454,-0.0283535,-0.0222221,-0.0141589,-0.0064289,-0.0007715,0.0021423,0.0027101,0.0019478,0.0008934,0.0001802,0.0000000,
848-0.0438697,-0.0401427,-0.0302624,-0.0175348,-0.0058001,0.0021423,0.0054435,0.0050956,0.0030280,0.0010071,-0.0000689,0.0000000,
849-0.0349361,-0.0318353,-0.0236374,-0.0131427,-0.0035858,0.0027101,0.0050956,0.0044802,0.0025380,0.0007613,-0.0001310,0.0000000,
850-0.0175894,-0.0159554,-0.0116456,-0.0061593,-0.0012204,0.0019478,0.0030280,0.0025380,0.0013935,0.0003943,-0.0000903,0.0000000,
851-0.0034180,-0.0030383,-0.0020445,-0.0008027,0.0002714,0.0008934,0.0010071,0.0007613,0.0003943,0.0001039,-0.0000288,0.0000000,
8520.0026655,0.0024746,0.0019643,0.0012941,0.0006519,0.0001802,-0.0000689,-0.0001310,-0.0000903,-0.0000288,0.0000065,0.0000000,
8530.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
8540.1296879,0.1167696,0.0828749,0.0402484,0.0028031,-0.0198868,-0.0258289,-0.0198082,-0.0095752,-0.0015321,0.0017471,0.0000000,
8550.1167696,0.1050786,0.0744136,0.0358774,0.0020773,-0.0183257,-0.0235546,-0.0179765,-0.0086344,-0.0013296,0.0016238,0.0000000,
8560.0828748,0.0744135,0.0522463,0.0244677,0.0002465,-0.0141589,-0.0175347,-0.0131427,-0.0061593,-0.0008027,0.0012941,0.0000000,
8570.0402484,0.0358773,0.0244676,0.0102938,-0.0018371,-0.0087088,-0.0098086,-0.0069817,-0.0030272,-0.0001546,0.0008601,0.0000000,
8580.0028031,0.0020773,0.0002465,-0.0018371,-0.0032679,-0.0035382,-0.0027388,-0.0014218,-0.0002427,0.0003867,0.0004426,0.0000000,
859-0.0198868,-0.0183257,-0.0141589,-0.0087088,-0.0035382,0.0001652,0.0019638,0.0021631,0.0014894,0.0006692,0.0001335,0.0000000,
860-0.0258289,-0.0235545,-0.0175347,-0.0098086,-0.0027388,0.0019638,0.0038001,0.0034089,0.0020002,0.0006705,-0.0000328,0.0000000,
861-0.0198081,-0.0179764,-0.0131427,-0.0069817,-0.0014218,0.0021631,0.0034089,0.0028848,0.0016103,0.0004821,-0.0000781,0.0000000,
862-0.0095751,-0.0086344,-0.0061593,-0.0030272,-0.0002427,0.0014894,0.0020002,0.0016103,0.0008674,0.0002423,-0.0000557,0.0000000,
863-0.0015321,-0.0013296,-0.0008027,-0.0001546,0.0003867,0.0006692,0.0006705,0.0004821,0.0002423,0.0000617,-0.0000183,0.0000000,
8640.0017471,0.0016238,0.0012941,0.0008601,0.0004426,0.0001335,-0.0000328,-0.0000781,-0.0000557,-0.0000183,0.0000038,0.0000000,
8650.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
8660.0245121,0.0212886,0.0129467,0.0028031,-0.0054675,-0.0095100,-0.0091333,-0.0059726,-0.0023192,0.0001220,0.0008690,0.0000000,
8670.0212885,0.0184299,0.0110400,0.0020773,-0.0051866,-0.0086670,-0.0082158,-0.0053132,-0.0020141,0.0001649,0.0008100,0.0000000,
8680.0129467,0.0110400,0.0061328,0.0002465,-0.0044008,-0.0064289,-0.0058001,-0.0035858,-0.0012204,0.0002714,0.0006519,0.0000000,
8690.0028031,0.0020773,0.0002465,-0.0018371,-0.0032679,-0.0035382,-0.0027388,-0.0014218,-0.0002427,0.0003867,0.0004426,0.0000000,
870-0.0054675,-0.0051866,-0.0044008,-0.0032679,-0.0020109,-0.0008640,-0.0000115,0.0004597,0.0005759,0.0004525,0.0002390,0.0000000,
871-0.0095100,-0.0086670,-0.0064289,-0.0035382,-0.0008640,0.0009487,0.0016881,0.0015623,0.0010061,0.0004364,0.0000852,0.0000000,
872-0.0091333,-0.0082158,-0.0058001,-0.0027388,-0.0000115,0.0016881,0.0021848,0.0017821,0.0010130,0.0003454,-0.0000014,0.0000000,
873-0.0059726,-0.0053132,-0.0035858,-0.0014218,0.0004596,0.0015623,0.0017821,0.0013691,0.0007327,0.0002175,-0.0000298,0.0000000,
874-0.0023192,-0.0020141,-0.0012204,-0.0002427,0.0005759,0.0010061,0.0010130,0.0007327,0.0003725,0.0000992,-0.0000240,0.0000000,
8750.0001220,0.0001649,0.0002714,0.0003867,0.0004525,0.0004364,0.0003454,0.0002175,0.0000992,0.0000219,-0.0000088,0.0000000,
8760.0008690,0.0008101,0.0006519,0.0004426,0.0002390,0.0000852,-0.0000014,-0.0000298,-0.0000240,-0.0000088,0.0000011,0.0000000,
8770.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
878-0.0418448,-0.0388006,-0.0306454,-0.0198868,-0.0095100,-0.0018223,0.0022627,0.0032331,0.0023989,0.0011163,0.0002267,0.0000000,
879-0.0388006,-0.0359603,-0.0283535,-0.0183257,-0.0086670,-0.0015308,0.0022350,0.0030945,0.0022783,0.0010565,0.0002142,0.0000000,
880-0.0306454,-0.0283535,-0.0222221,-0.0141589,-0.0064289,-0.0007715,0.0021423,0.0027101,0.0019478,0.0008934,0.0001802,0.0000000,
881-0.0198868,-0.0183257,-0.0141589,-0.0087088,-0.0035382,0.0001652,0.0019638,0.0021631,0.0014894,0.0006692,0.0001335,0.0000000,
882-0.0095100,-0.0086670,-0.0064289,-0.0035382,-0.0008640,0.0009487,0.0016881,0.0015623,0.0010061,0.0004364,0.0000852,0.0000000,
883-0.0018223,-0.0015308,-0.0007715,0.0001652,0.0009487,0.0013519,0.0013308,0.0010089,0.0005874,0.0002396,0.0000447,0.0000000,
8840.0022627,0.0022350,0.0021423,0.0019638,0.0016881,0.0013308,0.0009379,0.0005713,0.0002847,0.0001029,0.0000171,0.0000000,
8850.0032331,0.0030945,0.0027101,0.0021631,0.0015623,0.0010089,0.0005713,0.0002745,0.0001051,0.0000275,0.0000025,0.0000000,
8860.0023990,0.0022783,0.0019478,0.0014894,0.0010061,0.0005874,0.0002847,0.0001051,0.0000226,-0.0000019,-0.0000026,0.0000000,
8870.0011163,0.0010565,0.0008934,0.0006692,0.0004364,0.0002396,0.0001029,0.0000275,-0.0000019,-0.0000062,-0.0000025,0.0000000,
8880.0002267,0.0002142,0.0001802,0.0001335,0.0000852,0.0000447,0.0000171,0.0000025,-0.0000026,-0.0000025,-0.0000010,0.0000000,
8890.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
890-0.0631182,-0.0578490,-0.0438698,-0.0258289,-0.0091333,0.0022627,0.0071368,0.0068533,0.0041017,0.0013575,-0.0001089,0.0000000,
891-0.0578490,-0.0530014,-0.0401427,-0.0235546,-0.0082158,0.0022350,0.0066778,0.0063753,0.0038095,0.0012622,-0.0000978,0.0000000,
892-0.0438697,-0.0401427,-0.0302624,-0.0175348,-0.0058001,0.0021423,0.0054435,0.0050956,0.0030280,0.0010071,-0.0000689,0.0000000,
893-0.0258289,-0.0235545,-0.0175347,-0.0098086,-0.0027388,0.0019638,0.0038001,0.0034089,0.0020002,0.0006705,-0.0000328,0.0000000,
894-0.0091333,-0.0082158,-0.0058001,-0.0027388,-0.0000115,0.0016881,0.0021848,0.0017821,0.0010130,0.0003454,-0.0000014,0.0000000,
8950.0022627,0.0022350,0.0021423,0.0019638,0.0016881,0.0013308,0.0009379,0.0005713,0.0002847,0.0001029,0.0000171,0.0000000,
8960.0071368,0.0066778,0.0054435,0.0038001,0.0021848,0.0009379,0.0002000,-0.0000882,-0.0001034,-0.0000291,0.0000213,0.0000000,
8970.0068532,0.0063753,0.0050956,0.0034089,0.0017821,0.0005713,-0.0000882,-0.0002787,-0.0002039,-0.0000665,0.0000155,0.0000000,
8980.0041017,0.0038095,0.0030280,0.0020002,0.0010130,0.0002847,-0.0001034,-0.0002039,-0.0001440,-0.0000496,0.0000066,0.0000000,
8990.0013574,0.0012622,0.0010071,0.0006705,0.0003454,0.0001029,-0.0000291,-0.0000665,-0.0000496,-0.0000193,0.0000001,0.0000000,
900-0.0001089,-0.0000978,-0.0000689,-0.0000328,-0.0000014,0.0000171,0.0000213,0.0000155,0.0000066,0.0000001,-0.0000022,0.0000000,
9010.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
902-0.0511850,-0.0467286,-0.0349361,-0.0198082,-0.0059726,0.0032331,0.0068533,0.0061573,0.0035159,0.0010550,-0.0001880,0.0000000,
903-0.0467286,-0.0426434,-0.0318353,-0.0179765,-0.0053133,0.0030945,0.0063753,0.0057001,0.0032491,0.0009749,-0.0001724,0.0000000,
904-0.0349361,-0.0318353,-0.0236374,-0.0131427,-0.0035858,0.0027101,0.0050956,0.0044802,0.0025380,0.0007613,-0.0001310,0.0000000,
905-0.0198081,-0.0179764,-0.0131427,-0.0069817,-0.0014218,0.0021631,0.0034089,0.0028848,0.0016103,0.0004821,-0.0000781,0.0000000,
906-0.0059726,-0.0053132,-0.0035858,-0.0014218,0.0004596,0.0015623,0.0017821,0.0013691,0.0007327,0.0002175,-0.0000298,0.0000000,
9070.0032331,0.0030945,0.0027101,0.0021631,0.0015623,0.0010089,0.0005713,0.0002745,0.0001051,0.0000275,0.0000025,0.0000000,
9080.0068532,0.0063753,0.0050956,0.0034089,0.0017821,0.0005713,-0.0000882,-0.0002787,-0.0002039,-0.0000665,0.0000155,0.0000000,
9090.0061573,0.0057001,0.0044802,0.0028848,0.0013691,0.0002745,-0.0002787,-0.0003841,-0.0002511,-0.0000812,0.0000141,0.0000000,
9100.0035159,0.0032491,0.0025380,0.0016103,0.0007327,0.0001051,-0.0002039,-0.0002511,-0.0001613,-0.0000537,0.0000066,0.0000000,
9110.0010550,0.0009749,0.0007613,0.0004821,0.0002175,0.0000275,-0.0000665,-0.0000812,-0.0000537,-0.0000199,0.0000000,0.0000000,
912-0.0001880,-0.0001724,-0.0001310,-0.0000781,-0.0000298,0.0000025,0.0000155,0.0000141,0.0000066,0.0000000,-0.0000024,0.0000000,
9130.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
914-0.0262455,-0.0238678,-0.0175894,-0.0095752,-0.0023192,0.0023989,0.0041017,0.0035159,0.0019496,0.0005548,-0.0001273,0.0000000,
915-0.0238678,-0.0216940,-0.0159554,-0.0086344,-0.0020141,0.0022783,0.0038095,0.0032491,0.0017978,0.0005110,-0.0001172,0.0000000,
916-0.0175894,-0.0159554,-0.0116456,-0.0061593,-0.0012204,0.0019478,0.0030280,0.0025380,0.0013935,0.0003943,-0.0000903,0.0000000,
917-0.0095751,-0.0086344,-0.0061593,-0.0030272,-0.0002427,0.0014894,0.0020002,0.0016103,0.0008674,0.0002423,-0.0000557,0.0000000,
918-0.0023192,-0.0020141,-0.0012204,-0.0002427,0.0005759,0.0010061,0.0010130,0.0007327,0.0003725,0.0000992,-0.0000240,0.0000000,
9190.0023990,0.0022783,0.0019478,0.0014894,0.0010061,0.0005874,0.0002847,0.0001051,0.0000226,-0.0000019,-0.0000026,0.0000000,
9200.0041017,0.0038095,0.0030280,0.0020002,0.0010130,0.0002847,-0.0001034,-0.0002039,-0.0001440,-0.0000496,0.0000066,0.0000000,
9210.0035159,0.0032491,0.0025380,0.0016102,0.0007327,0.0001051,-0.0002039,-0.0002511,-0.0001613,-0.0000537,0.0000066,0.0000000,
9220.0019496,0.0017978,0.0013935,0.0008674,0.0003725,0.0000226,-0.0001440,-0.0001613,-0.0001020,-0.0000349,0.0000026,0.0000000,
9230.0005548,0.0005110,0.0003943,0.0002423,0.0000992,-0.0000019,-0.0000496,-0.0000537,-0.0000349,-0.0000136,-0.0000009,0.0000000,
924-0.0001273,-0.0001172,-0.0000903,-0.0000557,-0.0000240,-0.0000026,0.0000066,0.0000066,0.0000026,-0.0000009,-0.0000020,0.0000000,
9250.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
926-0.0054893,-0.0049178,-0.0034180,-0.0015321,0.0001220,0.0011163,0.0013575,0.0010550,0.0005548,0.0001483,-0.0000403,0.0000000,
927-0.0049178,-0.0043990,-0.0030383,-0.0013296,0.0001649,0.0010565,0.0012622,0.0009749,0.0005110,0.0001362,-0.0000371,0.0000000,
928-0.0034180,-0.0030383,-0.0020445,-0.0008027,0.0002714,0.0008934,0.0010071,0.0007613,0.0003943,0.0001039,-0.0000288,0.0000000,
929-0.0015321,-0.0013296,-0.0008027,-0.0001546,0.0003867,0.0006692,0.0006705,0.0004821,0.0002423,0.0000617,-0.0000183,0.0000000,
9300.0001220,0.0001649,0.0002714,0.0003867,0.0004525,0.0004364,0.0003454,0.0002175,0.0000992,0.0000219,-0.0000088,0.0000000,
9310.0011163,0.0010565,0.0008934,0.0006692,0.0004364,0.0002396,0.0001029,0.0000275,-0.0000019,-0.0000062,-0.0000025,0.0000000,
9320.0013574,0.0012622,0.0010071,0.0006705,0.0003454,0.0001029,-0.0000291,-0.0000665,-0.0000496,-0.0000193,0.0000001,0.0000000,
9330.0010550,0.0009749,0.0007613,0.0004821,0.0002175,0.0000275,-0.0000665,-0.0000812,-0.0000537,-0.0000199,0.0000000,0.0000000,
9340.0005548,0.0005110,0.0003943,0.0002423,0.0000992,-0.0000019,-0.0000496,-0.0000537,-0.0000349,-0.0000136,-0.0000009,0.0000000,
9350.0001483,0.0001362,0.0001039,0.0000617,0.0000219,-0.0000062,-0.0000193,-0.0000199,-0.0000136,-0.0000063,-0.0000015,0.0000000,
936-0.0000403,-0.0000371,-0.0000288,-0.0000183,-0.0000088,-0.0000025,0.0000001,0.0000000,-0.0000009,-0.0000015,-0.0000012,0.0000000,
9370.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
9380.0036277,0.0033657,0.0026655,0.0017471,0.0008690,0.0002267,-0.0001089,-0.0001880,-0.0001273,-0.0000403,0.0000092,0.0000000,
9390.0033657,0.0031230,0.0024746,0.0016238,0.0008100,0.0002142,-0.0000978,-0.0001724,-0.0001172,-0.0000371,0.0000085,0.0000000,
9400.0026655,0.0024746,0.0019643,0.0012941,0.0006519,0.0001802,-0.0000689,-0.0001310,-0.0000903,-0.0000288,0.0000065,0.0000000,
9410.0017471,0.0016238,0.0012941,0.0008601,0.0004426,0.0001335,-0.0000328,-0.0000781,-0.0000557,-0.0000183,0.0000038,0.0000000,
9420.0008690,0.0008101,0.0006519,0.0004426,0.0002390,0.0000852,-0.0000014,-0.0000298,-0.0000240,-0.0000088,0.0000011,0.0000000,
9430.0002267,0.0002142,0.0001802,0.0001335,0.0000852,0.0000447,0.0000171,0.0000025,-0.0000026,-0.0000025,-0.0000010,0.0000000,
944-0.0001089,-0.0000978,-0.0000689,-0.0000328,-0.0000014,0.0000171,0.0000213,0.0000155,0.0000066,0.0000001,-0.0000022,0.0000000,
945-0.0001880,-0.0001724,-0.0001310,-0.0000781,-0.0000298,0.0000025,0.0000155,0.0000141,0.0000066,0.0000000,-0.0000024,0.0000000,
946-0.0001273,-0.0001172,-0.0000903,-0.0000557,-0.0000240,-0.0000026,0.0000066,0.0000066,0.0000026,-0.0000009,-0.0000020,0.0000000,
947-0.0000403,-0.0000371,-0.0000288,-0.0000183,-0.0000088,-0.0000025,0.0000001,0.0000000,-0.0000009,-0.0000015,-0.0000012,0.0000000,
9480.0000092,0.0000085,0.0000065,0.0000038,0.0000011,-0.0000010,-0.0000022,-0.0000024,-0.0000020,-0.0000012,-0.0000005,0.0000000,
9490.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
9500.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
9510.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
9520.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
9530.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
9540.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
9550.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
9560.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
9570.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
9580.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
9590.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
9600.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,
9610.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000,0.0000000
962};
963
964bool FourierInserter3DMode11::insert_pixel(const float& xx, const float& yy, const float& zz, const std::complex<float> dt,const float& weight)
965{
966 int x0 = Util::fast_floor(xx-2.5);
967 int y0 = Util::fast_floor(yy-2.5);
968 int z0 = Util::fast_floor(zz-2.5);
969
970 if (subx0<0) { // normal full reconstruction
971 if (x0<-nx2-6 || y0<-ny2-6 || z0<-nz2-6 || x0>nx2+5 || y0>ny2+5 || z0>nz2+5 ) return false;
972
973 // no error checking on add_complex_fast, so we need to be careful here
974 int x1=x0+6;
975 int y1=y0+6;
976 int z1=z0+6;
977 if (x0<-nx2) x0=-nx2;
978 if (x1>nx2) x1=nx2;
979 if (y0<-ny2) y0=-ny2;
980 if (y1>ny2) y1=ny2;
981 if (z0<-nz2) z0=-nz2;
982 if (z1>nz2) z1=nz2;
983
984 float w=weight;
985
986 float r, gg;
987 for (int k = z0 ; k <= z1; k++) {
988 for (int j = y0 ; j <= y1; j++) {
989 for (int i = x0; i <= x1; i ++) {
990 gg=FourierInserter3DMode11::kernel[abs(Util::fast_floor((i-xx)*3.0f+0.5))][abs(Util::fast_floor((j-yy)*3.0f+0.5))][abs(Util::fast_floor((k-zz)*3.0f+0.5))];
991
992 size_t off;
993 off=data->add_complex_at_fast(i,j,k,dt*gg*w);
994 norm[off/2]+=w;
995// if (i==67&&j==19&&k==1) printf("%1.1f %1.1f %1.1f\t%d %d %d\t%d %d %d\t%f\t%f\t%f\t%f\t%f\t%f\n",xx,yy,zz,i,j,k,abs(Util::fast_floor((i-xx)*3.0f+0.5)),abs(Util::fast_floor((j-yy)*3.0f+0.5)),abs(Util::fast_floor((k-zz)*3.0f+0.5)),gg,w,dt.real(),dt.imag(),norm[off/2],data->get_value_at(off));
996
997 }
998 }
999 }
1000 return true;
1001 }
1002 printf("region writing not supported in mode \n");
1003 return false;
1004}
EMAN::Dict::has_key
bool has_key(const string &key) const
Ask the Dictionary if it as a particular key.
Definition: emobject.h:511
EMAN::EMConsts::I5G
static const float I5G
Definition: emobject.h:75
EMAN::EMConsts::I3G
static const float I3G
Definition: emobject.h:73
EMAN::EMConsts::I2G
static const float I2G
Definition: emobject.h:72
EMAN::FactoryBase::params
Dict params
This is the dictionary the stores the parameters of the object.
Definition: emobject.h:953
EMAN::Factory
Factory is used to store objects to create new instances.
Definition: emobject.h:725
EMAN::FourierInserter3DMode10::insert_pixel
virtual bool insert_pixel(const float &xx, const float &yy, const float &zz, const std::complex< float > dt, const float &weight=1.0)
Insert a complex pixel [dt[0]+dt[1]i] at (float) coordinate [xx,yy,zz] with weighting into a discrete...
Definition: reconstructor_tools.cpp:760
EMAN::FourierInserter3DMode11::kernel
static const float kernel[12][12][12]
Definition: reconstructor_tools.h:536
EMAN::FourierInserter3DMode11::insert_pixel
virtual bool insert_pixel(const float &xx, const float &yy, const float &zz, const std::complex< float > dt, const float &weight=1.0)
Insert a complex pixel [dt[0]+dt[1]i] at (float) coordinate [xx,yy,zz] with weighting into a discrete...
Definition: reconstructor_tools.cpp:964
EMAN::FourierInserter3DMode1::insert_pixel
virtual bool insert_pixel(const float &xx, const float &yy, const float &zz, const std::complex< float > dt, const float &weight=1.0)
Insert a complex pixel [dt[0]+dt[1]i] at (float) coordinate [xx,yy,zz] with weighting into a discrete...
Definition: reconstructor_tools.cpp:113
EMAN::FourierInserter3DMode2::insert_pixel
virtual bool insert_pixel(const float &xx, const float &yy, const float &zz, const std::complex< float > dt, const float &weight=1.0)
Insert a complex pixel [dt[0]+dt[1]i] at (float) coordinate [xx,yy,zz] with weighting into a discrete...
Definition: reconstructor_tools.cpp:128
EMAN::FourierInserter3DMode2l::insert_pixel
virtual bool insert_pixel(const float &xx, const float &yy, const float &zz, const std::complex< float > dt, const float &weight=1.0)
Insert a complex pixel [dt[0]+dt[1]i] at (float) coordinate [xx,yy,zz] with weighting into a discrete...
Definition: reconstructor_tools.cpp:191
EMAN::FourierInserter3DMode3::insert_pixel
virtual bool insert_pixel(const float &xx, const float &yy, const float &zz, const std::complex< float > dt, const float &weight=1.0)
Insert a complex pixel [dt[0]+dt[1]i] at (float) coordinate [xx,yy,zz] with weighting into a discrete...
Definition: reconstructor_tools.cpp:265
EMAN::FourierInserter3DMode5::insert_pixel
virtual bool insert_pixel(const float &xx, const float &yy, const float &zz, const std::complex< float > dt, const float &weight=1.0)
Insert a complex pixel [dt[0]+dt[1]i] at (float) coordinate [xx,yy,zz] with weighting into a discrete...
Definition: reconstructor_tools.cpp:315
EMAN::FourierInserter3DMode6::insert_pixel
virtual bool insert_pixel(const float &xx, const float &yy, const float &zz, const std::complex< float > dt, const float &weight=1.0)
Insert a complex pixel [dt[0]+dt[1]i] at (float) coordinate [xx,yy,zz] with weighting into a discrete...
Definition: reconstructor_tools.cpp:382
EMAN::FourierInserter3DMode7::insert_pixel
virtual bool insert_pixel(const float &xx, const float &yy, const float &zz, const std::complex< float > dt, const float &weight=1.0)
Insert a complex pixel [dt[0]+dt[1]i] at (float) coordinate [xx,yy,zz] with weighting into a discrete...
Definition: reconstructor_tools.cpp:539
EMAN::FourierInserter3DMode7::kernel
static const float kernel[9][9][9]
Definition: reconstructor_tools.h:380
EMAN::FourierInserter3DMode8::insert_pixel
virtual bool insert_pixel(const float &xx, const float &yy, const float &zz, const std::complex< float > dt, const float &weight=1.0)
Insert a complex pixel [dt[0]+dt[1]i] at (float) coordinate [xx,yy,zz] with weighting into a discrete...
Definition: reconstructor_tools.cpp:594
EMAN::FourierInserter3DMode9::insert_pixel
virtual bool insert_pixel(const float &xx, const float &yy, const float &zz, const std::complex< float > dt, const float &weight=1.0)
Insert a complex pixel [dt[0]+dt[1]i] at (float) coordinate [xx,yy,zz] with weighting into a discrete...
Definition: reconstructor_tools.cpp:714
EMAN::FourierPixelInserter3D::data
EMData * data
A pointer to the constructor argument real_data.
Definition: reconstructor_tools.h:139
EMAN::FourierPixelInserter3D::nx
int nx
Image volume data sizes a convenience variable used here and there.
Definition: reconstructor_tools.h:142
EMAN::FourierPixelInserter3D::norm
float * norm
A pointer to the constructor argument normalize_values.
Definition: reconstructor_tools.h:137
EMAN::Util::hypot3sq
static int hypot3sq(int x, int y, int z)
Euclidean distance function squared in 3D: f(x,y,z) = (x*x + y*y + z*z);.
Definition: util.h:805
EMAN::Util::square
static int square(int n)
Calculate a number's square.
Definition: util.h:736
EMAN::Util::fast_floor
static int fast_floor(float x)
A fast way to calculate a floor, which is largest integral value not greater than argument.
Definition: util.h:874
EMAN::Util::hypot3
static float hypot3(int x, int y, int z)
Euclidean distance function in 3D: f(x,y,z) = sqrt(x*x + y*y + z*z);.
Definition: util.h:827
EMAN::Util::get_min
static int get_min(int f1, int f2)
Get the minimum of 2 numbers.
Definition: util.h:922
EMAN::Util::fast_exp
static float fast_exp(const float &f)
Returns an approximate of exp(x) using a cached table uses actual exp(x) outside the cached range.
Definition: util.cpp:788
sqrt
EMData * sqrt() const
return square root of current image
Definition: emdata_core.cpp:1084
NotExistingObjectException
#define NotExistingObjectException(objname, desc)
Definition: exception.h:130
EMAN
E2Exception class.
Definition: aligner.h:40
reconstructor_tools.h
Generated by doxygen 1.9.2