Postprocess PIV results (simple calibration)#

import h5py
import numpy as np

from fluidimage import get_path_image_samples

from fluidimage.calibration import Calibration2DSimple

path_base = get_path_image_samples() / 'Milestone'
path_piv = path_base / 'piv_0000a-b.h5'

requesting all changes
adding changesets
adding manifests
adding file changes
added 1121 changesets with 5003 changes to 970 files
new changesets bed748d77948:e4bedcab9258
updating to branch default
420 files updated, 0 files merged, 0 files removed, 0 files unresolved

Create the calibration object#

with h5py.File(path_piv) as file:
    shape_im = file['couple/shape_images'][...]

point0 = (10, 20)
point1 = (10, 10)

calib = Calibration2DSimple(point0, point1, distance=0.01, aspect_ratio_pixel=1.0,
                            shape_image=shape_im, point_origin=[i//2 for i in shape_im])

print(calib.yphys1pixel)

0.001

Application to synthetic data#

We can try the calibration object on synthetic data:

ixs = np.linspace(0, 100, 11)
iys = np.linspace(0, 100, 11)
xphys, yphys = calib.pix2phys(ixs, iys)
print(xphys, yphys, sep='\n')

[-1.28 -1.27 -1.26 -1.25 -1.24 -1.23 -1.22 -1.21 -1.2  -1.19 -1.18]
[-1.08 -1.07 -1.06 -1.05 -1.04 -1.03 -1.02 -1.01 -1.   -0.99 -0.98]

dxs = np.random.randint(-5, 5, ixs.size)
dys = np.random.randint(-5, 5, ixs.size)
xphys, yphys, dxphys, dyphys = calib.displ2phys(ixs, iys, dxs, dys)
print(xphys, yphys, dxphys, dyphys, sep='\n')

[-1.28 -1.27 -1.26 -1.25 -1.24 -1.23 -1.22 -1.21 -1.2  -1.19 -1.18]
[-1.08 -1.07 -1.06 -1.05 -1.04 -1.03 -1.02 -1.01 -1.   -0.99 -0.98]
[-0.002 -0.002  0.004 -0.004 -0.003 -0.005 -0.002 -0.002  0.    -0.004
 -0.001]
[-0.004  0.002  0.     0.001 -0.001  0.002 -0.001  0.    -0.001  0.004
 -0.005]

Application to real data#

with h5py.File(path_piv) as file:
    deltaxs_final = file['piv1/deltaxs_final'][:]
    deltays_final = file['piv1/deltays_final'][:]
    ixvecs_final = file['piv1/ixvecs_final'][:]
    iyvecs_final = file['piv1/iyvecs_final'][:]

xphys, yphys, dxphys, dyphys = calib.displ2phys(
    ixvecs_final, iyvecs_final, deltaxs_final, deltays_final)

print(xphys, yphys, sep='\n')

[-1.264 -1.264 -1.264 ...  1.248  1.248  1.248]
[-1.064 -1.048 -1.032 ...  1.016  1.032  1.048]

2d grid#

We need to produce a good grid in the physical space to interpolate the data:

from fluidimage.postproc.piv import get_grid_pixel_from_piv_file

xs1d, ys1d = get_grid_pixel_from_piv_file(path_piv)

print(xs1d, ys1d, sep='\n')

WARNING  Module                                                                                                    
         /home/docs/checkouts/readthedocs.org/user_builds/fluidimage/envs/latest/lib/python3.11/site-packages/fluid
         image/calcul/fft.py has not been compiled for Transonic-Pythran

WARNING  Module                                                                                                    
         /home/docs/checkouts/readthedocs.org/user_builds/fluidimage/envs/latest/lib/python3.11/site-packages/fluid
         image/calcul/subpix.py has not been compiled for Transonic-Pythran

WARNING  Module                                                                                                    
         /home/docs/checkouts/readthedocs.org/user_builds/fluidimage/envs/latest/lib/python3.11/site-packages/fluid
         image/calcul/correl.py has not been compiled for Transonic-Pythran

WARNING  Module                                                                                                    
         /home/docs/checkouts/readthedocs.org/user_builds/fluidimage/envs/latest/lib/python3.11/site-packages/fluid
         image/calcul/mean_neighbors.py has not been compiled for Transonic-Pythran

WARNING  Module                                                                                                    
         /home/docs/checkouts/readthedocs.org/user_builds/fluidimage/envs/latest/lib/python3.11/site-packages/fluid
         image/calcul/interpolate/thin_plate_spline.py has not been compiled for Transonic-Pythran

Warning: function compute_tps_matrix_numpy not pythranized.
[  16   32   48   64   80   96  112  128  144  160  176  192  208  224
256  272  288  304  320  336  352  368  384  400  416  432  448
480  496  512  528  544  560  576  592  608  624  640  656  672
704  720  736  752  768  784  800  816  832  848  864  880  896
928  944  960  976  992 1008 1024 1040 1056 1072 1088 1104 1120
1152 1168 1184 1200 1216 1232 1248 1264 1280 1296 1312 1328 1344
1376 1392 1408 1424 1440 1456 1472 1488 1504 1520 1536 1552 1568
1600 1616 1632 1648 1664 1680 1696 1712 1728 1744 1760 1776 1792
1824 1840 1856 1872 1888 1904 1920 1936 1952 1968 1984 2000 2016
2048 2064 2080 2096 2112 2128 2144 2160 2176 2192 2208 2224 2240
2272 2288 2304 2320 2336 2352 2368 2384 2400 2416 2432 2448 2464
2496 2512 2528 2544]
[  16   32   48   64   80   96  112  128  144  160  176  192  208  224
256  272  288  304  320  336  352  368  384  400  416  432  448
480  496  512  528  544  560  576  592  608  624  640  656  672
704  720  736  752  768  784  800  816  832  848  864  880  896
928  944  960  976  992 1008 1024 1040 1056 1072 1088 1104 1120
1152 1168 1184 1200 1216 1232 1248 1264 1280 1296 1312 1328 1344
1376 1392 1408 1424 1440 1456 1472 1488 1504 1520 1536 1552 1568
1600 1616 1632 1648 1664 1680 1696 1712 1728 1744 1760 1776 1792
1824 1840 1856 1872 1888 1904 1920 1936 1952 1968 1984 2000 2016
2048 2064 2080 2096 2112 2128 2144]

print(len(xs1d), len(ys1d))

159 134

xs2d, ys2d = np.meshgrid(xs1d, ys1d)

assert xs2d.shape == ys2d.shape
print(xs2d.shape)

(134, 159)

print(xs2d)

[[  16   32   48 ... 2512 2528 2544]
 [  16   32   48 ... 2512 2528 2544]
 [  16   32   48 ... 2512 2528 2544]
 ...
 [  16   32   48 ... 2512 2528 2544]
 [  16   32   48 ... 2512 2528 2544]
 [  16   32   48 ... 2512 2528 2544]]

print(ys2d)

[[  16   16   16 ...   16   16   16]
 [  32   32   32 ...   32   32   32]
 [  48   48   48 ...   48   48   48]
 ...
 [2112 2112 2112 ... 2112 2112 2112]
 [2128 2128 2128 ... 2128 2128 2128]
 [2144 2144 2144 ... 2144 2144 2144]]

xs2dphys, ys2dphys = calib.pix2phys(xs2d, ys2d)

print(xs2dphys)

[[-1.264 -1.248 -1.232 ...  1.232  1.248  1.264]
 [-1.264 -1.248 -1.232 ...  1.232  1.248  1.264]
 [-1.264 -1.248 -1.232 ...  1.232  1.248  1.264]
 ...
 [-1.264 -1.248 -1.232 ...  1.232  1.248  1.264]
 [-1.264 -1.248 -1.232 ...  1.232  1.248  1.264]
 [-1.264 -1.248 -1.232 ...  1.232  1.248  1.264]]

print(ys2dphys)

[[-1.064 -1.064 -1.064 ... -1.064 -1.064 -1.064]
 [-1.048 -1.048 -1.048 ... -1.048 -1.048 -1.048]
 [-1.032 -1.032 -1.032 ... -1.032 -1.032 -1.032]
 ...
 [ 1.032  1.032  1.032 ...  1.032  1.032  1.032]
 [ 1.048  1.048  1.048 ...  1.048  1.048  1.048]
 [ 1.064  1.064  1.064 ...  1.064  1.064  1.064]]

We define a function to interpolate the data on the grid

from scipy.interpolate import griddata

def gridd(delta):
        delta_grid = griddata((yphys, xphys), delta, (ys2dphys, xs2dphys), method='cubic')
        return delta_grid.astype('float32')

dxphys_grid = gridd(dxphys)
dyphys_grid = gridd(dyphys)

print(dxphys_grid)

[[ 1.2309151e-04 -1.1506143e-05  1.4211003e-04 ...  1.5858485e-03
   1.6704550e-03            nan]
 [ 2.4544270e-04  3.8240830e-04  9.3949668e-04 ...  1.5839678e-03
   1.4983630e-03            nan]
 [ 2.8449911e-04  4.6133468e-04  1.1535143e-03 ...  1.4397623e-03
   1.5227608e-03            nan]
 ...
 [ 2.5880121e-04 -5.5452627e-03 -5.7651512e-03 ... -1.2556877e-04
  -6.2446838e-04            nan]
 [ 2.9158618e-04 -5.5306000e-03 -6.1151003e-03 ...  1.8470040e-05
  -7.5900904e-04            nan]
 [           nan            nan            nan ...            nan
             nan            nan]]

We can finally produce PIVOnGrid objects…