Skip to content

Commit

Permalink
Adding inferred time seg option (#12)
Browse files Browse the repository at this point in the history
  • Loading branch information
@acerjanic
acerjanic authored Jul 21, 2019
1 parent df03db5 commit 7ffdb48
Show file tree
Hide file tree
Showing 3 changed files with 72 additions and 61 deletions.
31 changes: 19 additions & 12 deletions PowerGrid/MPI/PowerGridPcSenseMPI_TS.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -76,18 +76,6 @@ int main(int argc, char** argv)
return 1;
}

if (TimeSegmentationInterp.compare("hanning") == 0) {
type = 1;
} else if (TimeSegmentationInterp.compare("minmax") == 0) {
type = 2;
} else if (TimeSegmentationInterp.compare("histo") == 0) {
type = 3;
} else {
std::cout << "Did not recognize temporal interpolator selection. " << std::endl
<< "Acceptable values are hanning or minmax." << std::endl;
return 1;
}

} catch (boost::program_options::error& e) {
std::cerr << "Error: " << e.what() << std::endl;
std::cout << desc << std::endl;
Expand Down Expand Up @@ -264,6 +252,25 @@ int main(int argc, char** argv)
getCompleteISMRMRDAcqData<float>(d, acqTrack, NSlice, NRep, NAvg, NEcho, NPhase, data, kx, ky,
kz, tvec);


// Deal with the number of time segments
if(!vm.count("TimeSegments")) {
switch(type) {
case 1:
L = ceil((arma::max(tvec) - arma::min(tvec))/2E-3);
break;
case 2:
L = ceil((arma::max(tvec) - arma::min(tvec))/3E-3);
break;
case 3:
L = ceil((arma::max(tvec) - arma::min(tvec))/3E-3);
break;
default:
L = 0;
}
std::cout << "Info: Setting L = " << L << " by default." << std::endl;
}

PMap = getISMRMRDCompletePhaseMap<float>(d, NSlice, NSet, NRep, NAvg, NPhase, NEcho, NSeg,
(uword)(Nx * Ny * Nz));
SMap = getISMRMRDCompleteSENSEMap<std::complex<float> >(d, sen, NSlice, (uword)(Nx * Ny * Nz));
Expand Down
51 changes: 32 additions & 19 deletions PowerGrid/PowerGridIsmrmrd.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -72,34 +72,28 @@ int main(int argc, char **argv) {
std::cout << desc << std::endl;
return 1;
}
/*
if(precisionString.compare("double") ==0) {
typedef double PGPrecision;
} else if(precisionString.compare("float") == 0) {
typedef float PGPrecision;
} else {
typedef double PGPrecision;
std::cout << "Did not recognize precision option. Defaulting to
double precision." << std::endl;
}
*/


if (FourierTrans.compare("DFT") == 0) {
FtType = 2;

} else if (FourierTrans.compare("NUFFT") == 0) {

FtType = 1;
if (TimeSegmentationInterp.compare("hanning") == 0) {
if(!vm.count("TimeSegmentationInterp")) {
type = 1;
} else if (TimeSegmentationInterp.compare("minmax") == 0) {
type = 2;
} else if (TimeSegmentationInterp.compare("histo") == 0) {
type = 3;
} else {
std::cout << "Did not recognize temporal interpolator selection. " << std::endl
<< "Acceptable values are hanning or minmax." << std::endl;
return 1;
if (TimeSegmentationInterp.compare("hanning") == 0) {
type = 1;
} else if (TimeSegmentationInterp.compare("minmax") == 0) {
type = 2;
} else if (TimeSegmentationInterp.compare("histo") == 0) {
type = 3;
} else {
std::cout << "Did not recognize temporal interpolator selection. " << std::endl
<< "Acceptable values are hanning or minmax." << std::endl;
return 1;
}
}
} else if (FourierTrans.compare("DFTGrads") == 0) {
FtType = 3;
Expand Down Expand Up @@ -220,6 +214,25 @@ int main(int argc, char **argv) {
getCompleteISMRMRDAcqData<float>(d, acqTrack, NSlice, NRep, NAvg, NEcho, NPhase, data, kx, ky,
kz, tvec);

// Deal with the number of time segments
if(!vm.count("TimeSegments")) {
switch(type) {
case 1:
L = ceil((arma::max(tvec) - arma::min(tvec))/2E-3);
break;
case 2:
L = ceil((arma::max(tvec) - arma::min(tvec))/3E-3);
break;
case 3:
L = ceil((arma::max(tvec) - arma::min(tvec))/3E-3);
break;
default:
L = 0;
}
std::cout << "Info: Setting L = " << L << " by default." << std::endl;
}


std::cout << "Number of elements in kx = " << kx.n_rows << std::endl;
std::cout << "Number of elements in ky = " << ky.n_rows << std::endl;
std::cout << "Number of elements in kz = " << kz.n_rows << std::endl;
Expand Down
51 changes: 21 additions & 30 deletions PowerGrid/PowerGridPcSenseTimeSeg.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -74,28 +74,6 @@ int main(int argc, char **argv)
std::cout << desc << std::endl;
return 1;
}
/*
if(precisionString.compare("double") ==0) {
typedef double PGPrecision;
} else if(precisionString.compare("float") == 0) {
typedef float PGPrecision;
} else {
typedef double PGPrecision;
std::cout << "Did not recognize precision option. Defaulting to
double precision." << std::endl;
}
*/
if (TimeSegmentationInterp.compare("hanning") == 0) {
type = 1;
} else if (TimeSegmentationInterp.compare("minmax") == 0) {
type = 2;
} else if (TimeSegmentationInterp.compare("histo") == 0) {
type = 3;
} else {
std::cout << "Did not recognize temporal interpolator selection. " << std::endl
<< "Acceptable values are hanning, minmax, or histo." << std::endl;
return 1;
}

}
catch (boost::program_options::error& e) {
Expand All @@ -114,9 +92,6 @@ int main(int argc, char **argv)
acqTracking* acqTrack;
processISMRMRDInput<float>(rawDataFilePath, d, hdr, FM, sen, acqTrack);

//savemat("testFM.mat", "FM", FM);
//savemat("testSen.mat", "sen", sen);

std::cout << "Number of elements in SENSE Map = " << sen.n_rows << std::endl;
std::cout << "Number of elements in Field Map = " << FM.n_rows << std::endl;

Expand Down Expand Up @@ -198,6 +173,26 @@ int main(int argc, char **argv)
getCompleteISMRMRDAcqData<float>(d, acqTrack, NSlice, NRep, NAvg, NEcho, NPhase, data, kx, ky,
kz, tvec);


// Deal with the number of time segments
if(!vm.count("TimeSegments")) {
switch(type) {
case 1:
L = ceil((arma::max(tvec) - arma::min(tvec))/2E-3);
break;
case 2:
L = ceil((arma::max(tvec) - arma::min(tvec))/3E-3);
break;
case 3:
L = ceil((arma::max(tvec) - arma::min(tvec))/3E-3);
break;
default:
L = 0;
}
std::cout << "Info: Setting L = " << L << " by default." << std::endl;
}


PMap = getISMRMRDCompletePhaseMap<float>(d, NSlice, NSet, NRep, NAvg, NPhase, NEcho, NSeg,
(uword) (Nx*Ny*Nz));
SMap = getISMRMRDCompleteSENSEMap<std::complex<float>>(d, sen, NSlice, (uword) (Nx*Ny*Nz));
Expand All @@ -212,18 +207,14 @@ int main(int argc, char **argv)

std::cout << "Number of columns in data = " << data.n_cols << std::endl;

//Gnufft<float> A(kx.n_rows, (float) 2.0, Nx, Ny, Nz, kx, ky, kz, ix,
// iy, iz);
//Gdft<float> A(kx.n_rows, Nx*Ny*Nz,kx,ky,kz,ix,iy,iz,FM,tvec);
pcSenseTimeSeg<float> S_DWI(kx, ky, kz, Nx, Ny, Nz, nc, tvec, L, type, SMap, FMap,
0-PMap);
//pcSENSE<float, Gnufft<float>> Sg(A, sen, kx.n_rows, Nx*Ny*Nz, nc);
QuadPenalty<float> R(Nx, Ny, Nz, beta, dims2penalize);

ImageTemp = reconSolve<float, pcSenseTimeSeg<float>, QuadPenalty<float>>(data, S_DWI, R, kx, ky, kz,
Nx,
Ny, Nz, tvec, NIter);
//writeISMRMRDImageData<float>(d, ImageTemp, Nx, Ny, Nz);

writeNiftiMagPhsImage<float>(filename,ImageTemp,Nx,Ny,Nz);
}
}
Expand Down

0 comments on commit 7ffdb48

Please sign in to comment.