contrib/mul/vimt3d/vimt3d_save.cxx Source File

Go to the documentation of this file.
00001 // This is mul/vimt3d/vimt3d_save.cxx
00002 #include "vimt3d_save.h"
00003 //:
00004 // \file
00005 // \author Ian Scott
00006 
00007 #include <mbl/mbl_log.h>
00008 #include <vil3d/vil3d_image_resource.h>
00009 #include <vil3d/vil3d_new.h>
00010 #include <vil3d/vil3d_save.h>
00011 #include <vil3d/file_formats/vil3d_meta_image_format.h>
00012 #include <vimt3d/vimt3d_transform_3d.h>
00013 #include <vimt3d/vimt3d_image_3d.h>
00014 #include <vimt3d/vimt3d_vil3d_v3i.h>
00015 #include <vimt3d/vimt3d_vil3d_v3m.h>
00016 
00017 
00018 static mbl_logger& logger()
00019 {
00020   static mbl_logger l("mul.vimt3d.save");
00021   return l;
00022 }
00023 
00024 //: Create a transform from the properties of image resource.
00025 // The transform will be from world co-ordinates in metres to image co-ordinates (or mm if requested).
00026 void vimt3d_save_transform(vil3d_image_resource_sptr &ir,
00027                            const vimt3d_transform_3d& trans, 
00028                            bool use_millimetres /*=false*/)
00029 {
00030   if (dynamic_cast<vimt3d_vil3d_v3i_image *>(ir.ptr()) || 
00031     dynamic_cast<vimt3d_vil3d_v3m_image *>(ir.ptr()) )
00032   {
00033     vgl_vector_3d<double> vox_per_mm = trans.delta(vgl_point_3d<double>(0,0,0), 
00034                                                    vgl_vector_3d<double>(1.0, 1.0, 1.0));
00035 
00036     // get the translation component
00037     double tx =  trans.matrix()(0,3);
00038     double ty =  trans.matrix()(1,3);
00039     double tz =  trans.matrix()(2,3);
00040 
00041     vimt3d_transform_3d tr;
00042     //const double units_scaling = use_millimetres ? 1000.0 : 1.0;
00043     tr.set_zoom_only (1000.0*vox_per_mm.x(),
00044                       1000.0*vox_per_mm.y(),
00045                       1000.0*vox_per_mm.z(), tx,ty,tz );
00046 
00047     if (dynamic_cast<vimt3d_vil3d_v3i_image *>(ir.ptr()))
00048       static_cast<vimt3d_vil3d_v3i_image &>(*ir).set_world2im(tr);
00049     else
00050       static_cast<vimt3d_vil3d_v3m_image &>(*ir).set_world2im(tr);
00051 
00052   }
00053   else if (dynamic_cast<vil3d_meta_image *>(ir.ptr()) )
00054   {
00055     vgl_vector_3d<double> vx_size = trans.delta(vgl_point_3d<double>(0.0,0.0,0.0),
00056                                                 vgl_vector_3d<double>(1.0,1.0,1.0));
00057     double ox = trans.matrix()(0,3);
00058     double oy = trans.matrix()(1,3);
00059     double oz = trans.matrix()(2,3);
00060     static_cast<vil3d_meta_image *>(ir.ptr())->set_offset(ox,oy,oz,
00061                                                     1/vx_size.x(),1/vx_size.y(),1/vx_size.z());
00062   }
00063   else
00064   {
00065     vimt3d_transform_3d i2w=trans.inverse();
00066     vgl_vector_3d<double> dp = i2w.delta(vgl_point_3d<double> (0,0,0),
00067                                          vgl_vector_3d<double> (1.0, 1.0, 1.0));
00068     float scale = use_millimetres ? 1000.0f : 1.0f;
00069     if (!ir->set_voxel_size(float(dp.x())/scale,float(dp.y())/scale,float(dp.z())/scale))
00070       MBL_LOG(WARN, logger(), "vimt3d_save_transform(): Unable to include voxel sizes:"
00071               <<dp.x()<<','<<dp.y()<<','<<dp.z() );
00072   }
00073 }
00074 
00075 
00076 bool vimt3d_save(const vcl_string& path,
00077                  const vimt3d_image_3d& image,
00078                  bool use_millimetres /*=false*/)
00079 {
00080   const vimt3d_image_3d & iv = image;
00081   const vil3d_image_view_base & ib = iv.image_base();
00082 
00083   vil3d_image_resource_sptr ir = vil3d_new_image_resource(
00084     path.c_str(), ib.ni(), ib.nj(), ib.nk(), ib.nplanes(), ib.pixel_format(),
00085     vil3d_save_guess_file_format(path.c_str()));
00086 
00087   if (!ir)
00088     return false;
00089 
00090   ir->put_view(ib);
00091 
00092   vimt3d_save_transform(ir, image.world2im(), use_millimetres);
00093   return true;
00094 }