- #1
NielsW
- 2
- 0
I have the following paper "Collecting Performance of a LiDAR Telescope at Short Distances":
http://earsel.org/wp-content/uploads/2016/11/3-3_03_Ohm.pdf
I am supposed to calculate the efficiency of the LiDAR, as shown in Fig. 4 in the paper. However, my graphs do not look at all as they do in the paper. I calculate b, B, M and P and with this AL as shown in the paper. Then I calculate Omega.
Then I integrate, as shown in Eq. (7). I integrate of the implicit variable r from 0 to R= beta *z / 2. The python code I wrote is shown below.
Reasons, I think I am doing something wrong:
• The Graphs look different
• There is a removable discontinuity when z = z0
• A lot of the graph is undefined, even when the graph is shown, parts of the makeup is undefined
http://earsel.org/wp-content/uploads/2016/11/3-3_03_Ohm.pdf
I am supposed to calculate the efficiency of the LiDAR, as shown in Fig. 4 in the paper. However, my graphs do not look at all as they do in the paper. I calculate b, B, M and P and with this AL as shown in the paper. Then I calculate Omega.
Then I integrate, as shown in Eq. (7). I integrate of the implicit variable r from 0 to R= beta *z / 2. The python code I wrote is shown below.
Reasons, I think I am doing something wrong:
• The Graphs look different
• There is a removable discontinuity when z = z0
• A lot of the graph is undefined, even when the graph is shown, parts of the makeup is undefined
Python:
import math
import matplotlib.pyplot as plt
import matplotlib.patches as mpatches
import numpy as np
from scipy.integrate import quad
f = 1.2 # focal length
beta = 2.4e-3 # laser beam divergence
L=0.1 # lense radius
d= 1.25 # distance of diaphrang
sign = 1 # in the plus minus part, the sign defines wether plus or minus
def R(z,beta):
return z * beta /2
def z0 (d,f): # particular depth
if (d-f)==0:
return None
return (d*f)/(d-f)
def Dcalc(f, d) : # radius of diaphrang
return R(z0(d,f), beta) * d / z0(d,f)
def b(z, f): # image distance
if z-f==0:
return None
return (z * f ) / (z - f )
def B (z, r, f): # image size
if z-f==0:
return None
return (r * f ) / (z - f )
def M(z, r,f,d) : # center of projection of diaphrgm to the plan of the lens
if b(z,f) == None or (b(z,f) - d ) == 0:
return None
return (B(z,r,f) * d ) / (b(z,f) - d )
def P(z,D,f,d) : # radius of the projection on lens plane
if b(z,f) == None or (b(z,f) - d ) == 0:
return None
return (D * b(z,f) ) / (b(z,f) - d )
def AL(z,r,f,d,D,L,sign): # Area
if ( M(z,r,f,d)==None or ((M(z,r,f,d) * L) == 0) or ((M(z,r,f,d) * P(z,D,f,d)) == 0 ) or P(z,D,f,d) ==None ):
return None
inAa = 1.1*(( M(z,r,f,d)**2 + L**2 - P(z,D,f,d)**2) / (2 * M(z,r,f,d) * L ))
inAb = 1.1*(( M(z,r,f,d)**2 - L**2 + P(z,D,f,d)**2) / (2 * M(z,r,f,d) * P(z,D,f,d)))
if not((inAa >-1 and inAa < 1) and (inAb >-1 and inAb < 1)):
return None
teilA = (L**2 * math.acos( inAa) + P(z,D,f,d)**2 * math.acos( inAb))
if sign >0:
teilB = 0.5 * (( 4 * L**2 * M(z,r,f,d)**2 + ( M(z,r,f,d)**2 + L**2 -P(z,D,f,d)**2)**2) )**(0.5)
else:
teilB = 0.5 * (( 4 * L**2 * M(z,r,f,d)**2 - ( M(z,r,f,d)**2 + L**2 -P(z,D,f,d)**2)**2) )**(0.5)
return teilA - teilB
def Omega(z,r,f,d,D,L,sign):
if AL(z,r,f,d,D,L,sign) == None or z==0:
return None
return AL(z,r,f,d,D,L,sign)/ z**2
def IntegrationOfOmega(z,r,f,d,D,L,sign):
if Omega(z,r,f,d,D,L,sign) ==None:
return 0
return Omega(z,r,f,d,D,L,sign)*r
def Sensitivity (f,d,D,L,sign, beta):
sens = []
zValue = []
ooz = []
for iii in np.arange( 2, 60, 0.5):
integrand = lambda r: IntegrationOfOmega(iii,r,f,d,D,L,sign)
zValue.append(iii)
if iii==0:
ooz.append(None)
else:
ooz.append(1.82*L**2/(iii**2))
if (R(iii,beta)) == 0 or R(iii,beta) == None:
sens.append (None)
elif (d-((iii*f)/(iii-f))==0):
sens.append (None)
else:
result, error = quad(integrand, 0, R(iii,beta))
sens.append(2 / (R(iii,beta) ** 2) * result)
return zValue , sens , ooz
D = Dcalc(f, d)
#'''
Solution25p = Sensitivity (f,1.2501, D,L, sign, beta)
Solution23p = Sensitivity (f,1.2329, D,L, sign, beta)
Solution22p = Sensitivity (f,1.2245, D,L, sign, beta)
#'''
Solution25m = Sensitivity (f,1.2501, D,L,-sign, beta)
Solution23m = Sensitivity (f,1.2329, D,L,-sign, beta)
Solution22m = Sensitivity (f,1.2245, D,L,-sign, beta)
#'''
plt.plot( Solution25p[0], Solution25p[1],label='z0=30m, sign=+' )
plt.plot( Solution23p[0], Solution23p[1],label='z0=45m, sign=+' )
plt.plot( Solution22p[0], Solution22p[1],label='z0=60m, sign=+')
#'''
plt.plot( Solution25m[0], Solution25m[1],label='z0=30m, sign=-' )
plt.plot( Solution23m[0], Solution23m[1],label='z0=30m, sign=-' )
plt.plot( Solution22m[0], Solution22m[1],label='z0=30m, sign=-' )
#'''
plt.legend()
plt.xlabel("depth z/m")
plt.ylabel("Sensitivity")
plt.show()
Last edited by a moderator: