day and night with automatic solar angle

DYAMOND_WINTER/day-night.py

@@ -141,6 +141,70 @@ rsdt.SetProjection = 'Spherical Projection'
 rsdt.LayerThickness = 50
 rsdt.MaskingValueVar = 'rsdt'
 
+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)
+cdftime = ct.utime(timeUnits)
+t = cdftime.num2date(currentTime)
+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')
+
+
+
 # create a new 'Texture Map to Sphere'
 night = TextureMaptoSphere(registrationName='night', Input=rsdt)
 night.PreventSeam = 0
@@ -268,7 +332,7 @@ 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/BlackMarble_2016_3km_for_pv.jpg')
+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')
@@ -377,7 +441,7 @@ if "SHiELD" in options["clouds"] :
 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]+'_white-dark.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])