X-Git-Url: http://git.kpe.io/?p=ctsim.git;a=blobdiff_plain;f=libctsim%2Fimagefile.cpp;h=b5198bdca173c0569f45749c49fe69cfed049eab;hp=4a49465065257ba2f02dd5c4df6c20e2b1cd4ebc;hb=69486c2e62f1919827fc6b0c5b65b8cd6fc5384a;hpb=ab78de6e33fff94c37847621458db3ad430a8bb3 diff --git a/libctsim/imagefile.cpp b/libctsim/imagefile.cpp index 4a49465..b5198bd 100644 --- a/libctsim/imagefile.cpp +++ b/libctsim/imagefile.cpp @@ -9,7 +9,7 @@ ** This is part of the CTSim program ** Copyright (c) 1983-2001 Kevin Rosenberg ** -** $Id: imagefile.cpp,v 1.46 2002/06/27 03:19:23 kevin Exp $ +** $Id$ ** ** This program is free software; you can redistribute it and/or modify ** it under the terms of the GNU General Public License (version 2) as @@ -29,6 +29,7 @@ #ifdef HAVE_CTN_DICOM #include "ctndicom.h" #endif +#include "interpolator.h" const double ImageFile::s_dRedGrayscaleFactor = 0.299; const double ImageFile::s_dGreenGrayscaleFactor = 0.587; @@ -50,27 +51,29 @@ const int ImageFile::EXPORT_FORMAT_RAW = 6; const char* ImageFile::s_aszExportFormatName[] = { - {"text"}, - {"pgm"}, - {"pgmascii"}, + "text", + "pgm", + "pgmascii", #ifdef HAVE_PNG - {"png"}, - {"png16"}, + "png", + "png16", #endif #ifdef HAVE_CTN_DICOM - {"dicom"}, + "dicom", #endif }; const char* ImageFile::s_aszExportFormatTitle[] = { - {"Text"}, - {"PGM"}, - {"PGM ASCII"}, - {"PNG"}, - {"PNG 16-bit"}, + "Text", + "PGM", + "PGM ASCII", +#ifdef HAVE_PNG + "PNG", + "PNG 16-bit", +#endif #ifdef HAVE_CTN_DICOM - {"Dicom"}, + "Dicom", #endif }; const int ImageFile::s_iExportFormatCount = sizeof(s_aszExportFormatName) / sizeof(const char*); @@ -88,21 +91,23 @@ const int ImageFile::IMPORT_FORMAT_DICOM = 2; const char* ImageFile::s_aszImportFormatName[] = { - {"ppm"}, + "ppm", #ifdef HAVE_PNG - {"png"}, + "png", #endif #ifdef HAVE_CTN_DICOM - {"dicom"}, + "dicom", #endif }; const char* ImageFile::s_aszImportFormatTitle[] = { - {"PPM"}, - {"PNG"}, + "PPM", +#ifdef HAVE_PNG + "PNG", +#endif #ifdef HAVE_CTN_DICOM - {"Dicom"}, + "Dicom", #endif }; const int ImageFile::s_iImportFormatCount = sizeof(s_aszImportFormatName) / sizeof(const char*); @@ -720,7 +725,7 @@ ImageFile::fft (ImageFile& result) const return false; } - fftw_complex* in = new fftw_complex [m_nx * m_ny]; + fftw_complex* in = static_cast (fftw_malloc (sizeof(fftw_complex) * m_nx * m_ny)); ImageFileArrayConst vReal = getArray(); ImageFileArrayConst vImag = getImaginaryArray(); @@ -729,17 +734,17 @@ ImageFile::fft (ImageFile& result) const unsigned int iArray = 0; for (ix = 0; ix < m_nx; ix++) { for (iy = 0; iy < m_ny; iy++) { - in[iArray].re = vReal[ix][iy]; + in[iArray][0] = vReal[ix][iy]; if (isComplex()) - in[iArray].im = vImag[ix][iy]; + in[iArray][1] = vImag[ix][iy]; else - in[iArray].im = 0; + in[iArray][1] = 0; iArray++; } } - fftwnd_plan plan = fftw2d_create_plan (m_nx, m_ny, FFTW_FORWARD, FFTW_IN_PLACE | FFTW_ESTIMATE | FFTW_USE_WISDOM); - fftwnd_one (plan, in, NULL); + fftw_plan plan = fftw_plan_dft_2d (m_nx, m_ny, in, in, FFTW_FORWARD, FFTW_ESTIMATE); + fftw_execute (plan); ImageFileArray vRealResult = result.getArray(); ImageFileArray vImagResult = result.getImaginaryArray(); @@ -747,13 +752,13 @@ ImageFile::fft (ImageFile& result) const unsigned int iScale = m_nx * m_ny; for (ix = 0; ix < m_nx; ix++) { for (iy = 0; iy < m_ny; iy++) { - vRealResult[ix][iy] = in[iArray].re / iScale; - vImagResult[ix][iy] = in[iArray].im / iScale; + vRealResult[ix][iy] = in[iArray][0] / iScale; + vImagResult[ix][iy] = in[iArray][1] / iScale; iArray++; } } - delete in; - fftwnd_destroy_plan (plan); + fftw_free(in); + fftw_destroy_plan (plan); Fourier::shuffleFourierToNaturalOrder (result); @@ -791,32 +796,31 @@ ImageFile::ifft (ImageFile& result) const Fourier::shuffleNaturalToFourierOrder (result); - fftw_complex* in = new fftw_complex [m_nx * m_ny]; + fftw_complex* in = static_cast(fftw_malloc(sizeof(fftw_complex) * m_nx * m_ny)); unsigned int iArray = 0; for (ix = 0; ix < m_nx; ix++) { for (iy = 0; iy < m_ny; iy++) { - in[iArray].re = vRealResult[ix][iy]; - in[iArray].im = vImagResult[ix][iy]; + in[iArray][0] = vRealResult[ix][iy]; + in[iArray][1] = vImagResult[ix][iy]; iArray++; } } - fftwnd_plan plan = fftw2d_create_plan (m_nx, m_ny, FFTW_BACKWARD, FFTW_IN_PLACE | FFTW_ESTIMATE | FFTW_USE_WISDOM); + fftw_plan plan = fftw_plan_dft_2d (m_nx, m_ny, in, in, FFTW_BACKWARD, FFTW_ESTIMATE); - fftwnd_one (plan, in, NULL); + fftw_execute (plan); iArray = 0; for (ix = 0; ix < m_nx; ix++) { for (iy = 0; iy < m_ny; iy++) { - vRealResult[ix][iy] = in[iArray].re; - vImagResult[ix][iy] = in[iArray].im; + vRealResult[ix][iy] = in[iArray][0]; + vImagResult[ix][iy] = in[iArray][1]; iArray++; } } - fftwnd_destroy_plan (plan); - - delete in; + fftw_destroy_plan (plan); + fftw_free(in); return true; } @@ -837,24 +841,24 @@ ImageFile::fftRows (ImageFile& result) const ImageFileArrayConst vReal = getArray(); ImageFileArrayConst vImag = getImaginaryArray(); - fftw_plan plan = fftw_create_plan (m_nx, FFTW_FORWARD, FFTW_IN_PLACE | FFTW_ESTIMATE | FFTW_USE_WISDOM); + fftw_complex* in = static_cast(fftw_malloc(sizeof(fftw_complex) * m_nx)); + fftw_plan plan = fftw_plan_dft_1d (m_nx, in, in, FFTW_FORWARD, FFTW_ESTIMATE); - fftw_complex* in = new fftw_complex [m_nx]; std::complex* pcRow = new std::complex [m_nx]; for (unsigned int iy = 0; iy < m_ny; iy++) { unsigned int ix; for (ix = 0; ix < m_nx; ix++) { - in[ix].re = vReal[ix][iy]; + in[ix][0] = vReal[ix][iy]; if (isComplex()) - in[ix].im = vImag[ix][iy]; + in[ix][1] = vImag[ix][iy]; else - in[ix].im = 0; + in[ix][1] = 0; } - fftw_one (plan, in, NULL); + fftw_execute (plan); for (ix = 0; ix < m_nx; ix++) - pcRow[ix] = std::complex(in[ix].re, in[ix].im); + pcRow[ix] = std::complex(in[ix][0], in[ix][1]); Fourier::shuffleFourierToNaturalOrder (pcRow, m_nx); for (ix = 0; ix < m_nx; ix++) { @@ -865,7 +869,7 @@ ImageFile::fftRows (ImageFile& result) const delete [] pcRow; fftw_destroy_plan (plan); - delete in; + fftw_free(in); return true; } @@ -883,12 +887,11 @@ ImageFile::ifftRows (ImageFile& result) const return false; } - fftw_complex* in = new fftw_complex [m_nx]; - ImageFileArrayConst vReal = getArray(); ImageFileArrayConst vImag = getImaginaryArray(); - - fftw_plan plan = fftw_create_plan (m_nx, FFTW_BACKWARD, FFTW_IN_PLACE | FFTW_ESTIMATE | FFTW_USE_WISDOM); + + fftw_complex* in = static_cast(fftw_malloc(sizeof(fftw_complex) * m_nx)); + fftw_plan plan = fftw_plan_dft_1d (m_nx, in, in, FFTW_BACKWARD, FFTW_ESTIMATE); std::complex* pcRow = new std::complex [m_nx]; unsigned int ix, iy; @@ -904,21 +907,21 @@ ImageFile::ifftRows (ImageFile& result) const Fourier::shuffleNaturalToFourierOrder (pcRow, m_nx); for (ix = 0; ix < m_nx; ix++) { - in[ix].re = pcRow[ix].real(); - in[ix].im = pcRow[ix].imag(); + in[ix][0] = pcRow[ix].real(); + in[ix][1] = pcRow[ix].imag(); } - fftw_one (plan, in, NULL); + fftw_execute (plan); for (ix = 0; ix < m_nx; ix++) { - vReal[ix][iy] = in[ix].re; - vImag[ix][iy] = in[ix].im; + vReal[ix][iy] = in[ix][0]; + vImag[ix][iy] = in[ix][1]; } } delete [] pcRow; fftw_destroy_plan (plan); - delete in; + fftw_free(in); return true; } @@ -939,24 +942,24 @@ ImageFile::fftCols (ImageFile& result) const ImageFileArrayConst vReal = getArray(); ImageFileArrayConst vImag = getImaginaryArray(); - fftw_plan plan = fftw_create_plan (m_ny, FFTW_FORWARD, FFTW_IN_PLACE | FFTW_ESTIMATE | FFTW_USE_WISDOM); + fftw_complex* in = static_cast(fftw_malloc(sizeof(fftw_complex) * m_ny)); + fftw_plan plan = fftw_plan_dft_1d (m_ny, in, in, FFTW_FORWARD, FFTW_ESTIMATE); std::complex* pcCol = new std::complex [m_ny]; - fftw_complex* in = new fftw_complex [m_ny]; for (unsigned int ix = 0; ix < m_nx; ix++) { unsigned int iy; for (iy = 0; iy < m_ny; iy++) { - in[iy].re = vReal[ix][iy]; + in[iy][0] = vReal[ix][iy]; if (isComplex()) - in[iy].im = vImag[ix][iy]; + in[iy][1] = vImag[ix][iy]; else - in[iy].im = 0; + in[iy][1] = 0; } - fftw_one (plan, in, NULL); + fftw_execute (plan); for (iy = 0; iy < m_ny; iy++) - pcCol[iy] = std::complex(in[iy].re, in[iy].im); + pcCol[iy] = std::complex(in[iy][0], in[iy][1]); Fourier::shuffleFourierToNaturalOrder (pcCol, m_ny); for (iy = 0; iy < m_ny; iy++) { @@ -967,7 +970,7 @@ ImageFile::fftCols (ImageFile& result) const delete [] pcCol; fftw_destroy_plan (plan); - delete in; + fftw_free(in); return true; } @@ -985,12 +988,11 @@ ImageFile::ifftCols (ImageFile& result) const return false; } - fftw_complex* in = new fftw_complex [m_ny]; - ImageFileArrayConst vReal = getArray(); ImageFileArrayConst vImag = getImaginaryArray(); - fftw_plan plan = fftw_create_plan (m_ny, FFTW_BACKWARD, FFTW_IN_PLACE | FFTW_ESTIMATE | FFTW_USE_WISDOM); + fftw_complex* in = static_cast(fftw_malloc(sizeof(fftw_complex) * m_ny)); + fftw_plan plan = fftw_plan_dft_1d (m_ny, in, in, FFTW_BACKWARD, FFTW_ESTIMATE); std::complex* pcCol = new std::complex [m_ny]; unsigned int ix, iy; @@ -1006,21 +1008,21 @@ ImageFile::ifftCols (ImageFile& result) const Fourier::shuffleNaturalToFourierOrder (pcCol, m_ny); for (iy = 0; iy < m_ny; iy++) { - in[iy].re = pcCol[iy].real(); - in[iy].im = pcCol[iy].imag(); + in[iy][0] = pcCol[iy].real(); + in[iy][1] = pcCol[iy].imag(); } - fftw_one (plan, in, NULL); + fftw_execute (plan); for (iy = 0; iy < m_ny; iy++) { - vReal[ix][iy] = in[iy].re; - vImag[ix][iy] = in[iy].im; + vReal[ix][iy] = in[iy][0]; + vImag[ix][iy] = in[iy][1]; } } delete [] pcCol; fftw_destroy_plan (plan); - delete in; + fftw_free(in); return true; }