Commit 9cebb2d3 authored by Florian Ziemen's avatar Florian Ziemen
Browse files

spheres with clouds

parent 73fa62d2
......@@ -63,9 +63,6 @@ options=parse_args()
# ----------------------------------------------------------------
# Create a new 'Light'
light1 = CreateLight()
light1.Intensity = 1.0
light1.Position = [50.0, -87.0, -28.6]
# Create a new 'Light'
light2 = CreateLight()
......@@ -102,7 +99,7 @@ renderView1.EyeAngle = 1.0
renderView1.Background = [1.0, 1.0, 1.0]
renderView1.BackEnd = 'OSPRay raycaster'
renderView1.EnvironmentNorth = [0.0, 0.0, 0.0]
renderView1.AdditionalLights = [light1, light2, light3]
renderView1.AdditionalLights = [light2, light3]
renderView1.OSPRayMaterialLibrary = materialLibrary1
SetActiveView(None)
......@@ -159,6 +156,8 @@ else:
cloud_psl.CellArrayStatus = options['cloud_vars']
cloud_psl.SetProjection = 'Spherical Projection'
cloud_psl.LayerThickness = 50
print (dir(cloud_psl))
cloud_psl.SkipGrid=True
if options.get("psl", False):
print ("Getting psl", file=sys.stderr)
......@@ -193,20 +192,8 @@ if not cloud_expr:
alpha.Function = ' (3/2 * ({cloud_expr}) *100)/ (3/2 * ({cloud_expr}) *100+7)'.format(cloud_expr=cloud_expr)
# create a new 'CDIReader'
print ("Getting light", file=sys.stderr)
rsdt_file = options["light"]
rsdt = CDIReader(registrationName='rsdt', FileNames=rsdt_file.split(","))
rsdt.Dimensions = '(clon, clat, toa)'
rsdt.CellArrayStatus = ['rsdt']
rsdt.SetProjection = 'Spherical Projection'
rsdt.LayerThickness = 50
rsdt.MaskingValueVar = 'rsdt'
print ("Done getting light", file=sys.stderr)
# create a new 'NetCDF Time Annotation'
netCDFTimeAnnotation1 = NetCDFTimeAnnotation(registrationName='NetCDFTimeAnnotation1', Input=rsdt)
netCDFTimeAnnotation1 = NetCDFTimeAnnotation(registrationName='NetCDFTimeAnnotation1', Input=cloud_psl)
netCDFTimeAnnotation1.Expression = '"%02i-%02i-%02i %02i:%02i" % (Date[0], Date[1], Date[2], Date[3], Date[4])'
if not options.get("no_time_string", False):
......@@ -218,81 +205,8 @@ if not options.get("no_time_string", False):
print ("Adjusting light", file=sys.stderr)
lightadjust = ProgrammableFilter(registrationName='light-adjust', Input=rsdt)
lightadjust.OutputDataSetType = 'vtkTable'
lightadjust.Script = """
import sys
sys.path.append ("/usr/local/Caskroom/miniconda/base/lib/python3.7/site-packages/")
from paraview.simple import GetActiveViewOrCreate, GetLight, AddLight
try:
import cftime as ct
except:
try:
import netcdftime as ct
except:
print("Need the python cftime (or the older netcdftime) module for the NetCDFTimeAnnotation plugin!", file=sys.stderr)
import datetime as dt
import numpy as np
inp = self.GetInputDataObject(0,0)
currentTime = inp.GetInformation().Get(vtk.vtkDataObject.DATA_TIME_STEP())
sdate=vtk.vtkStringArray()
timeUnitsArray= inp.GetFieldData().GetAbstractArray("time_units")
if timeUnitsArray:
timeUnits = timeUnitsArray.GetValue(0)
print ("Time units ", timeUnits)
cdftime = ct.utime(timeUnits)
t = cdftime.num2date(currentTime)
print ("currentTime ", currentTime, ", t = ", t)
first = ct.JulianDayFromDate(dt.datetime(t.year,1,1,0,0,0))
jd = (ct.JulianDayFromDate(t) - first + 1)
angle = 23.45 *np.pi/180 * np.sin(2 * np.pi * (284 + jd ) /365.25)
time = jd-np.floor(jd)
x = - np.cos(time * 2 * np.pi)
y = np.sin(time * 2 * np.pi)
z = np.tan(angle)
print (jd, time, x, y , z)
renderView1 = GetActiveViewOrCreate(\'RenderView\')
# Create a new \'Light\'
light = GetLight(0, view=renderView1)
if light is None:
light = AddLight(view=renderView1)
dist = 1e6
light.Position = [x*dist, y*dist, z*dist]
light.DiffuseColor = [1.0, 1.0, 1.0]"""
lightadjust.RequestInformationScript = ''
lightadjust.RequestUpdateExtentScript = ''
lightadjust.PythonPath = ''
# create a new 'Programmable Annotation'
dummylabel = ProgrammableAnnotation(registrationName='dummy-label', Input=lightadjust)
dummylabel.Script = """to = self.GetTableOutput()
arr = vtk.vtkStringArray()
arr.SetName("Text")
arr.SetNumberOfComponents(1)
arr.InsertNextValue("")
to.AddColumn(arr)"""
dummylabel.PythonPath = ''
dummylabelDisplay = Show(dummylabel, renderView1, 'TextSourceRepresentation')
print ("Done adjusting light", file=sys.stderr)
print ("Setting up textures", file=sys.stderr)
# create a new 'Texture Map to Sphere'
night = TextureMaptoSphere(registrationName='night', Input=rsdt)
night.PreventSeam = 0
# create a new 'Texture Map to Sphere'
day = TextureMaptoSphere(registrationName='day', Input=rsdt)
day = TextureMaptoSphere(registrationName='day', Input=cloud_psl)
day.PreventSeam = 0
# ----------------------------------------------------------------
......@@ -414,56 +328,8 @@ print ("Done setting up clouds", file=sys.stderr)
print ("Activating textures", file=sys.stderr)
# show data from night
nightDisplay = Show(night, renderView1, 'UnstructuredGridRepresentation')
# get color transfer function/color map for 'rsdt'
rsdtLUT = GetColorTransferFunction('rsdt')
rsdtLUT.AutomaticRescaleRangeMode = 'Never'
rsdtLUT.EnableOpacityMapping = 1
rsdtLUT.RGBPoints = [0.0, 1.0, 1.0, 1.0, 10.0, 1.0, 1.0, 1.0]
rsdtLUT.ScalarRangeInitialized = 1.0
# a texture
blackMarble_2016_3km_for_pv = CreateTexture('/home/k/k202134/Paraview/BlackMarble_2016_3km_for_pv-dark.png')
# get opacity transfer function/opacity map for 'rsdt'
rsdtPWF = GetOpacityTransferFunction('rsdt')
rsdtPWF.Points = [0.0, 1.0, 0.5, 0.0, 10.0, 0.0, 0.5, 0.0]
rsdtPWF.ScalarRangeInitialized = 1
# trace defaults for the display properties.
nightDisplay.Representation = 'Surface'
nightDisplay.ColorArrayName = ['CELLS', 'rsdt']
nightDisplay.LookupTable = rsdtLUT
nightDisplay.Specular = 0.99
nightDisplay.SpecularPower = 10.0
nightDisplay.Ambient = 0.33
nightDisplay.SelectTCoordArray = 'Texture Coordinates'
nightDisplay.SelectNormalArray = 'None'
nightDisplay.SelectTangentArray = 'None'
nightDisplay.Texture = blackMarble_2016_3km_for_pv
nightDisplay.OSPRayScaleArray = 'Texture Coordinates'
nightDisplay.OSPRayScaleFunction = 'PiecewiseFunction'
nightDisplay.SelectOrientationVectors = 'None'
nightDisplay.ScaleFactor = 40.0
nightDisplay.SelectScaleArray = 'alpha'
nightDisplay.GlyphType = 'Arrow'
nightDisplay.GlyphTableIndexArray = 'alpha'
nightDisplay.GaussianRadius = 2.0
nightDisplay.SetScaleArray = ['POINTS', 'Texture Coordinates']
nightDisplay.ScaleTransferFunction = 'PiecewiseFunction'
nightDisplay.OpacityArray = ['POINTS', 'Texture Coordinates']
nightDisplay.OpacityTransferFunction = 'PiecewiseFunction'
nightDisplay.DataAxesGrid = 'GridAxesRepresentation'
nightDisplay.PolarAxes = 'PolarAxesRepresentation'
nightDisplay.ScalarOpacityFunction = rsdtPWF
nightDisplay.ScalarOpacityUnitDistance = 10.049316464272104
nightDisplay.OpacityArrayName = ['CELLS', 'alpha']
nightDisplay.BumpMapInputDataArray = [None, '']
nightDisplay.ExtrusionInputDataArray = ['CELLS', 'alpha']
nightDisplay.SelectInputVectors = ['POINTS', 'Texture Coordinates']
nightDisplay.WriteLog = ''
print ("Done activating textures", file=sys.stderr)
......@@ -530,7 +396,7 @@ if options.get("psl", False):
# ----------------------------------------------------------------
# restore active source
SetActiveSource(rsdt)
SetActiveSource(cloud_psl)
# ----------------------------------------------------------------
......
......@@ -129,7 +129,7 @@ fn = options["clouds"]
if "UM" in options["clouds"] or "SAM" in options["clouds"] or "NICAM" in options["clouds"]:
if ":" in fn:
files = fn.split(":")
readers = [ NetCDFReader(registrationName='cloud_psl_%d'%n, FileName=[x]) for n, x in enumerate( files ) ]
readers = [ NetCDFReader(registrationName='cloud_psl_%d'%n, FileName=x.split(",")) for n, x in enumerate( files ) ]
for x in readers:
x.Dimensions = '(latitude, longitude)'
if "NICAM" in options["clouds"]:
......@@ -137,12 +137,12 @@ if "UM" in options["clouds"] or "SAM" in options["clouds"] or "NICAM" in options
x.VerticalScale = 200.0
cloud_psl = AppendAttributes(registrationName='cloud_psl', Input=readers )
else:
cloud_psl = NetCDFReader(registrationName='cloud_psl', FileName=[options["clouds"]])
cloud_psl = NetCDFReader(registrationName='cloud_psl', FileName=options["clouds"].split(","))
cloud_psl.VerticalScale = 200.0
else:
if ":" in fn:
files = fn.split(":")
readers = [ CDIReader(registrationName='cloud_psl_%d'%n, FileNames=[x]) for n, x in enumerate( files ) ]
readers = [ CDIReader(registrationName='cloud_psl_%d'%n, FileNames=x.split(",")) for n, x in enumerate( files ) ]
for n, x in enumerate(readers):
x.Dimensions = '(clon, clat, sfc)'
x.CellArrayStatus = [options['cloud_vars'][n]]
......@@ -152,7 +152,7 @@ else:
else:
# create a new 'CDIReader'
cloud_psl = CDIReader(registrationName='cloud_psl', FileNames=[options["clouds"]])
cloud_psl = CDIReader(registrationName='cloud_psl', FileNames=options["clouds"].split(","))
cloud_psl.Dimensions = '(clon, clat, sfc)'
cloud_psl.CellArrayStatus = options['cloud_vars']
cloud_psl.SetProjection = 'Spherical Projection'
......@@ -160,7 +160,7 @@ else:
if options.get("psl", False):
print ("Getting psl", file=sys.stderr)
psl = CDIReader(registrationName='psl', FileNames=[options["psl"]])
psl = CDIReader(registrationName='psl', FileNames=options["psl"].split(","))
psl.Dimensions = '(clon, clat, sfc)'
psl.CellArrayStatus = [options['psl_var']]
psl.SetProjection = 'Spherical Projection'
......@@ -194,7 +194,7 @@ alpha.Function = ' (3/2 * ({cloud_expr}) *100)/ (3/2 * ({cloud_expr}) *100+7)'.f
# create a new 'CDIReader'
print ("Getting light", file=sys.stderr)
rsdt_file = options["light"]
rsdt = CDIReader(registrationName='rsdt', FileNames=[rsdt_file])
rsdt = CDIReader(registrationName='rsdt', FileNames=rsdt_file.split(","))
rsdt.Dimensions = '(clon, clat, toa)'
rsdt.CellArrayStatus = ['rsdt']
rsdt.SetProjection = 'Spherical Projection'
......@@ -202,11 +202,6 @@ rsdt.LayerThickness = 50
rsdt.MaskingValueVar = 'rsdt'
print ("Done getting light", file=sys.stderr)
alpha = Calculator(registrationName='alpha', Input=rsdt)
alpha.AttributeType = 'Cell Data'
alpha.ResultArrayName = 'alpha'
alpha.Function = 'rsdt'
# create a new 'NetCDF Time Annotation'
netCDFTimeAnnotation1 = NetCDFTimeAnnotation(registrationName='NetCDFTimeAnnotation1', Input=rsdt)
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment