CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   ParaView (https://www.cfd-online.com/Forums/paraview/)
-   -   [General] Get data from Calculator filter in python script (https://www.cfd-online.com/Forums/paraview/124404-get-data-calculator-filter-python-script.html)

taalf October 4, 2013 08:20
Get data from Calculator filter in python script
 
Hi all,

New on python scripting with Paraview, I am trying for a while to do a (I hope) simple thing but without any success...

I have generated in my script a Calculator filter (let's call it Calculator1).

When I display it in a Paraview spreadsheet view, it looks like below:

___|__Point_ID__|___Points___|___Result___|
_0_|__0_________|__67.15...__|__4.35e+06__|

That is: a single row with a column called Result (among other columns).

I would like to get the value of this Result (4.35e+06 here) into a python variable.

Can anyone help me with this?...

Thank you very much!

Best regards,

William

wyldckat October 6, 2013 12:56
Greetings William,

Have a look into this section: http://www.paraview.org/Wiki/ParaVie...tiple_Datasets

In other words, you can get in into another Python calculator the "Result" field, for example:
Code:
inputs[0].PointData['Result']
But if you could provide a more detailed example, it would be easier to help you.

Best regards,
Bruno

taalf October 7, 2013 03:50
1 Attachment(s)
Hi Bruno,

Thank you very much for your help !

I cannot share my case but I enclose the python script I use... (and copy paste it in this post).

I have a wing surface with pressure data. I want to split it into several 2D profiles, and to calculate the lift coefficient on each of these profiles.

The last object "Calculator2" is the one from which I would like to extract the "Result" value...

I tried Calculator2.PointData['Result'] but it returns "Array: Result"...

Here is my script:



try:
paraview.simple

except:
from paraview.simple import *

# Get active boundary (actualy a 3D wing from a mesh with pressure data)
Boundary = GetActiveSource()

# Extract surface from this wing boundary (necessary to compute normal vectors)
ExtractSurface1 = ExtractSurface( Boundary )

# Compute normal vectors of the wing surface
GenerateSurfaceNormals1 = GenerateSurfaceNormals( ExtractSurface1 )

# For several positions on the z axis...
for z in range(500,5500,1000):
# Make a cut of the wing to have a local 2D profile
Slice1 = Slice( GenerateSurfaceNormals1 )
Slice1.SliceType = "Plane"
Slice1.SliceOffsetValues = [0.0]
Slice1.SliceType.Origin = [0.0, 0.0, z]
Slice1.SliceType.Normal = [0.0, 0.0, 1.0]

# Calculate the local "pressure" vectors on the 2D profile
Calculator1 = Calculator( Slice1 )
Calculator1.Function = 'Normals*pressure'

# Make the integral of pressure vectors to have the total force exerted on the 2D profile (a single vector)
IntegrateVariables1 = IntegrateVariables( Calculator1 )

# Project the force vector to get the lift coefficient (here, with angle of attack = 0 deg)
Calculator2 = Calculator( IntegrateVariables1 )
Calculator2.Function = 'Result.(0*iHat+1*jHat)'

# Print the local lift coefficient ! :-/
print 'C_L at z =',z/1000.0,'m'
print Calculator2.PointData['Result']# Doesn't work... :-|

taalf October 7, 2013 07:56
2 Attachment(s)
Hi,

I made a test case on which to run my script:

https://www.dropbox.com/s/tacsiz95qe...g_pressure.zip

It consist of a Ensight Gold case.

I just use an Extract Block filter to extract the "wing" boundary, then select the ExtractBlock1 object and open a Python shell and run the script.

Before: Attachment 25866

After: Attachment 25867

Best regards,

William

taalf October 7, 2013 09:11
I finally use a workaround consisting in exporting data into CSV files...

http://www.paraview.org/Wiki/ParaVie...rting_CSV_Data

Adapting to my script, it looks like:

writer = CreateWriter(str(z)+".csv", Calculator2)
writer.FieldAssociation = "Points"
writer.UpdatePipeline()
del writer

Like this I have separated files for each 2D profile. Not beautiful at all, but it works...

But I hope there is a way to avoid such CSV file export and to get values directly in the python script...

Best regards,

William

wyldckat October 13, 2013 06:33
Hi William,

Many thanks for sharing the information you've managed to figure out!

The idea I had given on the other post relied on using the filter "Python Calculator", which does provide access to values on arrays, but it works differently from the main Python scripting system that ParaView uses.

It took me a few minutes to figure this one out. The instructions are partially given here: http://paraview.org/Wiki/ParaView/Py...Source_Proxies
The rest I used the knowledge I have of VTK.

The example is as follows, with comments along the way:
Code:
Calculator2=GetActiveSource()

#get access to the VTK data
realCalculator2=servermanager.Fetch(Calculator2)

#The resulting calculation is in Point data format
rcPointData=realCalculator2.GetPointData()

#Get the array we want, which in my case is named "Result2", which is the result of "Result*iHat" from the previous integrated results
Result2Array=rcPointData.GetArray('Result2')

#Now you can access the values in two ways:
#1- Directly access all values on the array
manualVector = [Result2Array.GetValue(0), Result2Array.GetValue(1), Result2Array.GetValue(2)]

#2- Directly access the vector itself, without having to do it manually
automaticVector = Result2Array.GetTuple(0)
Best regards,
Bruno

taalf October 14, 2013 08:17
1 Attachment(s)
Thank you Bruno ! ! ! ! :D

I dreamed it, you did it ! :)

It works perfectly :)

Enclosed is my script fixed with your solution :)

Thank you again ! :D

Best regards,

William

taalf October 15, 2013 12:06
1 Attachment(s)
Hi all...

I am still struggling making my script.

I figured out that I have to use the Python Calculator to do what I need.

But after creating the Python Calculator object, I cannot use the solution above on it.

I put
Code:
data=servermanager.Fetch(areaCalculator)
But then if I do
Code:
pointData=data.GetPointData()
I have an error... Thought I can display the Python Calculator in a Paraview Spreadsheet and see the Point Data...

If someone can help.............. :(

Enclosed is my script (please, use the previous case and modify the path inside the python script).

Code:
try: paraview.simple
except: from paraview.simple import *

# Import data

wing_pressure_case = EnSightReader( CaseFileName='/home/tougeron/test/wing_pressure.case' )
wing_pressure_case.CellArrays = []
wing_pressure_case.PointArrays = ['pressure']

# Extract wing

wing=ExtractBlock(wing_pressure_case)
wing.BlockIndices = [2]

# Compute normals

wingSurface=ExtractSurface(wing)
wingSurfaceWithNormals=GenerateSurfaceNormals(wingSurface)

# Clip the wing

wingClip1=Clip(wingSurfaceWithNormals,ClipType="Plane")
wingClip1.Scalars = ['POINTS', 'pressure']
wingClip1.ClipType.Origin = [0.0, 0.0, 1000.0]
wingClip1.ClipType.Normal = [0.0, 0.0, 1.0]

wingClip2=Clip(wingClip1,ClipType="Plane")
wingClip2.Scalars = ['POINTS', 'pressure']
wingClip2.ClipType.Origin = [0.0, 0.0, 1010.0]
wingClip2.ClipType.Normal = [0.0, 0.0, -1.0]

# Get the local area

areaCalculator=PythonCalculator(wingClip2)
areaCalculator.ArrayAssociation='Point Data'
areaCalculator.ArrayName='area'
areaCalculator.Expression='area(inputs[0])'

# Get data

data=servermanager.Fetch(areaCalculator)

pointData=data.GetPointData() # <------- Doesn't work ! :-( ???

wyldckat October 15, 2013 16:18
Hi William,

I won't be able to look into the code before the weekend, but there is a quick tip I can give right now... actually, I already mentioned in the previous post, so I'll rephrase it:
  1. The Python shell in ParaView (accessible from the menu "Tools"), gives full access from an outside view of the Python environment. In other words, it has it's own workspace.
  2. The "Python Calculator", as well as the "Programmable Filter", use their own Python workspaces and do not provide the full feature-set that is available on the Python Shell.
Nonetheless, if my memory doesn't fail me, the "print" function works on both workspaces and should be able to print out to the Python Shell, therefore this is a way that you can use to diagnose how things are working inside the other workspace.


Another tip is from another thread and I quote:
Quote:
Originally Posted by mpeti (Post 455435)
I believe it is worth mentioning that it is possible to put the debugger breakpoints (pdb.set_trace()) in the progFilter.Script string and this way to figure out what exactly the input of the filter is like.
I've never used this myself, but it seems promising.

Best regards,
Bruno

taalf October 16, 2013 08:10
1 Attachment(s)
Hi @wyldckat,

Thank you for your time!

I read you, trust you, but I still don't see how to do...

But the thing which is very strange to me is that my script works perfectly with a Python Calculator made on a source Sphere, but not with the same Python Calculator made on my wing surface...

Here is my script:

Code:
try: paraview.simple
except: from paraview.simple import *

test_on='wing'
#test_on='sphere'

if test_on=='wing':
   
    # Import data
   
    wing_pressure_case = EnSightReader( CaseFileName='/home/tougeron/test/wing_pressure.case' )
    wing_pressure_case.CellArrays = []
    wing_pressure_case.PointArrays = ['pressure']
   
    # Extract wing
   
    wing=ExtractBlock(wing_pressure_case)
    wing.BlockIndices = [2]
   
    # Compute normals
   
    wingSurface=ExtractSurface(wing)
    wingSurfaceWithNormals=GenerateSurfaceNormals(wingSurface)
   
    # Clip the wing
   
    wingClip1=Clip(wingSurfaceWithNormals,ClipType="Plane")
    wingClip1.Scalars = ['POINTS', 'pressure']
    wingClip1.ClipType.Origin = [0.0, 0.0, 1000.0]
    wingClip1.ClipType.Normal = [0.0, 0.0, 1.0]
   
    wingClip2=Clip(wingClip1,ClipType="Plane")
    wingClip2.Scalars = ['POINTS', 'pressure']
    wingClip2.ClipType.Origin = [0.0, 0.0, 1010.0]
    wingClip2.ClipType.Normal = [0.0, 0.0, -1.0]
   
    # Get the local area
   
    areaCalculator=PythonCalculator(wingClip2)
   

if test_on=='sphere':
   
    # Create a sphere
   
    Sphere1 = Sphere()
    Sphere1.Radius = 5.0
    Sphere1.ThetaResolution = 36
    Sphere1.PhiResolution = 36
   
    # Get the local area
   
    areaCalculator=PythonCalculator(Sphere1)
   

areaCalculator.ArrayAssociation='Point Data'
areaCalculator.ArrayName='area'
areaCalculator.Expression='area(inputs[0])'

# Get data

data=servermanager.Fetch(areaCalculator)

nbPoints=data.GetNumberOfPoints()
print 'Found',nbPoints,'points'

pointData=data.GetPointData() # <------- Doesn't work with wing! :-( ???

areaArray=pointData.GetArray('area')

area=[]
for i in range(nbPoints): area.append(areaArray.GetTuple(i)[0])

print len(area),"values got successfully:"
print area[0]
print area[1]
print area[2]
print "..."
When the
Code:
#test_on='sphere'
line is uncommented, so the script can get the local area from the areaCalculator object which is a Python Calculator of which the source is a Sphere.

Output:

Quote:
Found 1226 points
1226 values got successfully:
0.0174589491329
0.0174589529311
0.0174589475408
...
When the same line is commented, the areaCalculator is made exactly by the same way but have the wing clip as a source. And in this case I get the following error:

Quote:
Found 354 points
Traceback (most recent call last):
File "<string>", line 69, in <module>
AttributeError: GetPointData
. . . . . .

:confused:

wyldckat October 16, 2013 14:01
Hi William,

A quick question: does it work if you work with it manually on ParaView?

Because I'm guessing that you don't know about the method "UpdatePipeline()": http://paraview.org/Wiki/ParaView/Py...ing_a_Pipeline
Quote:
At this point, if you are interested in getting some information about the output of the shrink filter, you can force it to update (which will also cause the execution of the cone source). For details about VTK's demand-driven pipeline model used by ParaView, see one of the VTK books.
Code:
>>> shrinkFilter.UpdatePipeline()
>>> shrinkFilter.GetDataInformation().GetNumberOfCells()
33L
>>> shrinkFilter.GetDataInformation().GetNumberOfPoints()
128L
As you can see, the data should only be available after you've forced an update. ;)

Best regards,
Bruno

taalf October 17, 2013 02:59
Yahaaaa !

Thank you very much @Wyldckat ! :D

with a
Code:
areaCalculator.UpdatePipeline()
line before the "Fetch" command it works perfectly! :D

I hope you don't think I don't read the documentation :)

Just that it is lots of info at a time and I don't have much time to develop my script...

Merci infiniment ! :D

Best regards,

William

taalf October 17, 2013 06:18
:(

Spoke to fast...

It doesn't change anything actually..........................

:(

wyldckat October 17, 2013 14:32
Hi William,

Try on the other variables before that part.
The biggest suspect is "wing_pressure_case" since it's the reader and seems to be a bit more tricky to use, given its configuration settings.

And if "UpdatePipeline()" doesn't work, try "Update()".

Best regards,
Bruno

taalf October 18, 2013 06:27
Hi @Wyldckat,

I gave up with this method.

I put an UpdatePipe() on all objects. I tried Update() without success.

I came back to the CSV file solution.

I export the data into a temporary CSV file, then import it with the csv python library and delete it. It takes one cent of second and it is very robust.

Maybe I will later understand what I did wrong.

Thank you very much for your precious help anyway ! ;)

Best regards,

William

wyldckat October 19, 2013 14:55
Hi William,

OK, I've figured out the problem. If you had checked what "data" was, it would tell you this:
Code:
>>> data
(vtkMultiBlockDataSet)0x5b143c0
"MultiBlocks" have to be dealt differently, since each one of the blocks contains a part of the geometry+data. In this case, it's likely just one block, but still it is packed inside a "MultiBlock".

To unwrap the data, you can use "MergeBlocks":
Code:
mergedAreaCalculator=MergeBlocks(areaCalculator)

# Get data

data=servermanager.Fetch(mergedAreaCalculator)
The rest of the code should now work as intended!

If I'm not mistaken, this filter can be applied directly to the reader, so that you won't have to be bothered by it later on in the code :).

Best regards,
Bruno

taalf October 21, 2013 03:38
Hello, @Wyldckat :)

I had seen that the wing data was under multi-block format and noticed this was the only one difference between it and the source sphere case.

But I really didn't know what to do with this! I tried lots of things, and finally gave up!

Now my script is ready to use with the CSV solution. Don't know if I will use the MergeBlock solution on it.

But I am happy to see the source of the mistake ! (I tested on my test script, it works)

Thank you very much! :)

Best regards,

William

Ouardia February 15, 2018 04:53
Using mergeBlock in Paraview
 
Hi all,

I have a question about using mergeBlock in Paraview ! I tried to save the data on a .vtk file, it works for the first iteration but I need to save the last iteration's data in my .vtk file.

Any one have an idea please ? :confused::confused:

Thanks


All times are GMT -4. The time now is 06:31.