1.6.0b

PyMieScatt.egg-info/PKG-INFO

@@ -1,13 +1,13 @@
 Metadata-Version: 1.1
 Name: PyMieScatt
-Version: 1.6.0a0
+Version: 1.6.0b0
 Summary: A collection of forward and inverse Mie solving routines based on Bohren and Huffman's Mie Theory derivations.
 Home-page: http://air.eece.wustl.edu/people/ben-sumlin/
 Author: Benjamin Sumlin
 Author-email: bsumlin@wustl.edu
 License: GPL
 Description-Content-Type: UNKNOWN
-Description: 1.6.0a - See `documentation <http://pymiescatt.readthedocs.io/>`_ for update notes.
+Description: 1.6.0b - See `documentation <http://pymiescatt.readthedocs.io/>`_ for update notes.
 Keywords: Mie Rayleigh scattering absorption extinction light refraction
 Platform: UNKNOWN
 Classifier: Development Status :: 5 - Production/Stable

PyMieScatt/Inverse.py

@@ -74,7 +74,7 @@ def Inversion_SD(Bsca,Babs,wavelength,dp,ndp,nMin=1,nMax=3,kMin=0,kMax=1,scatter
 
 def ContourIntersection(Qsca,Qabs,wavelength,diameter,Qback=None,n=None,k=None,nMin=1,nMax=3,kMin=0.00001,kMax=1,gridPoints=100,interpolationFactor=2,maxError=0.005,fig=None,ax=None,axisOption=0):
 #  http://pymiescatt.readthedocs.io/en/latest/inverse.html#ContourIntersectio
-  if Qabs == 0.0 or Qabs[0] == 0.0:
+  if (type(Qabs) == np.float64 and Qabs == 0.0) or (type(Qabs) in [list, tuple, np.ndarray] and Qabs[0]==0):
     k = 0.0
   if k == 0.0:
     kMin = -0.1
@@ -313,7 +313,7 @@ def ContourIntersection(Qsca,Qabs,wavelength,diameter,Qback=None,n=None,k=None,n
 
 def ContourIntersection_SD(Bsca,Babs,wavelength,dp,ndp,n=None,k=None,nMin=1,nMax=3,kMin=0.00001,kMax=1,Bback=None,gridPoints=60,interpolationFactor=2,maxError=0.005,fig=None,ax=None,axisOption=0):
 #  http://pymiescatt.readthedocs.io/en/latest/inverse.html#ContourIntersection_SD
-  if Babs == 0.0 or Babs[0] == 0.0:
+  if (type(Babs) == np.float64 and Babs == 0.0) or (type(Babs) in [list, tuple, np.ndarray] and Babs[0]==0):
     k = 0.0
   if k == 0.0:
     kMin = -0.1

PyMieScatt/_version.py

@@ -1 +1 @@
-__version__ = "1.6.0a"
+__version__ = "1.6.0b"

PyMieScatt/devmode/functionPrototyping.py

@@ -1,329 +1,10 @@
 # -*- coding: utf-8 -*-
 # http://pymiescatt.readthedocs.io/en/latest/inverse.html
-from PyMieScatt.Mie import Mie_ab
-import numpy as np
-import matplotlib.pyplot as plt
-from matplotlib.offsetbox import AnchoredText
-from matplotlib.contour import QuadContourSet
-from matplotlib.collections import LineCollection
-from scipy.ndimage import zoom
-from scipy.integrate import trapz
-from shapely import geometry
+import PyMieScatt as ps
+m,w,d = 1.5+0.5j,532,500
+dp = (15.7,16.3,16.8,17.5,18.1,18.8,19.5,20.2,20.9,21.7,22.5,23.3,24.1,25,25.9,26.9,27.9,28.9,30,31.1,32.2,33.4,34.6,35.9,37.2,38.5,40,41.4,42.9,44.5,46.1,47.8,49.6,51.4,53.3,55.2,57.3,59.4,61.5,63.8,66.1,68.5,71,73.7,76.4,79.1,82,85.1,88.2,91.4,94.7,98.2,101.8,105.5,109.4,113.4,117.6,121.9,126.3,131,135.8,140.7,145.9,151.2,156.8,162.5,168.5,174.7,181.1,187.7,194.6,201.7,209.1,216.7,224.7,232.9,241.4,250.3,259.5,269,278.8,289,299.6,310.6,322,333.8,346,358.7,371.8,385.4,399.5,414.2,429.4,445.1,461.4,478.3,495.8,514,532.8,552.3,572.5,593.5,615.3,637.8,661.2,685.4,710.5,736.5,763.5,791.5,820.5,850.5,881.7,914,947.5,982.2)
+ndp = (0,0,0,0,20.007,22.604,28.565,30.236,35.303,40.97,44.946,53.688,57.714,70.388,75.05,76.691,89.164,102.936,112.255,113.609,132.423,151.342,164.092,172.006,195.124,213.077,233.145,262.117,279.707,313.686,346.505,386.169,422.061,472.318,531.711,570.058,615.129,689.183,735.801,796.658,866.776,937.752,1003.532,1074.711,1141.764,1189.764,1288.912,1336.65,1418.396,1415.48,1543.714,1610.163,1644.419,1646.855,1686.363,1811.085,1831.266,1840.095,1852.131,1874.75,1884.19,1882.413,1854.246,1830.001,1808.497,1784.443,1753.237,1695.941,1630.566,1600.728,1526.502,1475.664,1408.608,1339.624,1257.692,1193.99,1123.253,1064.108,1002.481,937.285,865.517,799.852,743.185,704.676,640.034,586.251,547.975,495.346,458.38,427.772,389.496,356.098,332.884,302.85,270.998,250.737,223.463,212.376,178.743,157.359,144.926,127.777,111.063,98.037,87.379,77.376,65.592,58.967,46.849,43.639,35.027,31.979,25.907,20.387,16.343,14.932)
 
-def find_intersections(A,B):
-  X = []
-  Y = []
-  if isinstance(A, QuadContourSet):
-    pathOne = A.collections[0]
-  elif isinstance(A, LineCollection):
-    pathOne = A
-  if isinstance(B, QuadContourSet):
-    pathTwo = B.collections[0]
-  elif isinstance(B, LineCollection):
-    pathTwo = B
-  for p1 in pathOne.get_paths():
-    for p2 in pathTwo.get_paths():
-      v1 = p1.vertices
-      v2 = p2.vertices
-
-#      try:
-      poly1 = geometry.LineString(v1)
-      poly2 = geometry.LineString(v2)
-      sol = poly1.intersection(poly2)
-      if type(sol)==geometry.point.Point:
-        X.append(sol.x)
-        Y.append(sol.y)
-      else:
-        for s in sol:
-          X.append(s.x)
-          Y.append(s.y)
-#      except:
-#        pass
+Bsca,Babs,Bback = ps.fastMie_SD(m,w,dp,ndp)
 
-  return X,Y
-
-def fastMieQ(m, wavelength, diameter):
-#  http://pymiescatt.readthedocs.io/en/latest/inverse.html#fastMieQ
-  x = np.pi*diameter/wavelength
-  if x==0:
-    return 0, 0, 0
-  elif x>0:
-    nmax = np.round(2+x+4*(x**(1/3)))
-    n = np.arange(1,nmax+1)
-    n1 = 2*n+1
-    x2 = x**2
-    an,bn = Mie_ab(m,x)
-    qext = (2/x2)*np.sum(n1*(an.real+bn.real))
-    qsca = (2/x2)*np.sum(n1*(an.real**2+an.imag**2+bn.real**2+bn.imag**2))
-    qback = (1/x2)*(np.abs(np.sum(n1*((-1)**n)*(an-bn)))**2)
-    qabs = qext-qsca
-    return qsca, qabs, qback
-
-plt.close('all')
-
-m = 1.5+0.005j
-
-diameter = 500#np.random.uniform(200,1000)
-wavelength = 500#np.random.uniform(300,700)
-
-Q = fastMieQ(m,wavelength,diameter)
-
-Q = [[x,x*0.01] for x in Q]
-
-Qsca = Q[0]
-Qabs = Q[1]
-#Qback = Q[2]
-Qback=None
-n=None
-k=None
-nMin=1
-nMax=2
-kMin=0.00001
-kMax=1
-gridPoints=100
-interpolationFactor=2
-maxError=0.005
-fig=None
-ax=None
-axisOption=1
-
-def ContourIntersection(Qsca,Qabs,wavelength,diameter,Qback=None,n=None,k=None,nMin=1,nMax=3,kMin=0.00001,kMax=1,gridPoints=100,interpolationFactor=2,maxError=0.005,fig=None,ax=None,axisOption=0):
-#  http://pymiescatt.readthedocs.io/en/latest/inverse.html#ContourIntersectio
-  if Qabs == 0.0 or Qabs[0] == 0.0:
-    k = 0.0
-  if k == 0.0:
-    kMin = -0.1
-    axisOption = 1
-
-  error = lambda measured,calculated: np.abs((calculated-measured)/measured)
-  if Qback is not None:
-    if gridPoints*interpolationFactor<400:
-      gridPoints = 2*gridPoints
-  labels = []
-
-  incErrors = False
-  if type(Qsca) in [list, tuple, np.ndarray]:
-    incErrors = True
-    scaError = Qsca[1]
-    Qsca = Qsca[0]
-    labels.append("Qsca = {b:1.3f}±{e:1.3f}".format(b=Qsca,e=scaError))
-  else:
-    scaError = None
-    labels.append("Qsca = {b:1.3f}".format(b=Qsca))
-
-  if type(Qabs) in [list, tuple, np.ndarray]:
-    incErrors = True
-    absError = Qabs[1]
-    Qabs = Qabs[0]
-    labels.append("Qabs = {b:1.3f}±{e:1.3f}".format(b=Qabs,e=absError))
-  else:
-    absError = None
-    labels.append("Qabs = {b:1.3f}".format(b=Qabs))
-
-  if type(Qback) in [list, tuple, np.ndarray]:
-    backError = Qback[1]
-    Qback = Qback[0]
-    labels.append("Qback = {b:1.3f}±{e:1.3f}".format(b=Qback,e=backError))
-  elif Qback is not None:
-    backError = None
-    labels.append("Qback - {b:1.3f}".format(b=Qback))
-  else:
-    backError = None
-
-  nRange = np.linspace(nMin,nMax,gridPoints)
-  if k == 0.0:
-    kRange = np.linspace(kMin,kMax,gridPoints)
-  else:
-    kRange = np.logspace(np.log10(kMin),np.log10(kMax),gridPoints)
-  QscaList, QabsList, QbackList = [], [], []
-  for _n in nRange:
-    s, a, b = [], [], []
-    for _k in kRange:
-      m = _n+_k*1.0j
-      _Qsca,_Qabs,_Qback = fastMieQ(m,wavelength,diameter)
-      s.append(_Qsca)
-      a.append(_Qabs)
-      b.append(_Qback)
-    QscaList.append(s)
-    QabsList.append(a)
-    QbackList.append(b)
-  QscaList = zoom(np.transpose(np.array(QscaList)),interpolationFactor)
-  QabsList = zoom(np.transpose(np.array(QabsList)),interpolationFactor)
-  QbackList = zoom(np.transpose(np.array(QbackList)),interpolationFactor)
-
-  _n = zoom(nRange,interpolationFactor)
-  _k = zoom(kRange,interpolationFactor)
-
-  if fig is None and ax is None:
-    fig, ax = plt.subplots()
-  elif fig is None:
-    fig = ax.get_figure()
-  elif ax is None:
-    ax = fig.gca()
-
-  scaLevels = np.array([Qsca])
-  absLevels = np.array([Qabs])
-
-  if Qback is not None:
-    backLevels = np.array([Qback])
-    if backError is not None:
-      backErrorLevels = np.array([Qback+x for x in [-backError,backError]])
-
-  if n is None:
-    scaChart = ax.contour(_n,_k,QscaList,scaLevels,origin='lower',linestyles='dashdot',linewidths=1.5,colors=('red'))
-    if scaError is not None:
-      scaErrorLevels = np.array([Qsca+x for x in [-scaError, scaError]])
-      ax.contourf(_n,_k,QscaList,scaErrorLevels,origin='lower',colors=('red'),alpha=0.15)
-      ax.contour(_n,_k,QscaList,scaErrorLevels,origin='lower',linewidths=0.5,colors=('red'),alpha=0.5)
-  else:
-    if type(n) in [list, tuple, np.ndarray]:
-      scaErrorLevels = [n[0]*(1+x) for x in [-n[1],n[1]]]
-      scaChart = ax.vlines(n[0],kMin,kMax,linestyle='dashdot',linewidth=1.5,color='r')
-    else:
-      scaChart = ax.vlines(n,kMin,kMax,linestyle='dashdot',linewidth=1.5,color='r')
-
-  if k is None:
-    absChart = ax.contour(_n,_k,QabsList,absLevels,origin='lower',linewidths=1.5,colors=('blue'))
-    if absError is not None:
-      absErrorLevels = np.array([Qabs+x for x in [-absError, absError]])
-      ax.contourf(_n,_k,QabsList,absErrorLevels,origin='lower',colors=('blue'),alpha=0.15)
-      ax.contour(_n,_k,QabsList,absErrorLevels,origin='lower',linewidths=0.5,colors=('blue'),alpha=0.5)
-  else:
-    if type(k) in [list, tuple, np.ndarray]:
-      absErrorLevels = [k[0]*(1+x) for x in [-k[1],k[1]]]
-      absChart = ax.hlines(k[0],nMin,nMax,linestyle='solid',linewidth=1.5,color='b')
-    else:
-      absChart = ax.hlines(k,nMin,nMax,linestyle='solid',linewidth=1.5,color='b')
-
-  if Qback is not None:
-    backChart = ax.contour(_n,_k,QbackList,backLevels,origin='lower',linestyles='dotted',linewidths=1.5,colors=('green'))
-    if backError is not None:
-      backErrorLevels = np.array([Qback+x for x in [-backError, backError]])
-      ax.contourf(_n,_k,QbackList,backErrorLevels,origin='lower',colors=('green'),alpha=0.15)
-      ax.contour(_n,_k,QbackList,backErrorLevels,origin='lower',linewidths=0.5,colors=('green'),alpha=0.5)
-
-  m1 = find_intersections(scaChart,absChart)
-
-  if n is not None and type(n) in [list, tuple, np.ndarray]:
-    scaChart = ax.vlines(scaErrorLevels,kMin,kMax,linestyle='dashdot',linewidth=0.5,color='r',alpha=0.5)
-    ax.axvspan(scaErrorLevels[0],scaErrorLevels[1],alpha=0.15,color='r')
-  if k is not None and type(k) in [list, tuple, np.ndarray]:
-    absChart = ax.hlines(absErrorLevels,nMin,nMax,linestyle='solid',linewidth=0.5,color='b',alpha=0.5)
-    ax.axhspan(absErrorLevels[0],absErrorLevels[1],alpha=0.15,color='b')
-
-  if Qback is not None:
-    m2 = find_intersections(scaChart,backChart)
-    r1 = [np.round(x+y*1j,2) for x,y in zip(m1[0],m1[1])]
-    r2 = [np.round(x+y*1j,2) for x,y in zip(m2[0],m2[1])]
-    m_sol = list(set(r1).intersection(r2))
-    nSolution,kSolution = [xx.real for xx in m_sol],[xx.imag for xx in m_sol]
-  else:
-    nSolution,kSolution = m1[0],m1[1]
-
-  if type(nSolution)==np.float64:
-    solutionSet = [nSolution + (0+1j)*kSolution]
-  else:
-    solutionSet = [(x+y*1j) for x,y in zip(nSolution,kSolution)]
-
-  forwardCalculations = []
-  for s in solutionSet:
-    _s,_a,_ = fastMieQ(s,wavelength,diameter)
-    forwardCalculations.append([_s,_a])
-  solutionErrors = []
-  for f in forwardCalculations:
-    solutionErrors.append([error(f[0],Qsca),error(f[1],Qabs)])
-
-  solutionSet = np.array(solutionSet)
-  forwardCalculations = np.array(forwardCalculations)
-  solutionErrors = np.array(solutionErrors)
-
-  if n is None and k is None:
-    proper = solutionErrors <= maxError
-    solution = []
-    for r,c in proper:
-      if r and c:
-        solution.append(True)
-      else:
-        solution.append(False)
-
-    solutionSet = solutionSet[solution]
-    forwardCalculations = forwardCalculations[solution]
-    solutionErrors = solutionErrors[solution]
-    nSolutionsToPlot,kSolutionsToPlot = [x.real for x in solutionSet],[x.imag for x in solutionSet]
-  else:
-    nSolutionsToPlot,kSolutionsToPlot = m1[0],m1[1]
-
-  ax.scatter(nSolutionsToPlot,kSolutionsToPlot,marker='o',s=128,linewidth=1.5,edgecolor='k',facecolor='none',zorder=3)
-  ax.scatter(nSolutionsToPlot,kSolutionsToPlot,marker='o',s=128,linewidth=0,edgecolor='none',facecolor='c',zorder=1,alpha=0.25)
-
-  for x,y,s in zip(nSolutionsToPlot,kSolutionsToPlot,solutionErrors):
-    if n is not None:
-      ax.axhline(y,linewidth=0.5,alpha=0.5,zorder=0)
-    if k is not None:
-      ax.axvline(x,linewidth=0.5,alpha=0.5,zorder=0)
-
-  ax.set_xlabel('n',fontsize=16)
-  ax.set_ylabel('k',fontsize=16)
-
-  ax.set_xlim((np.min(nRange),np.max(nRange)))
-  ax.set_ylim((np.min(kRange),np.max(kRange)))
-  ax.tick_params(which='both',direction='in')
-
-  if axisOption == 0:
-    if max(kSolutionsToPlot) <= 0.5 or kMax <= 1:
-      ax.set_yscale('log')
-    else:
-      ax.set_yscale('linear')
-  elif axisOption == 1:
-    ax.set_xscale('linear')
-    ax.set_yscale('linear')
-  elif axisOption == 2:
-    ax.set_yscale('log')
-  elif axisOption == 3:
-    ax.set_xscale('log')
-  elif axisOption == 4:
-    ax.set_xscale('log')
-    ax.set_yscale('log')
-  else:
-    pass
-
-  _c = ax.get_children()
-  if Qback is None:
-    if incErrors:
-      # no Qback, with error bounds
-      graphElements = {'Qsca':_c[0],'Qabs':_c[1], # contours
-                       'QscaErrFill':_c[2],'QscaErrOutline1':_c[3],'QscaErrOutline2':_c[4],
-                       'QabsErrFill':_c[5],'QabsErrOutline1':_c[6],'QabsErrOutline2':_c[7],
-                       'SolMark':_c[8],'SolFill':_c[9], # the circly thingies at each solutions
-                       'CrosshairsH':_c[10:-10:2],'CrosshairsV':_c[11:-10:2], # solution crosshairs
-                       'LeftSpine':_c[-10],'RightSpine':_c[-9],'BottomSpine':_c[-8],'TopSpine':_c[-7], # spines
-                       'XAxis':_c[-6],'YAxis':_c[-5]} # the axes
-    else:
-      # no Qback, no error bounds
-      graphElements = {'Qsca':_c[0],'Qabs':_c[1], # contours
-                       'SolFill':_c[2],'SolMark':_c[3], # the circly thingies at each solutions
-                       'CrosshairsH':_c[4:-10:2],'CrosshairsV':_c[5:-10:2], # solution crosshairs
-                       'LeftSpine':_c[-10],'RightSpine':_c[-9],'BottomSpine':_c[-8],'TopSpine':_c[-7], # spines
-                       'XAxis':_c[-6],'YAxis':_c[-5]} # the axes
-
-  else:
-    if incErrors:
-      # with Qback and error bounds
-      graphElements = {'Qsca':_c[0],'Qabs':_c[1],'Qback':_c[2], # contours
-                       'QscaErrFill':_c[4],'QscaErrOutline1':_c[5],'QscaErrOutline2':_c[6],
-                       'QabsErrFill':_c[7],'QabsErrOutline1':_c[8],'QabsErrOutline2':_c[9],
-                       'SolMark':_c[10],'SolFill':_c[11], # the circly thingies at each solutions
-                       'CrosshairsH':_c[12:-10:2],'CrosshairsV':_c[13:-10:2], # solution crosshairs
-                       'LeftSpine':_c[-10],'RightSpine':_c[-9],'BottomSpine':_c[-8],'TopSpine':_c[-7], # spines
-                       'XAxis':_c[-6],'YAxis':_c[-5]} # the axes
-    else:
-      # with Qback, no error bounds
-      graphElements = {'Qsca':_c[0],'Qabs':_c[1],'Qback':_c[2], # contours
-                       'SolFill':_c[3],'SolMark':_c[4], # the circly thingies at each solution
-                       'CrosshairsH':_c[5:-10:2],'CrosshairsV':_c[6:-10:2], # solution crosshairs
-                       'LeftSpine':_c[-10],'RightSpine':_c[-9],'BottomSpine':_c[-8],'TopSpine':_c[-7], # spines
-                       'XAxis':_c[-6],'YAxis':_c[-5]} # the axes
-
-  return solutionSet,forwardCalculations,solutionErrors, fig, ax, graphElements
-
-inv = ContourIntersection(Qsca,Qabs,wavelength,diameter,n=1.5,k=0.5)
+inv = ps.ContourIntersection_SD(Bsca,Babs,w,dp,ndp,n=1.75)

build/lib/PyMieScatt/Inverse.py

@@ -74,7 +74,7 @@ def Inversion_SD(Bsca,Babs,wavelength,dp,ndp,nMin=1,nMax=3,kMin=0,kMax=1,scatter
 
 def ContourIntersection(Qsca,Qabs,wavelength,diameter,Qback=None,n=None,k=None,nMin=1,nMax=3,kMin=0.00001,kMax=1,gridPoints=100,interpolationFactor=2,maxError=0.005,fig=None,ax=None,axisOption=0):
 #  http://pymiescatt.readthedocs.io/en/latest/inverse.html#ContourIntersectio
-  if Qabs == 0.0 or Qabs[0] == 0.0:
+  if (type(Qabs) == np.float64 and Qabs == 0.0) or (type(Qabs) in [list, tuple, np.ndarray] and Qabs[0]==0):
     k = 0.0
   if k == 0.0:
     kMin = -0.1
@@ -313,7 +313,7 @@ def ContourIntersection(Qsca,Qabs,wavelength,diameter,Qback=None,n=None,k=None,n
 
 def ContourIntersection_SD(Bsca,Babs,wavelength,dp,ndp,n=None,k=None,nMin=1,nMax=3,kMin=0.00001,kMax=1,Bback=None,gridPoints=60,interpolationFactor=2,maxError=0.005,fig=None,ax=None,axisOption=0):
 #  http://pymiescatt.readthedocs.io/en/latest/inverse.html#ContourIntersection_SD
-  if Babs == 0.0 or Babs[0] == 0.0:
+  if (type(Babs) == np.float64 and Babs == 0.0) or (type(Babs) in [list, tuple, np.ndarray] and Babs[0]==0):
     k = 0.0
   if k == 0.0:
     kMin = -0.1

build/lib/PyMieScatt/_version.py

@@ -1 +1 @@
-__version__ = "1.6.0a"
+__version__ = "1.6.0b"