MNE-CPP  beta 1.0
mne_hemisphere.cpp
Go to the documentation of this file.
1 
2 //=============================================================================================================
37 //*************************************************************************************************************
38 //=============================================================================================================
39 // INCLUDES
40 //=============================================================================================================
41 
42 #include "mne_hemisphere.h"
43 
44 
45 //*************************************************************************************************************
46 //=============================================================================================================
47 // USED NAMESPACES
48 //=============================================================================================================
49 
50 using namespace MNELIB;
51 
52 
53 //*************************************************************************************************************
54 //=============================================================================================================
55 // DEFINE MEMBER METHODS
56 //=============================================================================================================
57 
59 : type(1)
60 , id(-1)
61 , np(-1)
62 , ntri(-1)
63 , coord_frame(-1)
64 , rr(MatrixX3f::Zero(0,3))
65 , nn(MatrixX3f::Zero(0,3))
66 , tris(MatrixX3i::Zero(0,3))
67 , nuse(-1)
68 , inuse(VectorXi::Zero(0))
69 , vertno(VectorXi::Zero(0))
70 , nuse_tri(-1)
71 , use_tris(MatrixX3i::Zero(0,3))
72 , nearest(VectorXi::Zero(0))
73 , nearest_dist(VectorXd::Zero(0))
74 , patch_inds(VectorXi::Zero(0))
75 , dist_limit(-1)
76 , dist(SparseMatrix<double>())
77 , tri_cent(MatrixX3d::Zero(0,3))
78 , tri_nn(MatrixX3d::Zero(0,3))
79 , tri_area(VectorXd::Zero(0))
80 , use_tri_cent(MatrixX3d::Zero(0,3))
81 , use_tri_nn(MatrixX3d::Zero(0,3))
82 , use_tri_area(VectorXd::Zero(0))
83 //, m_TriCoords()
84 //, m_pGeometryData(NULL)
85 {
86 }
87 
88 
89 //*************************************************************************************************************
90 
92 : type(p_MNEHemisphere.type)
93 , id(p_MNEHemisphere.id)
94 , np(p_MNEHemisphere.np)
95 , ntri(p_MNEHemisphere.ntri)
96 , coord_frame(p_MNEHemisphere.coord_frame)
97 , rr(p_MNEHemisphere.rr)
98 , nn(p_MNEHemisphere.nn)
99 , tris(p_MNEHemisphere.tris)
100 , nuse(p_MNEHemisphere.nuse)
101 , inuse(p_MNEHemisphere.inuse)
102 , vertno(p_MNEHemisphere.vertno)
103 , nuse_tri(p_MNEHemisphere.nuse_tri)
104 , use_tris(p_MNEHemisphere.use_tris)
105 , nearest(p_MNEHemisphere.nearest)
106 , nearest_dist(p_MNEHemisphere.nearest_dist)
107 , pinfo(p_MNEHemisphere.pinfo)
108 , patch_inds(p_MNEHemisphere.patch_inds)
109 , dist_limit(p_MNEHemisphere.dist_limit)
110 , dist(p_MNEHemisphere.dist)
111 , tri_cent(p_MNEHemisphere.tri_cent)
112 , tri_nn(p_MNEHemisphere.tri_nn)
113 , tri_area(p_MNEHemisphere.tri_area)
114 , use_tri_cent(p_MNEHemisphere.use_tri_cent)
115 , use_tri_nn(p_MNEHemisphere.use_tri_nn)
116 , use_tri_area(p_MNEHemisphere.use_tri_area)
117 , neighbor_tri(p_MNEHemisphere.neighbor_tri)
118 , neighbor_vert(p_MNEHemisphere.neighbor_vert)
119 , cluster_info(p_MNEHemisphere.cluster_info)
120 , m_TriCoords(p_MNEHemisphere.m_TriCoords)
121 {
122  //*m_pGeometryData = *p_MNEHemisphere.m_pGeometryData;
123 }
124 
125 
126 //*************************************************************************************************************
127 
129 {
130 
131 }
132 
133 
134 //*************************************************************************************************************
135 
137 {
138  type = 1;
139  id = -1;
140  np = -1;
141  ntri = -1;
142  coord_frame = -1;
143  rr = MatrixX3f::Zero(0,3);
144  nn = MatrixX3f::Zero(0,3);
145  tris = MatrixX3i::Zero(0,3);
146  nuse = -1;
147  inuse = VectorXi::Zero(0);
148  vertno = VectorXi::Zero(0);
149  nuse_tri = -1;
150  use_tris = MatrixX3i::Zero(0,3);
151  nearest = VectorXi::Zero(0);
152  nearest_dist = VectorXd::Zero(0);
153  pinfo.clear();
154  patch_inds = VectorXi::Zero(0);
155  dist_limit = -1;
156  dist = SparseMatrix<double>();
157  tri_cent = MatrixX3d::Zero(0,3);
158  tri_nn = MatrixX3d::Zero(0,3);
159  tri_area = VectorXd::Zero(0);
160  use_tri_cent = MatrixX3d::Zero(0,3);
161  use_tri_nn = MatrixX3d::Zero(0,3);
162  use_tri_area = VectorXd::Zero(0);
163 
164  neighbor_tri = QMap<int, QVector<int>>();
165  neighbor_vert = QMap<int, QVector<int>>();
166 
168 
169  m_TriCoords = MatrixXf();
170 }
171 
172 
173 //*************************************************************************************************************
174 
175 MatrixXf& MNEHemisphere::getTriCoords(float p_fScaling)
176 {
177  if(m_TriCoords.size() == 0)
178  {
179  m_TriCoords = MatrixXf(3,3*tris.rows());
180  for(qint32 i = 0; i < tris.rows(); ++i)
181  {
182  m_TriCoords.col(i*3) = rr.row( tris(i,0) ).transpose().cast<float>();
183  m_TriCoords.col(i*3+1) = rr.row( tris(i,1) ).transpose().cast<float>();
184  m_TriCoords.col(i*3+2) = rr.row( tris(i,2) ).transpose().cast<float>();
185  }
186  }
187 
188  m_TriCoords *= p_fScaling;
189 
190  return m_TriCoords;
191 }
192 
193 
194 //*************************************************************************************************************
195 
196 bool MNEHemisphere::transform_hemisphere_to(fiff_int_t dest, const FiffCoordTrans &p_Trans)
197 {
198  FiffCoordTrans trans(p_Trans);
199 
200  if (this->coord_frame == dest)
201  {
202 // res = src;
203  return true;
204  }
205 
206  if (trans.to == this->coord_frame && trans.from == dest)
207  trans.invert_transform();
208  else if(trans.from != this->coord_frame || trans.to != dest)
209  {
210  printf("Cannot transform the source space using this coordinate transformation");//Consider throw
211  return false;
212  }
213 
214  MatrixXf t = trans.trans.block(0,0,3,4);
215 // res = src;
216  this->coord_frame = dest;
217  MatrixXf t_rr = MatrixXf::Ones(this->np, 4);
218  t_rr.block(0, 0, this->np, 3) = this->rr;
219  MatrixXf t_nn = MatrixXf::Zero(this->np, 4);
220  t_nn.block(0, 0, this->np, 3) = this->nn;
221 
222  this->rr = (t*t_rr.transpose()).transpose();
223  this->nn = (t*t_nn.transpose()).transpose();
224 
225  return true;
226 }
227 
228 
229 //*************************************************************************************************************
230 //ToDo
232 {
233  if(this->type == 1 || this->type == 2)
234  p_pStream->write_int(FIFF_MNE_SOURCE_SPACE_TYPE, &this->type);
235  else
236  printf("Unknown source space type (%d)\n", this->type);
237  p_pStream->write_int(FIFF_MNE_SOURCE_SPACE_ID, &this->id);
238 
239 // data = this.get('subject_his_id', None)
240 // if data:
241 // write_string(fid, FIFF.FIFF_SUBJ_HIS_ID, data)
242  p_pStream->write_int(FIFF_MNE_COORD_FRAME, &this->coord_frame);
243 
244  if(this->type == 2) //2 = Vol
245  {
246  qDebug() << "ToDo: Write Volume not implemented jet!!!!!!!!";
247 // p_pStream->write_int(FIFF_MNE_SOURCE_SPACE_VOXEL_DIMS, this->shape)
248 // p_pStream->write_coord_trans(this->src_mri_t);
249 
250  p_pStream->start_block(FIFFB_MNE_PARENT_MRI_FILE);
251 // write_coord_trans(fid, this['vox_mri_t'])
252 
253 // write_coord_trans(fid, this['mri_ras_t'])
254 
255 // write_float_sparse_rcs(fid, FIFF.FIFF_MNE_SOURCE_SPACE_INTERPOLATOR,
256 // this['interpolator'])
257 
258 // if 'mri_file' in this and this['mri_file'] is not None:
259 // write_string(fid, FIFF.FIFF_MNE_SOURCE_SPACE_MRI_FILE,
260 // this['mri_file'])
261 
262 // write_int(fid, FIFF.FIFF_MRI_WIDTH, this['mri_width'])
263 // write_int(fid, FIFF.FIFF_MRI_HEIGHT, this['mri_height'])
264 // write_int(fid, FIFF.FIFF_MRI_DEPTH, this['mri_depth'])
265 
266  p_pStream->end_block(FIFFB_MNE_PARENT_MRI_FILE);
267  }
268 
269  p_pStream->write_int(FIFF_MNE_SOURCE_SPACE_NPOINTS, &this->np);
272 
273  // Which vertices are active
274  p_pStream->write_int(FIFF_MNE_SOURCE_SPACE_SELECTION, this->inuse.data(), this->inuse.size());
275  p_pStream->write_int(FIFF_MNE_SOURCE_SPACE_NUSE, &this->nuse);
276 
277  p_pStream->write_int(FIFF_MNE_SOURCE_SPACE_NTRI, &this->ntri);
278  if (this->ntri > 0)
279  p_pStream->write_int_matrix(FIFF_MNE_SOURCE_SPACE_TRIANGLES, this->tris.array() + 1);
280 
281  if (this->type != 2 && this->use_tris.rows() > 0)
282  {
283  // Use triangulation
285  p_pStream->write_int_matrix(FIFF_MNE_SOURCE_SPACE_USE_TRIANGLES, this->use_tris.array() + 1);
286  }
287 
288  // Patch-related information
289  if (this->nearest.size() > 0)
290  {
291  p_pStream->write_int(FIFF_MNE_SOURCE_SPACE_NEAREST, this->nearest.data(), this->nearest.size());
292  p_pStream->write_float_matrix(FIFF_MNE_SOURCE_SPACE_NEAREST_DIST, this->nearest_dist.cast<float>());
293  }
294 
295  // Distances
296  if (this->dist.rows() > 0)
297  {
298  // Save only upper triangular portion of the matrix
299  typedef Eigen::Triplet<float> T;
300  std::vector<T> tripletList;
301  tripletList.reserve(this->dist.nonZeros());
302  for (int k=0; k < this->dist.outerSize(); ++k)
303  for (SparseMatrix<double>::InnerIterator it(this->dist,k); it; ++it)
304  if(it.col() >= it.row())//only upper triangle -> todo iteration can be optimized
305  tripletList.push_back(T(it.row(), it.col(), (float)it.value()));
306  SparseMatrix<float> dists(this->dist.rows(), this->dist.cols());
307  dists.setFromTriplets(tripletList.begin(), tripletList.end());
308 
310  //ToDo check if write_float_matrix or write float is okay
311  p_pStream->write_float(FIFF_MNE_SOURCE_SPACE_DIST_LIMIT, &this->dist_limit); //p_pStream->write_float_matrix(FIFF_MNE_SOURCE_SPACE_DIST_LIMIT, this->dist_limit);
312  }
313 }
314 
315 
317 
318 //QGeometryData* MNEHemisphere::getGeometryData(float p_fScaling)
319 //{
320 // if(m_pGeometryData == NULL)
321 // {
322 // m_pGeometryData = new QGeometryData();
323 
324 // MatrixXd* triCoords = getTriCoords(p_fScaling);
325 
326 // m_pGeometryData->appendVertexArray(QArray<QVector3D>::fromRawData( reinterpret_cast<const QVector3D*>(triCoords->data()), triCoords->cols() ));
327 // }
328 
329 // return m_pGeometryData;
330 //}
void write_int_matrix(fiff_int_t kind, const MatrixXi &mat)
#define FIFF_MNE_SOURCE_SPACE_NPOINTS
#define FIFF_MNE_SOURCE_SPACE_NORMALS
#define FIFF_MNE_SOURCE_SPACE_NUSE
bool transform_hemisphere_to(fiff_int_t dest, const FiffCoordTrans &p_Trans)
#define FIFF_MNE_SOURCE_SPACE_NTRI
MatrixXf & getTriCoords(float p_fScaling=1.0f)
#define FIFF_MNE_SOURCE_SPACE_ID
#define FIFF_MNE_SOURCE_SPACE_TRIANGLES
MNEHemisphere class declaration.
#define FIFF_MNE_SOURCE_SPACE_POINTS
QMap< int, QVector< int > > neighbor_vert
MNEClusterInfo cluster_info
SparseMatrix< double > dist
void write_float(fiff_int_t kind, const float *data, fiff_int_t nel=1)
void write_int(fiff_int_t kind, const fiff_int_t *data, fiff_int_t nel=1)
void start_block(fiff_int_t kind)
#define FIFF_MNE_SOURCE_SPACE_NEAREST_DIST
Hemisphere provides geometry information.
#define FIFF_MNE_SOURCE_SPACE_USE_TRIANGLES
#define FIFF_MNE_COORD_FRAME
#define FIFF_MNE_SOURCE_SPACE_TYPE
#define FIFF_MNE_SOURCE_SPACE_NUSE_TRI
#define FIFF_MNE_SOURCE_SPACE_DIST
void write_float_matrix(fiff_int_t kind, const MatrixXf &mat)
Coordinate transformation description.
QMap< int, QVector< int > > neighbor_tri
void write_float_sparse_rcs(fiff_int_t kind, const SparseMatrix< float > &mat)
#define FIFF_MNE_SOURCE_SPACE_NEAREST
FIFF File I/O routines.
Definition: fiff_stream.h:129
void writeToStream(FiffStream *p_pStream)
void end_block(fiff_int_t kind)
#define FIFF_MNE_SOURCE_SPACE_SELECTION
#define FIFF_MNE_SOURCE_SPACE_DIST_LIMIT
QList< VectorXi > pinfo
Matrix< float, 4, 4, DontAlign > trans