#include <aligner.h>
Inheritance diagram for EMAN::Refine3DAlignerGrid:
Public Member Functions | |
virtual EMData * | align (EMData *this_img, EMData *to_img, const string &cmp_name="sqeuclidean", const Dict &cmp_params=Dict()) const |
To align 'this_img' with another image passed in through its parameters. | |
virtual EMData * | align (EMData *this_img, EMData *to_img) const |
virtual string | get_name () const |
Get the Aligner's name. | |
virtual string | get_desc () const |
virtual TypeDict | get_param_types () const |
Static Public Member Functions | |
Aligner * | NEW () |
Static Public Attributes | |
const string | NAME = "refine_3d_grid" |
Refines a preliminary 3D alignment using a sampling grid. This is a port from tomohunter, but the az sampling scheme is altered cuch that the points on the sphere are equidistant (Improves speed several hundered times). The distance between the points on the sphere is 'delta' and the range(distance from the pole, 0,0,0 position) is given as 'range'. IN general this refinement scheme is a bit slower than the Quaternion Simplex aligner, but perfroms better in the presence of noise(according to current dogma).
xform.align3d | The Transform from the previous course alignment. If unspecified the identity matrix is used | |
delta | The angluar distance bewteen points on the spehere used in the grid | |
range | The size of the grid. Measured in angluar distance from the north pole | |
dotrans | Do a translational search | |
search | The maximum length of the detectable translational shift - if you supply this parameter you can not supply the maxshiftx, maxshifty or maxshiftz parameters. Each approach is mutually exclusive | |
searchx | The maximum length of the detectable translational shift in the x direction- if you supply this parameter you can not supply the maxshift parameters | |
searchy | The maximum length of the detectable translational shift in the y direction- if you supply this parameter you can not supply the maxshift parameters | |
searchz | The maximum length of the detectable translational shift in the z direction- if you supply this parameter you can not supply the maxshift parameters | |
verbose | Turn this on to have useful information printed to standard out |
Definition at line 1384 of file aligner.h.
|
Implements EMAN::Aligner. Definition at line 1390 of file aligner.h. References align(). 01391 { 01392 return align(this_img, to_img, "sqeuclidean", Dict()); 01393 }
|
|
To align 'this_img' with another image passed in through its parameters. The alignment uses a user-given comparison method to compare the two images. If none is given, a default one is used.
Implements EMAN::Aligner. Definition at line 2263 of file aligner.cpp. References EMAN::EMData::calc_ccf(), EMAN::EMData::calc_max_location_wrap(), calc_max_location_wrap_cuda(), EMAN::Cmp::cmp(), data, EMAN::EMData::do_fft(), EMAN::Dict::end(), EMAN::EMData::get_ndim(), get_stats_cuda(), EMAN::EMData::get_value_at_wrap(), EMAN::EMData::get_xsize(), EMAN::EMData::get_ysize(), EMAN::EMData::get_zsize(), EMAN::Dict::has_key(), ImageDimensionException, InvalidParameterException, CudaPeakInfo::peak, phi, EMAN::EMData::process(), EMAN::EMData::process_inplace(), CudaPeakInfo::px, CudaPeakInfo::py, CudaPeakInfo::pz, EMAN::EMData::set_attr(), EMAN::Dict::set_default(), EMAN::Transform::set_trans(), sqrt(), and t. 02265 { 02266 if ( this_img->get_ndim() != 3 || to->get_ndim() != 3 ) { 02267 throw ImageDimensionException("This aligner only works for 3D images"); 02268 } 02269 02270 Transform* t; 02271 if (params.has_key("xform.align3d") ) { 02272 // Unlike the 2d refine aligner, this class doesn't require the starting transform's 02273 // parameters to form the starting guess. Instead the Transform itself 02274 // is perturbed carefully (using quaternion rotation) to overcome problems that arise 02275 // when you use orthogonally-based Euler angles 02276 t = params["xform.align3d"]; 02277 }else { 02278 t = new Transform(); // is the identity 02279 } 02280 02281 int searchx = 0; 02282 int searchy = 0; 02283 int searchz = 0; 02284 bool dotrans = params.set_default("dotrans",1); 02285 if (params.has_key("search")) { 02286 vector<string> check; 02287 check.push_back("searchx"); 02288 check.push_back("searchy"); 02289 check.push_back("searchz"); 02290 for(vector<string>::const_iterator cit = check.begin(); cit != check.end(); ++cit) { 02291 if (params.has_key(*cit)) throw InvalidParameterException("The search parameter is mutually exclusive of the searchx, searchy, and searchz parameters"); 02292 } 02293 int search = params["search"]; 02294 searchx = search; 02295 searchy = search; 02296 searchz = search; 02297 } else { 02298 searchx = params.set_default("searchx",3); 02299 searchy = params.set_default("searchy",3); 02300 searchz = params.set_default("searchz",3); 02301 } 02302 02303 float delta = params.set_default("delta",1.0f); 02304 float range = params.set_default("range",10.0f); 02305 bool verbose = params.set_default("verbose",false); 02306 02307 bool tomography = (cmp_name == "ccc.tomo") ? 1 : 0; 02308 EMData * tofft = 0; 02309 if(dotrans || tomography){ 02310 tofft = to->do_fft(); 02311 } 02312 02313 #ifdef EMAN2_USING_CUDA 02314 if(EMData::usecuda == 1) { 02315 if(!this_img->getcudarodata()) this_img->copy_to_cudaro(); // This is safer 02316 if(!to->getcudarwdata()) to->copy_to_cuda(); 02317 if(to->getcudarwdata()){if(tofft) tofft->copy_to_cuda();} 02318 } 02319 #endif 02320 02321 Dict d; 02322 d["type"] = "eman"; // d is used in the loop below 02323 Dict best; 02324 // best["score"] = numeric_limits<float>::infinity(); 02325 best["score"] = 1.0e37; 02326 bool use_cpu = true; 02327 Transform tran = Transform(); 02328 Cmp* c = Factory <Cmp>::get(cmp_name, cmp_params); 02329 02330 for ( float alt = 0; alt < range; alt += delta) { 02331 // now compute a sane az step size 02332 float saz = 360; 02333 if(alt != 0) saz = delta/sin(alt*M_PI/180.0f); // This gives consistent az step sizes(arc lengths) 02334 for ( float az = 0; az < 360; az += saz ){ 02335 if (verbose) { 02336 cout << "Trying angle alt " << alt << " az " << az << endl; 02337 } 02338 // account for any changes in az 02339 for( float phi = -range-az; phi < range-az; phi += delta ) { 02340 d["alt"] = alt; 02341 d["phi"] = phi; 02342 d["az"] = az; 02343 Transform tr(d); 02344 tr = tr*(*t); // compose transforms, this moves to the pole (aprox) 02345 02346 //EMData* transformed = this_img->process("xform",Dict("transform",&tr)); 02347 EMData* transformed; 02348 transformed = this_img->process("xform",Dict("transform",&tr)); 02349 02350 //need to do things a bit diffrent if we want to compare two tomos 02351 float score = 0.0f; 02352 if(dotrans || tomography){ 02353 EMData* ccf = transformed->calc_ccf(tofft); 02354 #ifdef EMAN2_USING_CUDA 02355 if(EMData::usecuda == 1){ 02356 use_cpu = false; 02357 CudaPeakInfo* data = calc_max_location_wrap_cuda(ccf->getcudarwdata(), ccf->get_xsize(), ccf->get_ysize(), ccf->get_zsize(), searchx, searchy, searchz); 02358 tran.set_trans((float)-data->px, (float)-data->py, (float)-data->pz); 02359 //CudaPeakInfoFloat* data = calc_max_location_wrap_intp_cuda(ccf->getcudarwdata(), ccf->get_xsize(), ccf->get_ysize(), ccf->get_zsize(), searchx, searchy, searchz); 02360 //tran.set_trans(-data->xintp, -data->yintp, -data->zintp); 02361 tr = tran*tr; // to reflect the fact that we have done a rotation first and THEN a transformation 02362 if (tomography) { 02363 float2 stats = get_stats_cuda(ccf->getcudarwdata(), ccf->get_xsize(), ccf->get_ysize(), ccf->get_zsize()); 02364 score = -(data->peak - stats.x)/sqrt(stats.y); // Normalize, this is better than calling the norm processor since we only need to normalize one point 02365 } else { 02366 score = -data->peak; 02367 } 02368 delete data; 02369 } 02370 #endif 02371 if(use_cpu){ 02372 if(tomography) ccf->process_inplace("normalize"); 02373 //vector<float> fpoint = ccf->calc_max_location_wrap_intp(searchx,searchy,searchz); 02374 //tran.set_trans(-fpoint[0], -fpoint[1], -fpoint[2]); 02375 //score = -fpoint[3]; 02376 IntPoint point = ccf->calc_max_location_wrap(searchx,searchy,searchz); 02377 tran.set_trans((float)-point[0], (float)-point[1], (float)-point[2]); 02378 score = -ccf->get_value_at_wrap(point[0], point[1], point[2]); 02379 tr = tran*tr;// to reflect the fact that we have done a rotation first and THEN a transformation 02380 02381 } 02382 delete ccf; ccf =0; 02383 delete transformed; transformed = 0;// this is to stop a mem leak 02384 } 02385 02386 if(!tomography){ 02387 if(!transformed) transformed = this_img->process("xform",Dict("transform",&tr)); // we are returning a map 02388 score = c->cmp(to,transformed); 02389 delete transformed; transformed = 0;// this is to stop a mem leak 02390 } 02391 02392 if(score < float(best["score"])) { 02393 // printf("%f\n",score); 02394 best["score"] = score; 02395 best["xform.align3d"] = &tr; // I wonder if this will cause a mem leak? 02396 } 02397 } 02398 } 02399 } 02400 02401 if(tofft) {delete tofft; tofft = 0;} 02402 if (c != 0) delete c; 02403 02404 //make aligned map; 02405 EMData* best_match = this_img->process("xform",Dict("transform", best["xform.align3d"])); // we are returning a map 02406 best_match->set_attr("xform.align3d", best["xform.align3d"]); 02407 best_match->set_attr("score", float(best["score"])); 02408 02409 return best_match; 02410 02411 }
|
|
Implements EMAN::Aligner. Definition at line 1400 of file aligner.h. 01401 { 01402 return "Refines a preliminary 3D alignment using a simplex algorithm. Subpixel precision."; 01403 }
|
|
Get the Aligner's name. Each Aligner is identified by a unique name.
Implements EMAN::Aligner. Definition at line 1395 of file aligner.h. 01396 {
01397 return NAME;
01398 }
|
|
Implements EMAN::Aligner. Definition at line 1410 of file aligner.h. References EMAN::TypeDict::put(). 01411 { 01412 TypeDict d; 01413 d.put("xform.align3d", EMObject::TRANSFORM,"The Transform storing the starting guess. If unspecified the identity matrix is used"); 01414 d.put("delta", EMObject::FLOAT, "The angular step size. Default is 1." ); 01415 d.put("range", EMObject::FLOAT, "The angular range size. Default is 10." ); 01416 d.put("dotrans", EMObject::BOOL,"Do a translational search. Default is True(1)"); 01417 d.put("search", EMObject::INT,"The maximum length of the detectable translational shift - if you supply this parameter you can not supply the maxshiftx, maxshifty or maxshiftz parameters. Each approach is mutually exclusive."); 01418 d.put("searchx", EMObject::INT,"The maximum length of the detectable translational shift in the x direction- if you supply this parameter you can not supply the maxshift parameters. Default is 3."); 01419 d.put("searchy", EMObject::INT,"The maximum length of the detectable translational shift in the y direction- if you supply this parameter you can not supply the maxshift parameters. Default is 3."); 01420 d.put("searchz", EMObject::INT,"The maximum length of the detectable translational shift in the z direction- if you supply this parameter you can not supply the maxshift parameters. Default is 3"); 01421 d.put("verbose", EMObject::BOOL,"Turn this on to have useful information printed to standard out."); 01422 return d; 01423 }
|
|
Definition at line 1405 of file aligner.h. 01406 { 01407 return new Refine3DAlignerGrid(); 01408 }
|
|
Definition at line 82 of file aligner.cpp. |