Efficient FFTShift-implementation for vectors/matrices
Dear all, I tried to implement the function fftshift from numpy (i.e. swap the half-spaces of all axis) for row vectors in a matrix by using the following code void fft_shift(Mat &fft_matrix) { PetscScalar *mat_ptr; MatDenseGetArray (fft_matrix, &mat_ptr); PetscInt r_0, r_1; MatGetOwnershipRange(fft_matrix, &r_0, &r_1); PetscInt local_row_num = r_1 - r_0; arma::cx_mat temp_mat(local_row_num, Ntime, arma::fill::zeros); for(int i = 0; i < Ntime; ++i) { const PetscInt row_shift = i * local_row_num; for(int j = 0; j < local_row_num; ++j) { const PetscInt cur_pos = j + row_shift; if(i < (int)(Ntime / 2)) temp_mat(j, i + int(Ntime / 2)) = *(mat_ptr + cur_pos); else temp_mat(j, i - int(Ntime / 2)) = *(mat_ptr + cur_pos); } } for(int i = 0; i < Ntime; ++i) { const PetscInt row_shift = i * local_row_num; for(int j = 0; j < local_row_num; ++j) { const PetscInt cur_pos = j + row_shift; *(mat_ptr + cur_pos) = temp_mat(j, i); } } MatDenseRestoreArray (fft_matrix, &mat_ptr); } but I do not like the approach of having a second matrix as temporary storage space. Are there more efficient approaches possible using PETSc-functions? Thanks! Regards, Roland Richter
I found a variety of things on the web, below. I don't understand this but for the even case it seems one simply modifies the input matrix before the FFT http://www.fftw.org/faq/section3.html#centerorigin https://stackoverflow.com/questions/5915125/fftshift-ifftshift-c-c-source-co... <https://stackoverflow.com/questions/5915125/fftshift-ifftshift-c-c-source-code> https://www.dsprelated.com/showthread/comp.dsp/20790-1.php <https://www.dsprelated.com/showthread/comp.dsp/20790-1.php>
On May 18, 2021, at 9:48 AM, Roland Richter <[email protected]> wrote:
Dear all,
I tried to implement the function fftshift from numpy (i.e. swap the half-spaces of all axis) for row vectors in a matrix by using the following code
void fft_shift(Mat &fft_matrix) { PetscScalar *mat_ptr; MatDenseGetArray (fft_matrix, &mat_ptr); PetscInt r_0, r_1; MatGetOwnershipRange(fft_matrix, &r_0, &r_1); PetscInt local_row_num = r_1 - r_0; arma::cx_mat temp_mat(local_row_num, Ntime, arma::fill::zeros); for(int i = 0; i < Ntime; ++i) { const PetscInt row_shift = i * local_row_num; for(int j = 0; j < local_row_num; ++j) { const PetscInt cur_pos = j + row_shift; if(i < (int)(Ntime / 2)) temp_mat(j, i + int(Ntime / 2)) = *(mat_ptr + cur_pos); else temp_mat(j, i - int(Ntime / 2)) = *(mat_ptr + cur_pos); } } for(int i = 0; i < Ntime; ++i) { const PetscInt row_shift = i * local_row_num; for(int j = 0; j < local_row_num; ++j) { const PetscInt cur_pos = j + row_shift; *(mat_ptr + cur_pos) = temp_mat(j, i); } } MatDenseRestoreArray (fft_matrix, &mat_ptr); }
but I do not like the approach of having a second matrix as temporary storage space. Are there more efficient approaches possible using PETSc-functions?
Thanks!
Regards,
Roland Richter
participants (2)
-
Barry Smith -
Roland Richter