Commit acb5409c authored by Florian Ziemen's avatar Florian Ziemen
Browse files

mistral version of day-night

parent f78c66de
......@@ -12,8 +12,9 @@ print (sys.argv)
from argparse import ArgumentParser
import glob
r = 126400.
def sphere_proj (lat, lon) :
r = 1264.
global r
theta = (90 - lat) * pi / 180.
phi = lon * pi / 180.
z = r * cos (theta)
......@@ -31,9 +32,11 @@ def parse_args():
parser.add_argument("--cloud_expr")
parser.add_argument("--cloud_vars")
parser.add_argument("--psl",)
parser.add_argument("--psl_var",)
parser.add_argument("--light")
parser.add_argument("--lonlat", nargs = 2, type=float)
parser.add_argument("--output", "-o")
parser.add_argument("--no_time_string", action="store_true")
op = parser.parse_args()
options = vars(op)
print (options)
......@@ -44,6 +47,9 @@ def parse_args():
if not options.get("output", False) :
options["output"] = options["clouds"][:-3]
if not options.get("psl_var", False) :
options["psl_var"] = "psl"
print (options)
return options
......@@ -64,7 +70,7 @@ light1.Position = [50.0, -87.0, -28.6]
# Create a new 'Light'
light2 = CreateLight()
light2.Coords = 'Ambient'
light2.Intensity = 0.25
light2.Intensity = 0.3
light2.Position = [0.0, -100.0, 0.0]
light2.DiffuseColor = [0.7529411764705882, 0.8196078431372549, 1.0]
......@@ -89,7 +95,7 @@ renderView1.KeyLightIntensity = 0.0
renderView1.StereoType = 'Crystal Eyes'
renderView1.CameraPosition = [1144.0024360966738, -682.5869061407763, 142.4868048724324]
renderView1.CameraViewUp = [0, 0., 1]
renderView1.CameraViewAngle = 20.0
renderView1.CameraViewAngle = 20. * 1264. / r
renderView1.CameraFocalDisk = 1.0
renderView1.CameraParallelScale = 564.2333216713644
renderView1.EyeAngle = 1.0
......@@ -120,17 +126,18 @@ SetActiveView(renderView1)
# ----------------------------------------------------------------
fn = options["clouds"]
if "UM" in options["clouds"] or "SAM" in 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 ) ]
for x in readers:
x.Dimensions = '(latitude, longitude)'
if "NICAM" in options["clouds"]:
x.Dimensions = '(lat, lon)'
x.VerticalScale = 200.0
cloud_psl = AppendAttributes(registrationName='cloud_psl', Input=readers )
else:
cloud_psl = NetCDFReader(registrationName='cloud_psl', FileName=[options["clouds"]])
cloud_psl.Dimensions = '(latitude, longitude)'
cloud_psl.VerticalScale = 200.0
else:
if ":" in fn:
......@@ -155,39 +162,22 @@ if options.get("psl", False):
print ("Getting psl", file=sys.stderr)
psl = CDIReader(registrationName='psl', FileNames=[options["psl"]])
psl.Dimensions = '(clon, clat, sfc)'
psl.CellArrayStatus = ['psl']
psl.CellArrayStatus = [options['psl_var']]
psl.SetProjection = 'Spherical Projection'
psl.LayerThickness = 50
psl.MaskingValueVar = 'psl'
# create a new 'Cell Data to Point Data'
psl_to_point = CellDatatoPointData(registrationName='psl_to_point', Input=psl)
psl_to_point.ProcessAllArrays = 0
psl_to_point.CellDataArraytoprocess = ['psl']
psl_to_point.CellDataArraytoprocess = [options['psl_var']]
# create a new 'Contour'
psl_contour = Contour(registrationName='psl_contour', Input=psl_to_point)
psl_contour.ContourBy = ['POINTS', 'psl']
psl_contour.ContourBy = ['POINTS', options['psl_var'] ]
psl_contour.Isosurfaces = [99000.0, 100500.0, 102000.0, 103500.0]
psl_contour.PointMergeMethod = 'Uniform Binning'
print ("Done getting psl", file=sys.stderr)
# create a new 'NetCDF Time Annotation'
netCDFTimeAnnotation1 = NetCDFTimeAnnotation(registrationName='NetCDFTimeAnnotation1', Input=cloud_psl)
netCDFTimeAnnotation1.Expression = '"%02i-%02i-%02i %02i:%02i" % (Date[0], Date[1], Date[2], Date[3], Date[4])'
# show data in view
netCDFTimeAnnotation1Display = Show(netCDFTimeAnnotation1, renderView1, 'TextSourceRepresentation')
# Properties modified on netCDFTimeAnnotation1Display
netCDFTimeAnnotation1Display.FontSize = 96
# Properties modified on netCDFTimeAnnotation1
# Properties modified on netCDFTimeAnnotation1Display
netCDFTimeAnnotation1Display.FontFamily = 'Courier'
......@@ -198,19 +188,40 @@ alpha.ResultArrayName = 'alpha'
cloud_expr = options.get("cloud_expr", False)
if not cloud_expr:
cloud_expr = " ( 5* clivi + clwvi ) "
alpha.Function = ' (3/2 * {cloud_expr} *100)/ (3/2 * {cloud_expr} *100+7)'.format(cloud_expr=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])
rsdt.Dimensions = '(clon, clat, sfc)'
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)
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)
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):
# show data in view
netCDFTimeAnnotation1Display = Show(netCDFTimeAnnotation1, renderView1, 'TextSourceRepresentation')
netCDFTimeAnnotation1Display.FontSize = 96
netCDFTimeAnnotation1Display.FontFamily = 'Courier'
print ("Adjusting light", file=sys.stderr)
lightadjust = ProgrammableFilter(registrationName='light-adjust', Input=rsdt)
lightadjust.OutputDataSetType = 'vtkTable'
......@@ -236,8 +247,10 @@ 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)
......@@ -465,37 +478,37 @@ if options.get("psl", False):
print ("... initialized ... ", file=sys.stderr)
# get color transfer function/color map for 'psl'
pslLUT = GetColorTransferFunction('psl')
pslLUT = GetColorTransferFunction(options['psl_var'])
pslLUT.AutomaticRescaleRangeMode = 'Never'
pslLUT.RGBPoints = [100000.0, 0.23137254902, 0.298039215686, 0.752941176471, 101250.0, 0.865, 0.865, 0.865, 102500.0, 0.705882352941, 0.0156862745098, 0.149019607843]
pslLUT.ScalarRangeInitialized = 1.0
# trace defaults for the display properties.
psl_contourDisplay.Representation = 'Surface'
psl_contourDisplay.ColorArrayName = ['POINTS', 'psl']
psl_contourDisplay.ColorArrayName = ['POINTS', options['psl_var']]
psl_contourDisplay.LookupTable = pslLUT
psl_contourDisplay.LineWidth = 4.0
psl_contourDisplay.SelectTCoordArray = 'None'
psl_contourDisplay.SelectNormalArray = 'None'
psl_contourDisplay.SelectTangentArray = 'None'
psl_contourDisplay.Scale = [1.001, 1.001, 1.001]
psl_contourDisplay.OSPRayScaleArray = 'psl'
psl_contourDisplay.OSPRayScaleArray = options['psl_var']
psl_contourDisplay.OSPRayScaleFunction = 'PiecewiseFunction'
psl_contourDisplay.SelectOrientationVectors = 'None'
psl_contourDisplay.ScaleFactor = 39.99499206542969
psl_contourDisplay.SelectScaleArray = 'psl'
psl_contourDisplay.SelectScaleArray = options['psl_var']
psl_contourDisplay.GlyphType = 'Arrow'
psl_contourDisplay.GlyphTableIndexArray = 'psl'
psl_contourDisplay.GlyphTableIndexArray = options['psl_var']
psl_contourDisplay.GaussianRadius = 1.9997496032714843
psl_contourDisplay.SetScaleArray = ['POINTS', 'psl']
psl_contourDisplay.SetScaleArray = ['POINTS', options['psl_var']]
psl_contourDisplay.ScaleTransferFunction = 'PiecewiseFunction'
psl_contourDisplay.OpacityArray = ['POINTS', 'psl']
psl_contourDisplay.OpacityArray = ['POINTS', options['psl_var']]
psl_contourDisplay.OpacityTransferFunction = 'PiecewiseFunction'
psl_contourDisplay.DataAxesGrid = 'GridAxesRepresentation'
psl_contourDisplay.PolarAxes = 'PolarAxesRepresentation'
psl_contourDisplay.BumpMapInputDataArray = ['POINTS', 'psl']
psl_contourDisplay.ExtrusionInputDataArray = ['POINTS', 'psl']
psl_contourDisplay.SelectInputVectors = ['POINTS', 'psl']
psl_contourDisplay.BumpMapInputDataArray = ['POINTS', options['psl_var']]
psl_contourDisplay.ExtrusionInputDataArray = ['POINTS', options['psl_var']]
psl_contourDisplay.SelectInputVectors = ['POINTS', options['psl_var']]
psl_contourDisplay.WriteLog = ''
# init the 'PiecewiseFunction' selected for 'ScaleTransferFunction'
......@@ -512,8 +525,8 @@ if options.get("psl", False):
# note: the Get..() functions create a new object, if needed
# ----------------------------------------------------------------
# get opacity transfer function/opacity map for 'psl'
pslPWF = GetOpacityTransferFunction('psl')
# get opacity transfer function/opacity map for options['psl_var']
pslPWF = GetOpacityTransferFunction(options['psl_var'])
pslPWF.Points = [100000.0, 0.0, 0.5, 0.0, 102500.0, 1.0, 0.5, 0.0]
pslPWF.ScalarRangeInitialized = 1
print ("Done displaying psl", file=sys.stderr)
......@@ -528,7 +541,8 @@ if "SHiELD" in options["clouds"] :
psl_contour.Isosurfaces = [ x / 100. for x in psl_contour.Isosurfaces ]
pslLUT.RGBPoints = [1000.0, 0.23137254902, 0.298039215686, 0.752941176471, 1012.500, 0.865, 0.865, 0.865, 1025.0, 0.705882352941, 0.0156862745098, 0.149019607843]
alphaLUT.RescaleTransferFunction(0.0, 1.0)
alphaPWF.RescaleTransferFunction(0.0, 1.0)
if __name__ == '__main__':
view = GetActiveView()
......@@ -543,7 +557,8 @@ if __name__ == '__main__':
return "{fn}_{atts}_{n:03d}.png".format(fn=fn, n=n, atts=attstr)
for n,t in enumerate (tsteps):
print ("rendering for time %f"%t)
view.ViewTime = t
if len (tsteps) > 1:
view.ViewTime = t
# generate extracts
# renderView1.Background = [1.0, 1.0, 1.0]
......@@ -558,7 +573,8 @@ if __name__ == '__main__':
lonlat = [-30, 0]
view.CameraPosition = sphere_proj(lon=lonlat[0], lat=lonlat[1])
ll_string="%.0fN_%.0fE"%(lonlat[1],lonlat[0])
view.CameraParallelScale = 120
# view.CameraParallelScale = 120
SaveScreenshot(genname (options["output"], n, [ll_string, "transparent", "4k"]), renderView1, ImageResolution=[3840, 2160], TransparentBackground=1)
view.CameraPosition = oldpos
ll_string=""
......@@ -568,7 +584,7 @@ if __name__ == '__main__':
lonlat = options["lonlat"]
view.CameraPosition = sphere_proj(lon=lonlat[0], lat=lonlat[1])
ll_string="%.0fN_%.0fE"%(lonlat[1],lonlat[0])
view.CameraParallelScale = 120
# view.CameraParallelScale = 120
SaveScreenshot(genname (options["output"], n, [ll_string, "transparent", "4k"]), renderView1, ImageResolution=[3840, 2160], TransparentBackground=1)
view.CameraPosition = oldpos
ll_string=""
......
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