[Mehrez]

Hello,
I'm running a transient simulation of flow around cylinder and storing the pressure and velocity fields in files corresponding to different times.
I apply then on the pressure the ParaView filter "FFT Of Selection Over Time".
I'm not sure what is this filter doing exactly !? I didn't find deep documentation on the Fourier transform with ParaView. Do you have an idea on how does this filter choose the point on which it applies the FFT ?
Any docs are welcome.

Thanks

Mehrez

[babakflame]

Greetings Mehrez

I have been through this question as well. At first I thought maybe this filter plots energy spectrum along a specified line. I tried to search some links to find out what does this filter do exactly:

http://www.itk.org/Wiki/ParaView/Use...ters#Table_FFT
http://comments.gmane.org/gmane.comp...iew.user/22835

The above links are a bit vague. However some interesting points are mentioned. I would be very pleased if somebody describe this filter with an example. Since I have some doubts on the way of using this filter. For instance should we use it on a slice or line or why it does not accept velocity components ( I mean a separated component not a vector array) ?


Best,
Bobi

[Mehrez]

Originally Posted by babakflame

Greetings Mehrez

I have been through this question as well. At first I thought maybe this filter plots energy spectrum along a specified line. I tried to search some links to find out what does this filter do exactly:

http://www.itk.org/Wiki/ParaView/Use...ters#Table_FFT
http://comments.gmane.org/gmane.comp...iew.user/22835

The above links are a bit vague. However some interesting points are mentioned. I would be very pleased if somebody describe this filter with an example. Since I have some doubts on the way of using this filter. For instance should we use it on a slice or line or why it does not accept velocity components ( I mean a separated component not a vector array) ?


Best,
Bobi

Hi dear Bobi,
I have tried to undertand what is really done in the FFT filter but didn't get any convincing answer.
I have written my own python script that reads the OpenFoam files and plots the FFT.
Please tell me if you are interested, I can share it.
Best regards,
Mehrez

[babakflame]

Greetings Mehrez

Thanks for your reply. In deed, I am interested in deriving Energy Spectrum (1D or 3D) from LES simulations in OpenFOAM. I hope that your python script do that. It would be my very pleasure if you share it with me.


Best,
Bobi

[Mehrez]

Here is the python script:



# -*- coding:Latin-1 -*-


from math import *
import numpy
import re
import pylab
import math
import os
import glob


# Script to read files of a given variable from an transient OpenFoam case
# It performs an fft on the selected variable from the different time fields (variable may be "p" or "U"...)
# WARNING ! time steps must be equidistant
# MUST BE GIVEN : path and a reference point through which the analysis will be done


path='/home/agnaou/Bureau/Square/porosity75%/Re27000/' # path of the OpenFoam case
ref_point=5000 # an integer comprised between 0 and the number of mesh cells - 1
variable="p" # variable on which the FFT will be performed


################################################## ################################################## ########################################
# Read the pressure fields
################################################## ################################################## ########################################


folders=[x[0] for x in os.walk(path)] # search the different folders contained in the path directory

files_number=0 # number of variable fields to be treated or number of time steps
min_t=1e+50 # minimum time of the simulation
max_t=0 # maximum time of the simulation
for i in range (len(folders)):
os.chdir(folders[i])
for file in glob.glob(variable):
files_number=files_number+1
t=float(folders[i][len(folders[0]):len(folders[i])])
if (t<min_t):
min_t=t
if (t>max_t):
max_t=t
time_step=(float(max_t-min_t))/(float(files_number-1)) # time step of the simulation

VARIABLE=numpy.empty((files_number),dtype='object' ) # fields of the selected variable stored and ordered as a function of time

for i in range (len(folders)):
os.chdir(folders[i])
for file in glob.glob(variable):

TakeValue=open(folders[i]+"/"+variable)
Lines=TakeValue.readlines()
mesh=int(Lines[20])
t=float(folders[i][len(folders[0]):len(folders[i])])
VARIABLE[int((t-min_t)/time_step)]=numpy.empty((mesh),dtype='object')

for j in range (mesh):
VARIABLE[int((t-min_t)/time_step)][j]=float(Lines[j+22])

point=numpy.empty((files_number),dtype='float64')
for i in range (files_number):
point[i]=VARIABLE[i][ref_point]


################################################## ################################################## ########################################
# perform the ffts
################################################## ################################################## ########################################


fft=numpy.fft.fft(point)
frequency= numpy.fft.fftfreq(files_number, d=time_step)

fftr=abs(fft.real)
ffti=abs(fft.imag)
fftb=numpy.sqrt(fft.imag**2+fft.real**2)
frequency= numpy.fft.fftfreq(files_number, d=time_step)


################################################## ################################################## ########################################
# Plot
################################################## ################################################## ########################################


pylab.subplot(221)
pylab.title("Original Data")
pylab.grid()
pylab.plot(numpy.arange(files_number)*time_step,po int,'r-',alpha=1)
pylab.xlabel("Time")
pylab.ylabel(variable)

pylab.subplot(222)
pylab.title("Real FFT")
pylab.xlabel("Frequency")
pylab.ylabel("Power")
pylab.grid()
#pylab.plot(frequency,fftr,'b-',alpha=1)
pylab.scatter(frequency,fftr, c='b')

pylab.subplot(223)
pylab.title("Imaginary FFT")
pylab.xlabel("Frequency")
pylab.ylabel("Power")
pylab.grid()
pylab.plot(frequency,ffti,'g-',alpha=1)

pylab.subplot(224)
pylab.title("Real+Imaginary FFT")
pylab.xlabel("Frequency")
pylab.ylabel("Power")
pylab.grid()
pylab.plot(frequency,fftb,'k-',alpha=1)

pylab.show()

[babakflame]

Greetings Mehrez

Many thanks for sharing your work. I will try it and may contact you in case of questions.

Best,
Bobi

[babakflame]

Greetings Mehrez

I have activated Programmable filter on my paraview 3.12.0. Then I copied your script in the empty script space. I modified the path and cell number in your script. However, after applying your Filter nothing happens. Actually, I am interested in energy spectrum along the axis of my case. Is it possible to hint me to efficiently use your Filter?

I can send you my case so that probably it would be easier to get help.


Best,
Bobi

[Mehrez]

Hi dear Bobi,
You can use the script directly on a python shell without using Paraview.
You should have installed the loaded modules (matplotlib...) if it is not done yet.
You are right, you have to change the path, put a cell number and put a variable (replace 'p' in the original script by 'k' or 'eps' ...)
Please send me your case to take a look.
Best regards,
Mehrez.

Originally Posted by babakflame

Greetings Mehrez

I have activated Programmable filter on my paraview 3.12.0. Then I copied your script in the empty script space. I modified the path and cell number in your script. However, after applying your Filter nothing happens. Actually, I am interested in energy spectrum along the axis of my case. Is it possible to hint me to efficiently use your Filter?

I can send you my case so that probably it would be easier to get help.


Best,
Bobi

[babakflame]

Greetings Mehrez

I have installed matplotlib via the instructions in the following link:
http://matplotlib.org/users/installing.html

However, I have still my problem with using it. Is it possible to give step by step hints to use matplotlib to draw power spectrum with your code?

PS: I sent you the case with file2send, For instance I need Power Spectrum for a point in the shear layer.

PS: How do you determine cell number in your script?

Best,
Bobi

[Mehrez]

Hi dear Bobi,
I'm sorry to respond so late. The link you sent me by email didn't work (may be because I opened it late :-/).
I you are still interested, please resend me your files.
Hope that my python script helped you.

Best regards,

Mehrez

Originally Posted by babakflame

Greetings Mehrez

I have installed matplotlib via the instructions in the following link:
http://matplotlib.org/users/installing.html

However, I have still my problem with using it. Is it possible to give step by step hints to use matplotlib to draw power spectrum with your code?

PS: I sent you the case with file2send, For instance I need Power Spectrum for a point in the shear layer.

PS: How do you determine cell number in your script?

Best,
Bobi

[Kossivi]

Originally Posted by Mehrez

Here is the python script:



# -*- coding:Latin-1 -*-


from math import *
import numpy
import re
import pylab
import math
import os
import glob


# Script to read files of a given variable from an transient OpenFoam case
# It performs an fft on the selected variable from the different time fields (variable may be "p" or "U"...)
# WARNING ! time steps must be equidistant
# MUST BE GIVEN : path and a reference point through which the analysis will be done


path='/home/agnaou/Bureau/Square/porosity75%/Re27000/' # path of the OpenFoam case
ref_point=5000 # an integer comprised between 0 and the number of mesh cells - 1
variable="p" # variable on which the FFT will be performed


################################################## ################################################## ########################################
# Read the pressure fields
################################################## ################################################## ########################################


folders=[x[0] for x in os.walk(path)] # search the different folders contained in the path directory

files_number=0 # number of variable fields to be treated or number of time steps
min_t=1e+50 # minimum time of the simulation
max_t=0 # maximum time of the simulation
for i in range (len(folders)):
os.chdir(folders[i])
for file in glob.glob(variable):
files_number=files_number+1
t=float(folders[i][len(folders[0]):len(folders[i])])
if (t<min_t):
min_t=t
if (t>max_t):
max_t=t
time_step=(float(max_t-min_t))/(float(files_number-1)) # time step of the simulation

VARIABLE=numpy.empty((files_number),dtype='object' ) # fields of the selected variable stored and ordered as a function of time

for i in range (len(folders)):
os.chdir(folders[i])
for file in glob.glob(variable):

TakeValue=open(folders[i]+"/"+variable)
Lines=TakeValue.readlines()
mesh=int(Lines[20])
t=float(folders[i][len(folders[0]):len(folders[i])])
VARIABLE[int((t-min_t)/time_step)]=numpy.empty((mesh),dtype='object')

for j in range (mesh):
VARIABLE[int((t-min_t)/time_step)][j]=float(Lines[j+22])

point=numpy.empty((files_number),dtype='float64')
for i in range (files_number):
point[i]=VARIABLE[i][ref_point]


################################################## ################################################## ########################################
# perform the ffts
################################################## ################################################## ########################################


fft=numpy.fft.fft(point)
frequency= numpy.fft.fftfreq(files_number, d=time_step)

fftr=abs(fft.real)
ffti=abs(fft.imag)
fftb=numpy.sqrt(fft.imag**2+fft.real**2)
frequency= numpy.fft.fftfreq(files_number, d=time_step)


################################################## ################################################## ########################################
# Plot
################################################## ################################################## ########################################


pylab.subplot(221)
pylab.title("Original Data")
pylab.grid()
pylab.plot(numpy.arange(files_number)*time_step,po int,'r-',alpha=1)
pylab.xlabel("Time")
pylab.ylabel(variable)

pylab.subplot(222)
pylab.title("Real FFT")
pylab.xlabel("Frequency")
pylab.ylabel("Power")
pylab.grid()
#pylab.plot(frequency,fftr,'b-',alpha=1)
pylab.scatter(frequency,fftr, c='b')

pylab.subplot(223)
pylab.title("Imaginary FFT")
pylab.xlabel("Frequency")
pylab.ylabel("Power")
pylab.grid()
pylab.plot(frequency,ffti,'g-',alpha=1)

pylab.subplot(224)
pylab.title("Real+Imaginary FFT")
pylab.xlabel("Frequency")
pylab.ylabel("Power")
pylab.grid()
pylab.plot(frequency,fftb,'k-',alpha=1)

pylab.show()

Hello Mehrez,

I'm really interested in this code. I copy it and try run it but nothing happens.
Please can you send it to me through email and an example file on which to run which work for so that I figure out how it works (I'm newbe in python...).

Thanks.

Kossivi.

[Mehrez]

Hi,

Can you please show us the error message that you got ?
Thanks
Mehrez

Originally Posted by Kossivi

Hello Mehrez,

I'm really interested in this code. I copy it and try run it but nothing happens.
Please can you send it to me through email and an example file on which to run which work for so that I figure out how it works (I'm newbe in python...).

Thanks.

Kossivi.

[sadsid]

Originally Posted by Mehrez

Hi,

Can you please show us the error message that you got ?
Thanks
Mehrez

Hello.

I have few confusions. For example, my OF output folders are 0, 30, 60, and 90. Each folder has a pressure file.

Now, this code will collect pressure from each of these folders? If we want to prepare the file manually instead of this code and apply FFT directly then in this case I want to be sure about the information that your code collects from the output folders.