messing with readers -- all else should be fine

DYAMOND_WINTER/day-night.py

@@ -12,6 +12,15 @@ print (sys.argv)
 from argparse import ArgumentParser
 import glob
 
+def sphere_proj (lat, lon) :
+    r = 1400
+    theta = (90 - lat) * pi / 180.
+    phi = lon * pi / 180.
+    z = r * cos (theta)
+    x = r * sin(theta) * cos (phi)
+    y = r * sin(theta) * sin (phi)
+    return (x, y, z)
+
 options = {}
 def parse_args():
     '''Parses the command line arguments'''
@@ -19,15 +28,25 @@ def parse_args():
     parser = ArgumentParser()
     parser.description = "Plot icon files"
     parser.add_argument("--clouds")
+    parser.add_argument("--cloud_expr")
+    parser.add_argument("--cloud_vars")
     parser.add_argument("--psl",)
     parser.add_argument("--light")
+    parser.add_argument("--lonlat", nargs = 2, type=float)
     op = parser.parse_args()
     options = vars(op)
     print (options)
+    if not options.get("cloud_vars", False) :
+        options["cloud_vars"] = "cloud"
+    options["cloud_vars"] = options["cloud_vars"].split(",")
+    print (options)
     return options
 
 options=parse_args()
 
+
+
+
 # ----------------------------------------------------------------
 # setup views used in the visualization
 # ----------------------------------------------------------------
@@ -95,20 +114,28 @@ SetActiveView(renderView1)
 # setup the data processing pipelines
 # ----------------------------------------------------------------
 
-if "UM" in options["clouds"]:
-    cloud_psl = NetCDFReader(registrationName='cloud_psl', FileName=[options["clouds"]])
-    cloud_psl.Dimensions = '(latitude, longitude)'
-    cloud_psl.VerticalScale = 200.0
+if "UM" in options["clouds"] or "SAM" in options["clouds"] :
+    fn = options["clouds"]
+    if ":" in fn:
+        files = fn.split(":")
+        readers = [ NetCDFReader(registrationName='cloud_psl_%d'%n, FileName=[x]) for n, x in enumerate( files) ]
+        [ x.Dimensions = '(latitude, longitude)' for x in readers]
+        [ x.VerticalScale = 200.0 for x in readers ]
+        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:
     # create a new 'CDIReader'
     cloud_psl = CDIReader(registrationName='cloud_psl', FileNames=[options["clouds"]])
     cloud_psl.Dimensions = '(clon, clat, sfc)'
-    cloud_psl.CellArrayStatus = ['cloud']
+    cloud_psl.CellArrayStatus = options['cloud_vars']
     cloud_psl.SetProjection = 'Spherical Projection'
     cloud_psl.LayerThickness = 50
-    cloud_psl.MaskingValueVar = 'cloud'
 
 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']
@@ -125,6 +152,7 @@ if options.get("psl", False):
     psl_contour.ContourBy = ['POINTS', 'psl']
     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)
@@ -134,7 +162,10 @@ netCDFTimeAnnotation1.Expression = '"On %02i.%02i.%02i at %02i:%02i" % (Date[0],
 alpha = Calculator(registrationName='alpha', Input=cloud_psl)
 alpha.AttributeType = 'Cell Data'
 alpha.ResultArrayName = 'alpha'
-alpha.Function = ' (3/2 * cloud *100)/ (3/2 * cloud *100+7)'
+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)
 
 
 # create a new 'CDIReader'
@@ -378,8 +409,10 @@ nightDisplay.SelectInputVectors = ['POINTS', 'Texture Coordinates']
 nightDisplay.WriteLog = ''
 
 if options.get("psl", False):
+    print ("Displaying psl", file=sys.stderr)
     # show data from psl_contour
     psl_contourDisplay = Show(psl_contour, renderView1, 'GeometryRepresentation')
+    print ("... initialized ... ", file=sys.stderr)
 
     # get color transfer function/color map for 'psl'
     pslLUT = GetColorTransferFunction('psl')
@@ -433,10 +466,11 @@ if options.get("psl", False):
     pslPWF = GetOpacityTransferFunction('psl')
     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)
 
 # ----------------------------------------------------------------
 # restore active source
-SetActiveSource(alpha)
+SetActiveSource(cloud_psl)
 # ----------------------------------------------------------------
 
 
@@ -445,13 +479,34 @@ if "SHiELD" in options["clouds"] :
     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]
 
 if __name__ == '__main__':
-    # generate extracts
-    renderView1.Background = [1.0, 1.0, 1.0]
-    SaveScreenshot(options["clouds"][:-3]+'_white.png', renderView1, ImageResolution=[1920, 1080])
-    SaveScreenshot(options["clouds"][:-3]+'_white4k.png', renderView1, ImageResolution=[3840, 2160])  # TransparentBackground=1
-    renderView1.Background = [0., 0., 0.]
-    SaveScreenshot(options["clouds"][:-3]+'_black.png', renderView1, ImageResolution=[1920, 1080])
-    SaveScreenshot(options["clouds"][:-3]+'_black4k.png', renderView1, ImageResolution=[3840, 2160])
-    #
-    SaveScreenshot(options["clouds"][:-3]+'_transparent.png', renderView1, ImageResolution=[1920, 1080], TransparentBackground=1)  
-    SaveScreenshot(options["clouds"][:-3]+'_transparent4k.png', renderView1, ImageResolution=[3840, 2160], TransparentBackground=1)  
+    view = GetActiveView()
+    if options.get("lonlat", False ):
+        lonlat = options["lonlat"]
+        view.CameraPosition = sphere_proj(lon=lonlat[0], lat=lonlat[1])
+        ll_string="%.0fN_%.0fE"%(lonlat[1],lonlat[0])
+        view.CameraParallelScale = 120
+
+    else:
+        ll_string=""
+    print ('in main', file=sys.stderr)
+    reader=GetActiveSource()
+    tsteps=reader.TimestepValues
+    print ("Timesteps: ", ", ".join((str (x) for x in tsteps)))
+    def genname(fn, n, atts):
+        attstr = "_".join((x for x in atts if x))
+        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
+
+        # generate extracts
+        renderView1.Background = [1.0, 1.0, 1.0]
+        SaveScreenshot(genname (options["clouds"][:-3], n, [ll_string, "white", "HD"]), renderView1, ImageResolution=[1920, 1080])
+        SaveScreenshot(genname (options["clouds"][:-3], n, [ll_string, "white", "4k"]), renderView1, ImageResolution=[3840, 2160])  # TransparentBackground=1
+        renderView1.Background = [0., 0., 0.]
+        SaveScreenshot(genname (options["clouds"][:-3], n, [ll_string, "black", "HD"]), renderView1, ImageResolution=[1920, 1080])
+        SaveScreenshot(genname (options["clouds"][:-3], n, [ll_string, "black", "4k"]), renderView1, ImageResolution=[3840, 2160])  # TransparentBackground=1
+        #
+        SaveScreenshot(genname (options["clouds"][:-3], n, [ll_string, "transparent", "HD"]), renderView1, ImageResolution=[1920, 1080], TransparentBackground=1)
+        SaveScreenshot(genname (options["clouds"][:-3], n, [ll_string, "transparent", "4k"]), renderView1, ImageResolution=[3840, 2160], TransparentBackground=1)
+