1 /**********
2 *
3 * This file is part of the TINA Open Source Image Analysis Environment
4 * henceforth known as TINA
5 *
6 * TINA is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as
8 * published by the Free Software Foundation.
9 *
10 * TINA is distributed in the hope that it will be useful,
11 * but WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 * GNU General Public License for more details.
14 *
15 * You should have received a copy of the GNU General Public License
16 * along with TINA; if not, write to the Free Software Foundation, Inc.,
17 * 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
18 *
19 * ANY users of TINA who require exemption from the existing licence must
20 * negotiate a new licence with Dr. Neil.A.Thacker, the sole agent for
21 * the University of Manchester.
22 *
23 **********
24 *
25 * Program : TINA
26 * File :
27 * Date :
28 * Version :
29 * CVS Id :
30 *
31 * Author : mjs
32 *
33 * Notes :
34 *
35 *********
36 */
37
38 /**
39 * @file
40 * @brief Functions to perform Cerebral Blood Flow calculations
41 *
42 * The functions in this file perform Cerebral Blood Flow calculations (as detailed in
43 * Tina Memo 2002-001) for a volume of images, to produce normalised square-rooted cerebral blood
44 * volume (sqrt RCBV), 3D median filtered time of arrival (TTM) , net mean transit time (NMTT),
45 * and net cerebral blood flow (NCBF).
46 * In addition, there is also the capability to produce regional estimates of these parameters
47 * based on the Talairach stereotaxic atlas.
48 */
49
50
51 #include "tlmedFlow_tool.h"
52
53 #if HAVE_CONFIG_H
54 #include <config.h>
55 #endif
56
57 #include <stdio.h>
58 #include <sys/param.h>
59 #include <string.h>
60 #include <math.h>
61
62 #include <tina/sys/sysDef.h>
63 #include <tina/sys/sysPro.h>
64 #include <tina/math/mathDef.h>
65 #include <tina/math/mathPro.h>
66 #include <tina/image/imgDef.h>
67 #include <tina/image/imgPro.h>
68 #include <tina/geometry/geomDef.h>
69 #include <tina/geometry/geomPro.h>
70 #include <tina/medical/medDef.h>
71 #include <tina/medical/medPro.h>
72 #include <tinatool/draw/drawDef.h>
73 #include <tinatool/draw/drawPro.h>
74 #include <tinatool/wdgts/wdgtsDef.h>
75 #include <tinatool/wdgts/wdgtsPro.h>
76 #include <tinatool/tlvision/tlvisEdge_epi_mouse.h>
77 #include <tinatool/tlbase/tlbaseDef.h>
78 #include <tinatool/tlbase/tlbasePro.h>
79 #include <tinatool/tlmedical/tlmedPro.h>
80 #include <tinatool/tlmedical/tlmedDef.h>
81
82 static Sequence *cbv_seq = NULL;
83 static Sequence *ttm_seq = NULL;
84 static Sequence *flow_seq = NULL;
85 static Sequence *nmtt_seq = NULL;
86 static Sequence *smoothed_seq = NULL;
87 static Sequence *nmtt_mask_seq = NULL;
88
89 static shistogram **shist_vec = NULL, **shist_vec_left = NULL, **shist_vec_right = NULL;
90 static int hist_no = 1;
91 static float *imptrs = NULL;
92 static void *ptalno = NULL;
93 static int hist_choice = 0;
94
95
96 /** @brief Assigns the current sequence as the CBV volume
97 * @param void
98 * @return void
99 *
100 * Copies the current sequence, as shown in the sequence tool to the
101 * CBV sequence. If a CBV sequence already exists, this is deleted first.
102 * The images in the current sequence should be the "raw" CBV images, as would be obtained
103 * using the TINA NMR analysis tool.
104 */
105
106 static void get_cbv_vol_proc(void)
107 {
108 Sequence *seq = NULL;
109 List *lptr = NULL;
110
111 if ((seq = seq_get_current()) == NULL)
112 {
113 format("No sequence available\n");
114 return;
115 }
116
117 if (cbv_seq != NULL)
118 {
119 seq_rm(cbv_seq);
120 cbv_seq = NULL;
121 }
122 cbv_seq = seq_copy(seq);
123 }
124
125
126 /** @brief Assigns the current sequence to be the TTM sequence
127 * @param void
128 * @return void
129 *
130 * Copies the current sequence, as shown in the sequence tool to the
131 * TTM sequence. If a TTM sequence already exists, the TTM, NMTT, and 3D median smoothed TTM
132 * sequences are all deleted first.
133 * The images in the current sequence should be the raw TTM images as would be obtained using the
134 * TINA NMR analysis software
135 */
136 static void get_ttm_vol_proc(void)
137 {
138 Sequence *seq = NULL;
139 List *lptr = NULL;
140
141 if ((seq = seq_get_current()) == NULL)
142 {
143 format("No sequence available\n");
144 return;
145 }
146
147 if (ttm_seq != NULL)
148 {
149 seq_rm(ttm_seq);
150 seq_rm(nmtt_seq);
151 seq_rm(smoothed_seq);
152 ttm_seq = NULL;
153 nmtt_seq = NULL;
154 smoothed_seq = NULL;
155 }
156
157
158 ttm_seq = seq_copy(seq);
159 }
160
161 /** @brief Processes the CBV sequence to obtain the normalised square-rooted RCBV maps
162 * @param void
163 * @return void
164 *
165 * The function initially replaces each image in the sequence with its square-root, and then
166 * produces a maximum intensity projection (MIP) of the images. The MIP undergoes tangential
167 * smoothing, and the maximum voxel value in the image found. This voxel is assumed to be
168 * 100% blood, so the square-root CBV volume is normalised to this voxel value to give % maps of
169 * blood volume.
170 * The MIP is pushed onto the stack (visible on Imcalc TV) and the current sequence is set to the
171 * sqrt RCBV maps (visible on the Sequence Tv).
172 */
173 static void process_cbv_proc(void)
174 {
175 List *lptr = NULL;
176 Imrect *mip_im = NULL;
177 Imrect *mip_im2 = NULL;
178 double max_val = 0.0;
179 int x, y;
180
181 if (cbv_seq == NULL)
182 {
183 format("No CBV sequence\n");
184 return;
185 }
186
187 for (lptr = cbv_seq->start; lptr != NULL; lptr = lptr->next)
188 {
189 Imrect *seq_im = NULL, *seq_im2 = NULL;
190
191 seq_im = lptr->to;
192 seq_im2 = im_sqrt(seq_im);
193 lptr->to = seq_im2;
194 im_free(seq_im);
195
196 if (lptr == cbv_seq->start)
197 {
198 mip_im = im_maxsel(seq_im2, seq_im2);
199 }
200 else
201 {
202 mip_im2 = im_maxsel(mip_im, seq_im2);
203 im_free(mip_im);
204 mip_im = im_copy(mip_im2);
205 im_free(mip_im2);
206 mip_im2 = NULL;
207 }
208 }
209
210 mip_im2 = im_tsmooth(mip_im);
211 im_free(mip_im);
212 max_val = im_locate_max(mip_im2, &y, &x);
213 stack_push((void *)mip_im2, IMRECT, im_free);
214 format("max value %5.3f \n",max_val);
215
216 if (max_val > 2000.0)
217 {
218 format("Warning: max sqrt(CBV) value is greater than 2000,\n check there are no outliers in the data\n");
219 }
220
221
222 for (lptr = cbv_seq->start; lptr != NULL; lptr = lptr->next)
223 {
224 Imrect *seq_im = NULL, *seq_im2 = NULL;
225
226 seq_im = lptr->to;
227 seq_im2 = im_times((1.0/max_val), seq_im);
228 lptr->to = seq_im2;
229 im_free(seq_im);
230
231 }
232 seq_set_current(cbv_seq);
233
234
235 }
236
237 /** @brief Processes the TTM maps to produce the NMTT maps
238 * @param void
239 * @return void
240 *
241 * Initially the TTM images are thresholded at 0.0 to remove spurious negative values.
242 * The TTM maps are then 3d median smoothed, using the slice of interest and the adjacent two slices.
243 * The smoothed images are maintained in a separate sequence (smoothed_seq).
244 * The NMTT maps are then produced using the spatial differentiation in 3D of the smoothed
245 * TTM maps. Note that the NMTT volume is 2 slices shorter than the TTM volumes.
246 */
247 static void process_ttm_proc(void)
248 {
249 /* 3d median smooths entire volume
250 nb, need to maintain two sequences here, and delete original
251 at end of conversion
252 calculate 3d spatial derivative.nb, two less (? or is it 4) images
253 in sequence than before.
254 */
255
256 List *lptr = NULL;
257 List *smoothed_images = NULL, *nmtt_images = NULL;
258 double mtt_scale = 0.0;
259
260 if (ttm_seq == NULL)
261 {
262 format("No TTM sequence\n");
263 return;
264 }
265
266 /* threshold images first */
267
268 for (lptr = ttm_seq->start; lptr != NULL; lptr = lptr->next)
269 {
270 Imrect *thresh_im = NULL, *temp = NULL;
271
272 temp = lptr->to;
273 thresh_im = im_thresh(0.0, temp);
274 lptr->to = thresh_im;
275 im_free(temp);
276
277 }
278
279
280 /* 3d median smooth - produces a new sequence */
281
282 for (lptr = ttm_seq->start; lptr != NULL; lptr = lptr->next)
283 {
284 Imrect *im_new = NULL, *im_bef = NULL, *im_aft = NULL, *im_cen = NULL;
285 List *smoothed_link = NULL;
286
287 im_cen = lptr->to;
288 /* The following if statement deals with the cases at the start and end of the sequence
289 where there are no below and above slices respectively. In these cases the slice of
290 interest replaces the missing slice */
291 if ((lptr->last != NULL) && (lptr->last->to != NULL))
292 {
293 im_bef = lptr->last->to;
294 }
295 else
296 {
297 im_bef = lptr->to;
298 }
299 if ((lptr->next != NULL) && (lptr->next->to != NULL))
300 {
301 im_aft = lptr->next->to;
302 }
303 else
304 {
305 im_aft = lptr->to;
306 }
307 im_new = im_median_3D(im_bef, im_cen, im_aft);
308
309
310 smoothed_link = link_alloc(im_new, IMRECT);
311 smoothed_images = list_addtoend(smoothed_images, smoothed_link);
312 }
313 /* can be deleted later, just for debugging */
314 smoothed_seq = seq_alloc();
315 smoothed_seq->start = smoothed_images;
316 smoothed_seq->current = smoothed_images;
317 smoothed_seq->end = list_get_end(smoothed_images);
318
319 /* seq_set_current(smoothed_seq);*/
320
321 /* calculate 3d spatial derivative (creates a new sequence of NMTT maps)*/
322 /* should also create the relevant masks */
323 /* 1.79688/(3*2) */
324
325 mtt_scale = (ttm_seq->dim[0])/(2*ttm_seq->dim[2]);
326
327 if (smoothed_images->next->to == NULL)
328 {
329 format("2nd image in seq is null");
330 return;
331 }
332
333 for (lptr = smoothed_images->next; lptr->next != NULL; lptr = lptr->next)
334 {
335 /* bear in mind the sequence will now be two images shorter (2-24)*/
336 Imrect * im_cen = NULL, *im_bef = NULL, *im_aft = NULL;
337 Imrect *im_grad_h = NULL, *im_grad_v = NULL, *im_diff_z = NULL, *im_z = NULL;
338 Imrect *im_gradsq = NULL, *im_grad_z_sq = NULL, *im_gradsq_tot = NULL;
339 Imrect *im_nmtt = NULL;
340 List *nmtt_link = NULL;
341
342 im_cen = lptr->to;
343 im_bef = lptr->last->to;
344 im_aft = lptr->next->to;
345
346 im_grad_h = imf_grad_h(im_cen);
347 im_grad_v = imf_grad_v(im_cen);
348 im_gradsq = imf_sumsq(im_grad_h, im_grad_v);
349
350 im_diff_z = im_diff(im_aft, im_bef);
351 im_z = im_times(mtt_scale, im_diff_z);
352 im_grad_z_sq = imf_prod(im_z, im_z);
353
354 im_gradsq_tot = im_sum(im_gradsq, im_grad_z_sq);
355 im_nmtt = im_sqrt(im_gradsq_tot);
356 nmtt_link = link_alloc(im_nmtt, IMRECT);
357 nmtt_images = list_addtoend(nmtt_images, nmtt_link);
358
359 im_free(im_grad_h);
360 im_free(im_grad_v);
361 im_free(im_gradsq);
362 im_free(im_diff_z);
363 im_free(im_z);
364 im_free(im_grad_z_sq);
365 im_free(im_gradsq_tot);
366
367
368 }
369
370 nmtt_seq = seq_alloc();
371 nmtt_seq->start = nmtt_images;
372 nmtt_seq->current = nmtt_images;
373 nmtt_seq->end = list_get_end(nmtt_images);
374 nmtt_seq->dim[0] = ttm_seq->dim[0];
375 nmtt_seq->dim[1] = ttm_seq->dim[1];
376 nmtt_seq->dim[2] = ttm_seq->dim[2];
377 nmtt_seq->dim[3] = ttm_seq->dim[3];
378 nmtt_seq->offset = ttm_seq->offset+1;
379 nmtt_seq->stride = 1;
380 seq_set_current(nmtt_seq);
381
382
383 }
384
385 /** @brief Produces masks for applying to the NMTT maps to remove edge artefacts
386 * @param void
387 * @return void
388 *
389 * The NMTT maps have edge artefacts particularly at the boundary between the brain
390 * and background. Masks to remove these artefacts can be created using a convolution of
391 * the binary masks of the smoothed TTM maps used to create the NMTT maps.
392 * More specifically, the TTM map at the slice of interest is binarised and barrel shifted
393 * by one pixel in all directions to produce 4 masks. These are convolved, along with masks
394 * slice above and below the slice of interest.
395 * The NMTT masks are set to be the current sequence (visible in the Sequence Tool Tv).
396 * NB, the volume of NMTT masks is 2 images shorter than the volume of TTM images used in
397 * their creation.
398 */
399 static void get_nmtt_masks_proc(void)
400 {
401 List *lptr = NULL, *nmtt_mask_list = NULL;
402
403
404 if (nmtt_mask_seq != NULL)
405 {
406 seq_rm(nmtt_mask_seq);
407 nmtt_mask_seq = NULL;
408 }
409
410 if (smoothed_seq == NULL)
411 return;
412
413 if (smoothed_seq->start->next == NULL)
414 return;
415
416 for (lptr=smoothed_seq->start->next; lptr->next != NULL; lptr = lptr->next)
417 {
418 List *nmtt_mask_link = NULL;
419 Imrect *befa = NULL, *befb = NULL, *afta = NULL, *aftb = NULL;
420 Imrect *left = NULL, *right = NULL;
421 Imrect *up = NULL, *down = NULL;
422 Imrect *cena = NULL, *cenb = NULL;
423 Imrect *temp1=NULL, *temp2=NULL, *nmtt_mask = NULL;
424
425 befa = im_bthresh(0.0, lptr->last->to);
426 afta = im_bthresh(0.0, lptr->next->to);
427 cena = im_bthresh(0.0, lptr->to);
428
429 befb = im_cast(befa, uchar_v);
430 aftb = im_cast(afta, uchar_v);
431 cenb = im_cast(cena, uchar_v);
432
433 left = im_bshift(cenb, -1, 0);
434 right = im_bshift(cenb, 1, 0);
435 up = im_bshift(cenb, 0, -1);
436 down = im_bshift(cenb, 0, 1);
437
438 temp1 = im_prod(befb, aftb);
439 temp2 = im_prod(temp1, left); im_free(temp1);
440 temp1 = im_prod(temp2, right); im_free(temp2);
441 temp2 = im_prod(temp1, up); im_free(temp1);
442 nmtt_mask = im_prod(temp2, down); im_free(temp2);
443
444 im_free(befa); im_free(befb);
445 im_free(afta); im_free(aftb);
446 im_free(cena); im_free(cenb);
447 im_free(left); im_free(right);
448 im_free(up); im_free(down);
449
450 nmtt_mask_link = link_alloc(nmtt_mask, IMRECT);
451 nmtt_mask_list = list_addtoend(nmtt_mask_list, nmtt_mask_link);
452
453 }
454
455 nmtt_mask_seq = seq_alloc();
456 nmtt_mask_seq->start = nmtt_mask_list;
457 nmtt_mask_seq->current = nmtt_mask_list;
458 nmtt_mask_seq->end = list_get_end(nmtt_mask_list);
459 nmtt_mask_seq->dim[0] = ttm_seq->dim[0];
460 nmtt_mask_seq->dim[1] = ttm_seq->dim[1];
461 nmtt_mask_seq->dim[2] = ttm_seq->dim[2];
462 nmtt_mask_seq->dim[3] = ttm_seq->dim[3];
463 nmtt_mask_seq->offset = ttm_seq->offset+1;
464 nmtt_mask_seq->stride = 1;
465 seq_set_current(nmtt_mask_seq);
466
467 }
468
469 /** @brief Produces the log Net Cerebral Blood Flow images
470 * @param void
471 * @return void
472 *
473 * Calculates the log(NCBF) as described in Tina Memo 2002-001.
474 * The log(NCBF) sequence is set to the current sequence as visible in
475 * the Sequence tool Tv.
476 */
477 static void get_flow_vol_proc(void)
478 {
479
480 List *nmtt_lptr = NULL, *cbv_lptr = NULL, *mask_lptr = NULL;
481 double scale_factor = 0.0;
482 double voxel_length = 0.0;
483 List *flow_list = NULL;
484
485
486 if (flow_seq != NULL)
487 {
488 seq_rm(flow_seq);
489 flow_seq = NULL;
490 }
491 if ((nmtt_seq == NULL) || (cbv_seq == NULL) || (nmtt_mask_seq == NULL))
492 return;
493
494 /* scale factor is 60secs*voxel length* (M/rho)^2/3 */
495 voxel_length = ttm_seq->dim[0]/10.0;
496 scale_factor = 60*voxel_length*20.4684;
497
498 for (nmtt_lptr = nmtt_seq->start, mask_lptr = nmtt_mask_seq->start, cbv_lptr = cbv_seq->start->next; nmtt_lptr != NULL; nmtt_lptr = nmtt_lptr->next, mask_lptr = mask_lptr->next, cbv_lptr = cbv_lptr->next)
499 {
500 List *flow_link = NULL;
501 Imrect *flow_im = NULL, *nmtt_im = NULL, *cbv_im = NULL, *mask_im = NULL;
502 Imrect *temp1 = NULL, *temp2 = NULL, *cbv_sq = NULL, *cbv_sq_mask = NULL;
503 Imrect *log_flow_im = NULL, *thresh_log_flow_im = NULL;
504 double noise = 0.0, sigma_x = 0.0, sigma_y=0.0;
505
506 nmtt_im = nmtt_lptr->to;
507 cbv_im = cbv_lptr->to;
508
509 /* misnomer in nomenclature, the cbv_im is actually the sqrt cbv image */
510 mask_im = mask_lptr->to;
511 cbv_sq = im_prod(cbv_im, cbv_im);
512
513 cbv_sq_mask = im_prod(cbv_sq, mask_im);
514 temp1 = im_prod(nmtt_im, mask_im);
515
516 sigma_x = imf_diffx_noise(temp1, nmtt_im->region);
517 sigma_y = imf_diffy_noise(temp1, nmtt_im->region);
518 noise = (sigma_x + sigma_y)/2;
519
520 temp2 = im_div(cbv_sq_mask, temp1, noise, noise);
521 flow_im = im_times(scale_factor, temp2);
522 log_flow_im = im_log(flow_im);
523 thresh_log_flow_im = im_thresh(-10.0, log_flow_im);
524
525 im_free(temp1);
526 im_free(temp2);
527 im_free(cbv_sq);
528 im_free(cbv_sq_mask);
529 im_free(flow_im);
530 im_free(log_flow_im);
531 flow_link = link_alloc(thresh_log_flow_im, IMRECT);
532 flow_list = list_addtoend(flow_list, flow_link);
533
534 }
535 flow_seq = seq_alloc();
536 flow_seq->start = flow_list;
537 flow_seq->current = flow_list;
538 flow_seq->end = list_get_end(flow_list);
539 flow_seq->dim[0] = ttm_seq->dim[0];
540 flow_seq->dim[1] = ttm_seq->dim[1];
541 flow_seq->dim[2] = ttm_seq->dim[2];
542 flow_seq->dim[3] = ttm_seq->dim[3];
543 flow_seq->offset = ttm_seq->offset+1;
544 flow_seq->stride = 1;
545 seq_set_current(flow_seq);
546
547 }
548
549 /** @brief Creates (initialises) histograms for each Talairach region
550 * @param void
551 * @return void
552 *
553 * Initialises the histograms that will contain the regional (Taliarach) values
554 * of the flow maps. There are three sets of histograms corresponding to the whole brain
555 * and the left and right hemispheres.
556 *
557 */
558 static void init_tal_proc(void)
559 {
560
561
562 int i;
563
564 if (shist_vec != NULL)
565 {
566 for (i=0; i<153; i++)
567 {
568 hfree(shist_vec[i]);
569 }
570
571 }
572
573 if (shist_vec_left != NULL)
574 {
575 for (i=0; i<153; i++)
576 {
577 hfree(shist_vec_left[i]);
578 }
579
580 }
581
582 if (shist_vec_right != NULL)
583 {
584 for (i=0; i<153; i++)
585 {
586 hfree(shist_vec_right[i]);
587 }
588
589 }
590
591 shist_vec = (shistogram **)pvector_alloc(0, 153);
592 shist_vec_left = (shistogram **)pvector_alloc(0, 153);
593 shist_vec_right = (shistogram **)pvector_alloc(0, 153);
594
595 for (i=0; i<153; i++)
596 {
597 shistogram *tal_hist = NULL;
598
599 tal_hist = hbook1(i, pLexRecords[i], -10, 10, 50);
600 shist_vec[i] = tal_hist;
601 /*format("%d %s\n", i, pLexRecords[i]);*/
602 }
603
604 for (i=0; i<153; i++)
605 {
606 shistogram *tal_hist_left = NULL;
607 shistogram *tal_hist_right = NULL;
608
609 tal_hist_left = hbook1(i, pLexRecords[i], -10, 10, 50);
610 tal_hist_right = hbook1(i, pLexRecords[i], -10, 10, 50);
611 shist_vec_left[i] = tal_hist_left;
612 shist_vec_right[i] = tal_hist_right;
613
614 }
615
616 /* ideally also want something here to initialise flow to tal atlas in coreg tool if not already performed. */
617
618 return;
619
620
621
622 }
623
624 /** @brief Fills previously created regional histograms with log(NCBF) values
625 * @param void
626 * @return void
627 *
628 * Assumes that the volume of interest has already been registered to the Talairach
629 * template, or that the correct AIR file has been loaded. Requires the Coreg tool to
630 * be opened.
631 * The function iterates over all the voxels in the log(NCBF) volume and projects the positions
632 * onto the Talairach atlas to determine which region each voxel falls in. If the voxel is
633 * within the brain (not background) its value is added to the appropriate histograms. There
634 * are histograms covering a hierarchy of subdivisions, from lobar down through to the level of
635 * Brodmann area, and both left and right hemisphere and whole-brain histograms exist for all areas.
636 */
637 static void hist_flow_proc(void)
638 {
639 /* scan entire volumes
640 take flow value and input into correct histogram*/
641 float vox_val = 0.0;
642 Vec3 talpos, start_pos, tal3;
643 ListEntry *list_entry = NULL;
644 int i, j, k;
645 int seq_end = 0;
646 Imrect *im = NULL;
647 Imregion *imreg = NULL;
648 List *lptr = NULL;
649 int lx, ux, ly, uy;
650
651 if (flow_seq == NULL)
652 {
653 format("There's no flow sequence \n");
654 return;
655 }
656
657 format("Please ensure you have registered the flow volume to the Taliarach atlas or loaded in the appropriate air file.\n");
658
659 seq_end = get_end_frame(flow_seq);
660 lptr = flow_seq->start;
661 im = lptr->to;
662 imreg = im->region;
663 lx = imreg->lx; ly = imreg->ly; ux = imreg->ux; uy = imreg->uy;
664
665 for (k = flow_seq->offset, lptr = flow_seq->start; k < seq_end+1, lptr != NULL; k++, lptr = lptr->next)
666 {
667 im = lptr->to;
668
669 for (i = lx; i < ux; i++)
670 {
671
672 for (j = ly; j < uy; j++)
673 {
674
675 talpos = vec3(i, j, k); /*where ij and k are the voxel positions might need to add 0.5 to all */
676 vox_val = im_get_pixf(im, j, i); /* I HOPE :) */
677 /* I think because we've got the actual position that we don't need to add 0.5? */
678 tal3 = coreg_proj(talpos.el[0], talpos.el[1], talpos.el[2]+0.5);
679 /*!! this might prove problematic. need to do coreg_init or something, but I'm not sure which data set it needs to be done on */
680 tal3 = tal_convert(tal3);
681 list_entry = hist_fill_search(tal3);
682
683 if ((list_entry != NULL) && (vox_val != 0.0000))
684 {
685 int gross = 0, lobe = 0, sublobar = 0, record = 0;
686 int sublob1 = 0, sublob2 = 0;
687 gross = list_entry->field1;
688 hfill1(shist_vec[gross], vox_val, 1.0);
689
690 lobe = list_entry->field2;
691 hfill1(shist_vec[lobe], vox_val, 1.0);
692
693 sublobar = list_entry->field3;
694 hfill1(shist_vec[sublobar], vox_val, 1.0);
695
696 sublob1 = list_entry->field4;
697 hfill1(shist_vec[sublob1], vox_val, 1.0);
698
699 sublob2 = list_entry->field5;
700 hfill1(shist_vec[sublob2], vox_val, 1.0);
701
702 if (lobe >7 && lobe <53 && sublobar > 7 && sublobar <53)
703 {
704 if (gross == 2) /* ie, in left cerebrum */
705 {
706 hfill1(shist_vec_left[lobe], vox_val, 1.0);
707 hfill1(shist_vec_left[sublobar], vox_val, 1.0);
708 hfill1(shist_vec_left[sublob1], vox_val, 1.0);
709 hfill1(shist_vec_left[sublob2], vox_val, 1.0);
710 }
711 if (gross == 3) /*ie, in right cerebrum */
712 {
713 hfill1(shist_vec_right[lobe], vox_val, 1.0);
714 hfill1(shist_vec_right[sublobar], vox_val, 1.0);
715 hfill1(shist_vec_right[sublob1], vox_val, 1.0);
716 hfill1(shist_vec_right[sublob2], vox_val, 1.0);
717 }
718 }
719
720 }
721 else
722 /*format("no desc here\n");*/
723 {}
724 }
725 }
726 }
727 return;
728 }
729
730 /** @brief Sets the value of the global hist_choice variable
731 * @param int val The value determining histogram type
732 * @return void
733 *
734 * Sets the value (0, 1, or 2) of the global variable determining histogram type
735 * (whole-brain, left or right hemisphere).
736 */
737 static void hist_choice_proc(int val)
738 {
739 hist_choice = val;
740
741 }
742
743 /** @brief Plots the user-determined regional histogram of log(NCBF) values
744 * @param void
745 * @return void
746 *
747 * Takes the histopgram type (as specified by "Hist: Whole/Left/Right" on the user
748 * interface) and histogram number (as specified by "Hist no :" on the user interface)
749 * and plots the corresponding histogram on the Imcalc Graph Tv. Also provides the
750 * region name and histogram mean in the text box of the TinaTool.
751 */
752 static void hist_plot_proc(void)
753 {
754 Tv *graph_tv = NULL;
755 float hist_mean = 0.0;
756 shistogram **shist_choice = NULL;
757
758 if ((graph_tv = imcalc_graph_tv_get()) == NULL)
759 {
760 format("Open and install Graph Tv first\n");
761 return;
762 }
763
764 if ((shist_vec == NULL) || (shist_vec_left == NULL) || (shist_vec_right == NULL))
765 {
766 format("No histograms exist\n");
767 return;
768 }
769
770 if ((hist_no >152) ||(hist_no < 1))
771 {
772 format("No such histogram %d\n", hist_no);
773 return;
774 }
775
776 if (hist_no >=8 && hist_no <153)
777 {
778 switch (hist_choice)
779 {
780 case 0:
781 {
782 shist_choice = shist_vec;
783 break;
784 }
785 case 1:
786 {
787 shist_choice = shist_vec_left;
788 break;
789 }
790 case 2:
791 {
792 shist_choice = shist_vec_right;
793 break;
794 }
795 default:
796 shist_choice = shist_vec;
797 }
798 }
799
800 else
801 {
802 shist_choice = shist_vec;
803 }
804 graph_hfit(graph_tv, shist_choice[hist_no]);
805 hist_mean = shist_choice[hist_no]->mean;
806 format("Histogram of: %s ", pLexRecords[hist_no]);
807 format("mean: %f\n", hist_mean);
808
809 return;
810 }
811
812
813 /** @brief Code for the flow tool interface
814 * @param int x tool x-dimension
815 * @param int y tool y-dimension
816 * @return void
817 *
818 *
819 */
820 void flow_tool(int x, int y)
821 {
822 static void *tool = NULL;
823
824 if (tool)
825 {
826 tw_show_tool(tool);
827 return;
828 }
829
830 /*tv_proc();*/
831
832 tool = (void *)tw_tool("Flow Calculation Tool", x, y);
833
834 tw_label("Initialise Image Volumes");
835 tw_newrow();
836 tw_button("CBV volume", get_cbv_vol_proc, NULL);
837 tw_button("TTM volume", get_ttm_vol_proc, NULL);
838 tw_newrow();
839 tw_label("Perform Flow Calculations");
840 tw_newrow();
841 tw_button("Process CBV", process_cbv_proc, NULL);
842 tw_label(" -> ");
843 tw_button("process TTM", process_ttm_proc, NULL);
844 tw_label(" -> ");
845 tw_button("NMTT Masks", get_nmtt_masks_proc, NULL);
846 tw_label(" -> ");
847 tw_button("Get Flow Vol", get_flow_vol_proc, NULL);
848 tw_newrow();
849 tw_label("Initialise Regional Histograms");
850 tw_newrow();
851 tw_button("Init Hist.", init_tal_proc, NULL);
852 tw_button("Calc Flow Hists", hist_flow_proc, NULL);
853 tw_newrow();
854 tw_choice("Hist: ", hist_choice_proc, 0, "Whole", "Left", "Right", NULL);
855 ptalno = (void *)tw_iglobal(" Hist no: ", &hist_no, 10);
856 tw_button("Plot Hist", hist_plot_proc, NULL);
857
858 tw_end_tool();
859 }
860
861
862
863
864
865
866
This page was automatically generated by the
LXR engine.
Visit the LXR main site for more
information.