f function of Sphere class in profile.py

[Pinnacleman]

Hi,

I'm trying to add a new type of surface profile to rayoptics. I know I need to implement the sag function and the surface derivative function for the new surface type I want to add. So I'm studying the profile.py module and particularly the Sphere class since it is the simplest surface profile. I'm having hard time to understand the f function in Sphere class, which I re-posted below

def f(self, p):
return p[2] - 0.5self.cv(np.dot(p, p))

This f function should be Eq. (16) of Spencer's paper, which is basically the sag function of a sphere . But I cannot derive it into the implemented form. Is there something wrong with my understanding?

Thanks,

John

[mjhoptics]

Hello @Pinnacleman,
As it turns out, the f() implementation in Spherical is not used in the ray trace; the intersect() method is overridden directly in Spherical. I vaguely recall some issue that might have been related to what you've found; implementing intersect() made those issues moot (although I should either fix it or remove it to prevent others from befalling your fate).
As far as your original question, I would look at the implementations of the aspheric part of RadialPolynomial. The class uses the sag() function to implement f() but this is not required. f() maybe implemented directly if that makes more sense for your profile.
Feel free to ask further questions.
Sorry for the delay in responding.
Mike Hayford

[mjhoptics]

Hello @Pinnacleman,

I've looked into the code more closely and I believe it is correct. It is an implementation of Spencer's eq 25, for conic surfaces. The virtue of this form of the equation is that a sqrt() isn't required. Back in the day, specializing to avoid a sqrt() delivered significant speed gains in the ray trace. Spherical is by far the most frequently used profile, and so paying attention to algorithm efficiency has some benefit here.

Starting with the code:

F = p[2] - 0.5*self.cv*(np.dot(p, p))

Expanded out, the expression is:

F = p[2] - 0.5*self.cv*(p[0]*p[0] + p[1]*p[1] + p[2]*p[2])

Spencer's notation corresponds to:

rho**2 = p[0]*p[0] + p[1]*p[1]
Z = p[2]

Substituting notations, the code becomes

F = Z - 0.5*cv*(rho**2 + Z**2)

which is Spencer's eq 25 for kappa=1.

I've added comments to f() in Spherical and Conic to document this.

Regards,

Mike Hayford

[Pinnacleman]

Hi, Mike:

Thanks for your explanation. As you pointed out, I avoided this problem by overriding the intersect() function in my new surface type class.

John