1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
|
#include "filter_resize.h"
void FilterResize::init() {
for(int i=0; i<3; i++) {
newsize[i].set_description(STD_string(dataDimLabel[1+i])+"-size");
append_arg(newsize[i],"newsize"+itos(i));
}
}
bool FilterResize::process(Data<float,4>& data, Protocol& prot) const {
TinyVector<int,4> oldshape(data.shape());
TinyVector<int,4> newshape(data.shape());
for(int i=0; i<3; i++) newshape(1+i)=newsize[i];
data.congrid(newshape);
// update stuff in the protocol after regridding
prot.seqpars.set_MatrixSize(phaseDirection,newshape(phaseDim));
prot.seqpars.set_MatrixSize(readDirection,newshape(readDim));
if(prot.geometry.get_Mode()==slicepack) {
prot.geometry.set_nSlices(newshape(sliceDim));
float distfactor=secureDivision(oldshape(sliceDim), newshape(sliceDim));
prot.geometry.set_sliceDistance(distfactor*prot.geometry.get_sliceDistance());
} else {
prot.seqpars.set_MatrixSize(sliceDirection,newshape(sliceDim));
}
return true;
}
////////////////////////////////////////////////////////////////////////////////
void FilterResample::init() {
newsize.set_description("new size");
append_arg(newsize,"newsize");
}
bool FilterResample::process(Data<float,4>& data, Protocol& prot) const {
double totaldur=data.extent(timeDim)*prot.seqpars.get_RepetitionTime();
TinyVector<int,4> newshape(data.shape());
newshape(0)=newsize;
data.congrid(newshape);
// update stuff in the protocol after regridding
prot.seqpars.set_NumOfRepetitions(newsize);
prot.seqpars.set_RepetitionTime(secureDivision(totaldur,newsize));
return true;
}
////////////////////////////////////////////////////////////////////////////////
void FilterIsotrop::init() {
size=0;
size.set_description("voxelsize").set_unit(ODIN_SPAT_UNIT);
append_arg(size,"voxelsize");
}
bool FilterIsotrop::process(Data<float,4>& data, Protocol& prot) const {
Log<Filter> odinlog(c_label(),"process");
Geometry &geo=prot.geometry;
TinyVector<int,4> shape(data.shape()),newshape=shape;
TinyVector<float,3> vox_mult;
TinyVector<float,3> vox_size(
FileFormat::voxel_extent(geo,sliceDirection, shape(sliceDim)),
FileFormat::voxel_extent(geo,phaseDirection, shape(phaseDim)),
FileFormat::voxel_extent(geo,readDirection, shape(readDim))
);
float m;
if(size>0) m=size;
else m=min(vox_size);
vox_mult=vox_size/m;
ODINLOG(odinlog,normalDebug) << "Regridding to isotropic " << m << ODIN_SPAT_UNIT << STD_endl;
ODINLOG(odinlog,normalDebug) << "current voxel_size: " << vox_size << STD_endl;
ODINLOG(odinlog,normalDebug) << "shape scaling vector: " << vox_mult << STD_endl;
for(dataDim d=readDim;d>timeDim;d=dataDim(d-1))
newshape(d)=(int)(vox_mult(d-1)*shape(d));
ODINLOG(odinlog,normalDebug) << "new shape: " << newshape << STD_endl;
data.congrid(newshape);
//prevent update (called by set_nSlices) from overriding FOVslice
if(geo.get_Mode()==slicepack) {
geo.set_sliceThickness(m);
geo.set_sliceDistance(m);
}
if(geo.get_Mode()==voxel_3d) {
geo.set_FOV(sliceDirection,m*newshape(sliceDim));
}
// update stuff in the protocol after regridding
geo.set_nSlices(newshape(sliceDim));
prot.seqpars.set_MatrixSize(phaseDirection,newshape(phaseDim));
prot.seqpars.set_MatrixSize(readDirection,newshape(readDim));
ODINLOG(odinlog,normalDebug) << "new voxel_size: " << TinyVector<float,3>(
FileFormat::voxel_extent(geo,sliceDirection, data.shape()(sliceDim)),
FileFormat::voxel_extent(geo,phaseDirection, data.shape()(phaseDim)),
FileFormat::voxel_extent(geo,readDirection, data.shape()(readDim))) << STD_endl;
return true;
}
|
