Will check if anti-aliased data exist.

#

# trace generated using paraview version 5.8.1

#

# To ensure correct image size when batch processing, please search 

# for and uncomment the line `# renderView*.ViewSize = [*,*]`

#### import the simple module from the paraview

#

# import the simple module from the paraview

from paraview.simple import *

import glob

import os

import numpy as np

import vtk as vtk

ResetSession()

print(100*'\n')

##def updateFlow():

##    extractSelection1 = FindSource('ExtractSelection1')

##    Delete(extractSelection1)

##    del extractSelection1

##

##    extractSelection1 = ExtractSelection(Input=programmableFilter2)

##    myExtract = FindSource('ExtractSelection1')

##    try:

##        plotDataOverTime1Display.CompositeDataSetIndex =\

##                            FindSource('ExtractSelection1').Selection.IDs[1::2]

##    except:

##        print("Error 1 in updateFlow() function.")

##    try:

##        plotDataOverTime

##        plotDataOverTime1Display.SeriesVisibility = ['Flux (id=' + \

##            str(myID-1) +')' for myID in FindSource('ExtractSelection1').Selection.IDs[1::2]]

##    except:

##        print("Error 2 in updateFlow() function.")

def globalMaximum(temporalStatistics_object):        

    minrange  =  servermanager.Fetch(temporalStatistics).\

                 GetPointData().GetArray('Flux_minimum').GetValueRange()

    maxrange  =  servermanager.Fetch(temporalStatistics).\

                 GetPointData().GetArray('Flux_maximum').GetValueRange()

def globalMaximum(temporalStatistics_object,  quantity):

    """Computes  global  maximum  over  all  time  frames."""

    minrange  =  (servermanager.Fetch(temporalStatistics_object).

                 GetPointData().GetArray(quantity  +  '_minimum').

                  GetValueRange())

    maxrange  =  (servermanager.Fetch(temporalStatistics_object).

                 GetPointData().GetArray(quantity  +  '_maximum').

                  GetValueRange())

    return  np.max(np.abs([minrange,maxrange]))

# paths and filenames

base = r'C:\\Users\\galr\\Desktop\\pacienti_rozpracovani\\G318_patient18_bodyPart1_visit1_manual_autom'.replace('\\','\\\\')  +  2*os.sep

ResetSession()

print(100*'\n')

print('-------------------  STARTING  A  NEW  SESSION -----------------------')

##base = (r'C:\\Users\\galr\\Desktop\\pacienti_rozpracovani\\G318_patient18_bodyPart1_visit1_manual_autom_correct_orientation'.replace('\\','\\\\')  +  2*os.sep)

base = (r'C:\\Users\\galr\\Desktop\\pacienti_rozpracovani\\G318_patient14_bodyPart1_visit1_manual_autom_failed_correct_orientation'.replace('\\','\\\\')  +  2*os.sep)

##base  =  r'C:\Users\galr\Desktop\pacienti_rozpracovani\G318_patient62_bodyPart2_visit0_autom_manual'.replace('\\','\\\\')  +  2*os.sep

##base  =  (os.sep.join(os.path.normpath(base).split(os.path.sep)  +  list(os.sep))).decoce(unicode-escape)

##base  = os.path.normpath(r'C:\\Users\\galr\\Desktop\\pacienti_rozpracovani\\G318_patient62_bodyPart2_visit0_autom_manual')  +  2*os.sep 

##base = "C:\\Users\\galr\\Desktop\\pacienti_rozpracovani\\HeartNavJC20210817_"  +  2*os.sep  

##base = "C:\\Users\\galr\\Desktop\\pacienti_rozpracovani\\HeartNavMM20210817_"  +  2*os.sep

base = "C:\\Users\\galr\\Desktop\\pacienti_rozpracovani\\HeartNavAD20210824_"  +  2*os.sep 

meas_info_filename = "measurement_info" # series numbers and VENCs

sampling_percent = 1 # sometimes probably due to memory reasons only low values (e.g. 17) are possible 

##base = "C:\\Users\\galr\\Desktop\\pacienti_rozpracovani\\HeartNavAD20210824_"  +  2*os.sep 

meas_info_filename = "measurement_info"  # Import  series numbers and VENCs.

sampling_percent = 100 # sometimes probably due to memory reasons only low values (e.g. 17) are possible 

beginning_only = 0

selectFlag = 1;

automFlag  =  1;

selectFlag = 1

automFlag  =  1

my_color_map = 'jet'

if  automFlag ==  1:

    sampling_percent = 100

if (automFlag ==  1):

    sampling_percent = 100  #

#load measurement_file 

"""Loads Series numbers, VENCs and subtype of sequence

# get active view ASI MOZNO ZAKOMENTOVAT PRO PVBATCH

renderView1 = GetActiveViewOrCreate('RenderView')

print('Loading  data  in: ' + base)

vti_folders = sorted(glob.glob(os.sep.join([base,'antialiasing','vti_[0-9]*'])))

if vti_folders == []:

    vti_folders = sorted(glob.glob(os.sep.join([base,'vti_[0-9]*'])))

    print('Possibly  velocity  aliased  data loaded.')

else:

    print('Velcity  anti-aliased data loaded.')

magn_folder  =  [entry  for entry  in  vti_folders if  entry.endswith('_0')]

if  magn_folder:

# RESAMpLE MESH to 4D DATA, create a new 'Resample With Dataset'

resampleWithDataset1 = ResampleWithDataset(SourceDataArrays=legacyVTKReader1,

    DestinationMesh=xMLImageDataReader1)

resampleWithDataset1 = ResampleWithDataset(

        SourceDataArrays=legacyVTKReader1,

        DestinationMesh=xMLImageDataReader1

        )

# show data in view ASI MOZNO ZAKOMENTOVAT PRO PVBATCH

resampleWithDataset1Display = Show(resampleWithDataset1, renderView1, 'UniformGridRepresentation')

resampleWithDataset1Display = Show(

        resampleWithDataset1,

        renderView1,

        'UniformGridRepresentation'

        )

# COMBINE all 4D components, create a new 'Append Attributes'

appendAttributes1 = AppendAttributes(Input=[xMLImageDataReader1, xMLImageDataReader2, xMLImageDataReader3, xMLImageDataReader4])

appendAttributes1 = AppendAttributes(Input=[xMLImageDataReader1,

                                            xMLImageDataReader2,

                                            xMLImageDataReader3,

                                            xMLImageDataReader4])

# show data in view ASI MOZNO ZAKOMENTOVAT PRO PVBATCH

appendAttributes1Display = Show(appendAttributes1, renderView1, 'UniformGridRepresentation')

appendAttributes1Display = Show(appendAttributes1,

                                renderView1,

                                'UniformGridRepresentation')

# COMPUTE the velocity vector field, create a new 'Calculator'

calculator1 = Calculator(Input=appendAttributes1)

# Properties modified on calculator1

calculator1.ResultArrayName = 'velocity_vectors'

calculator1.Function = 'scalars*' + components_in_calc1[0] + \

			'+scalars_input_1*' + components_in_calc1[1] + \

			'+scalars_input_2*' + components_in_calc1[2]

calculator1.Function = ('scalars*' + components_in_calc1[0]

			+  '+scalars_input_1*' + components_in_calc1[1] 

			+  '+scalars_input_2*' + components_in_calc1[2])

# show data in view ASI MOZNO ZAKOMENTOVAT PRO PVBATCH

calculator1Display = Show(calculator1, renderView1, 'UniformGridRepresentation')

calculator1Display = Show(calculator1,

                          renderView1,

                          'UniformGridRepresentation')

# ADD the vector field to previous data, create a new 'Append Attributes'

appendAttributes2 = AppendAttributes(Input=[calculator1, resampleWithDataset1])

# show data in view ASI MOZNO ZAKOMENTOVAT PRO PVBATCH

appendAttributes2Display = Show(appendAttributes2, renderView1, 'UniformGridRepresentation')

appendAttributes2Display = Show(appendAttributes2,

                                renderView1,

                                'UniformGridRepresentation')

# MASK the vector field by MESH, create a new 'Calculator'

calculator2 = Calculator(Input=appendAttributes2)

calculator2Display = Show(calculator2, renderView1, 'UniformGridRepresentation')

# LOAD centerline, create a new 'Legacy VTK Reader'

#legacyVTKReader2 = LegacyVTKReader(FileNames=[centerline_path])

legacyVTKReader2 = LegacyVTKReader(FileNames=[centerline_path for ii in range(len(vti_files[0]))])

list_of_centerline_copies =  [centerline_path]*(len(vti_files[0]))

legacyVTKReader2 = LegacyVTKReader(FileNames=list_of_centerline_copies)

# show data in view ASI MOZNO ZAKOMENTOVAT PRO PVBATCH

legacyVTKReader2Display = Show(legacyVTKReader2, renderView1, 'UnstructuredGridRepresentation')

legacyVTKReader2Display = Show(legacyVTKReader2,

                               renderView1,

                               'UnstructuredGridRepresentation')

# DIFFERENTIATE centerline for Tangents, create a new 'Programmable Filter'

programmableFilter1 = ProgrammableFilter(Input=legacyVTKReader2)

"""

# show data in view ASI MOZNO ZAKOMENTOVAT PRO PVBATCH

programmableFilter1Display = Show(programmableFilter1, renderView1, 'UnstructuredGridRepresentation')

programmableFilter1Display = Show(programmableFilter1,

                                  renderView1,

                                  'UnstructuredGridRepresentation')

# DEFINE a ball around a point to use later in integration, create a new 'Calculator'

calculator3 = Calculator(Input=calculator2)

calculator3.Function = 'masked_velocity*(floor((1+exp(-2*20*(10-mag(coords-(166*iHat+145*jHat+95*kHat)))))^(-1)+0.1))'

#the constants with value 20 and 2 determine the sharpness of a Step function

#the constant with value 10 determines the radius of the ball

calculator3.Function = ('masked_velocity*(floor((1+exp(-2*20*(10-mag(coords-'

                        +'(166*iHat+145*jHat+95*kHat)))))^(-1)+0.1))')

#the constants 20 and 2 determine the sharpness of the Step function

#the constant 10 determines the radius of the ball

#the 0.1 shift ensures step close to 0

calculator3.ResultArrayName = 'cropped_vectors'

#nacteni centerline_skeleton

centerline_skeleton = paraview.servermanager.Fetch(programmableFilter1)

points = np.array(centerline_skeleton.GetPointData().GetArray('Coords'))

radii = np.array(paraview.servermanager.Fetch(legacyVTKReader2).GetPointData().GetArray('Radius'))

radii = (np.array(paraview.servermanager.Fetch(legacyVTKReader2).

                 GetPointData().GetArray('Radius')))

diffPoints = np.array(centerline_skeleton.GetPointData().GetArray('Tangent'))

# CREATE a new 'Slice' to iterate along the centerline

# create a new 'Generate Surface Normals'

generateSurfaceNormals1 = GenerateSurfaceNormals(Input=slice1)

# show data in view ASI MOZNO ZAKOMENTOVAT PRO PVBATCH

generateSurfaceNormals1Display = Show(generateSurfaceNormals1, renderView1, 'GeometryRepresentation')

generateSurfaceNormals1Display = Show(generateSurfaceNormals1,

                                      renderView1,

                                      'GeometryRepresentation')

# COMPUTE the normal projection of velocity vectos, create a new 'Calculator'

calculator4 = Calculator(Input=generateSurfaceNormals1)

pointDatatoCellData1 = PointDatatoCellData(Input=calculator4)

pointDatatoCellData1.PointDataArraytoprocess = ['projected_vectors']

# show data in view ASI MOZNO ZAKOMENTOVAT PRO PVBATCH

pointDatatoCellData1Display = Show(pointDatatoCellData1, renderView1, 'UnstructuredGridRepresentation')

pointDatatoCellData1Display = Show(pointDatatoCellData1,

                                   renderView1,

                                   'UnstructuredGridRepresentation')

npoints = len(points)

                    np.save(os.sep.join([base,'flux_sampPerc' + str(sampling_percent),'flux' + str(tt) + '.npy']),flux)   

    if  automFlag  ==  0:

    if (automFlag  ==  0):

        # LOAD flux data, create a new 'Programmable Filter'

        programmableFilter2 = ProgrammableFilter(Input=programmableFilter1)

        # set active view

    temporalStatistics.ComputeStandardDeviation = 0

    temporalStatistics.UpdatePipeline()

    globMax  =  globalMaximum(temporalStatistics)

    globMax  =  globalMaximum(temporalStatistics,'Flux')

    # Rescale transfer function

    fluxLUT.RescaleTransferFunction(-globMax, globMax)

        elif  automFlagPython  ==  1:

            secID_array = np.array([int(_) for _ in  secID_quant[0]])

            print(secID_array)

            ##print(secID_array)

            ID_list  =  list(set(secID_array))

            print(ID_list)

            ##print(ID_list)

            ID_indices =  []

            for  ID  in  ID_list:

                ID_indices.append(np.where(secID_array  ==  ID))

                print(ID_indices[ID])

                print(type(ID_indices[ID])) ##tuple

                ##print(ID_indices[ID])

                ##print(type(ID_indices[ID])) ##tuple

                dist_array  =  np.array(dist_quant)

                ##print(type(coords_array))

                ##cc  =  coords_array[0,list(ID_indices[ID])]

## alternatively, if you want to write images, you can use SaveScreenshot(...).

print('Select points on the centerline (hit the "Select Points On" button or the "d" key) and hit "FlowInTime" button in order to plot flow rate in time.')