00002 REM BATHOS 00004 REM OPTICAL DESIGN PROGRAMME * SEPTEMBER 1978 00006 REM BRIAN BLANDFORD, 25 NEW QUEBEC ST, LONDON W.1. 00010 DIM P$(11),P(11,7),V(3,10),Q(4,12),X(13,10),F(13) 00012 DATA C,D,N,V,H,U,L,HBAR,UBAR,LBAR,S,SBAR 00014 FOR J1=0 TO 11 00016 READ P$(J1) 00020 NEXT J1 00022 READ N$,K9,P9,P(2,0),P(3,0) 00024 FOR J9=1 TO K9 00026 FOR I9=0 TO 3 00028 READ P(I9,J9) 00030 NEXT I9 00032 NEXT J9 00034 READ P(4,1),P(5,0),P(7,1),P(8,0) 00035 REM LENS DATA LINES 35 TO 50 00036 REM FORMAT: 00037 REM "LENS TITLE",NO OF SURFACES, WAVELENGTH 00038 REM INDEX, V-VALUE IN OBJECT SPACE 00039 REM FOR EACH SURFACE: 00040 REM CURVATURE, SEPN, INDEX, V IN SPACE FOLLOWING 00041 REM PARAXIAL MARGINAL H,1 U,0 AND PUPIL HBAR,1 UBAR,0 00042 REM WELFORD M=U SIGN CONVENTION. 00043 REM EXAMPLE: 00044 DATA "DOUBLET, CONRADY P.316",3,.0000231,1,1E9 00045 DATA -0.18,0.1,1.622368,36.03 00046 DATA 0.067,0.3,1.526419,60.0 00047 DATA -0.399,0,1,1E9 00048 DATA 0.3,0,0,0.32492 00049 DATA 00050 DATA 00051 O0=0.7231584E-7 00052 GOSUB 00244 00053 PRINT 00054 PRINT "WOULD YOU LIKE SYSTEM DATA PRINTOUT:"; 00056 INPUT R$ 00058 IF R$="NO" THEN 00061 00060 GOSUB 00226 00061 PRINT 00062 PRINT"WOULD YOU LIKE PARAXIAL AND SEIDELS DATA:"; 00064 INPUT R$ 00066 IF R$="NO" THEN 00069 00068 GOSUB 00386 00069 PRINT 00070 PRINT "DO YOU WISH TO DEFINE A NEW LIST OF VARIABLES:"; 00072 INPUT R$ 00074 IF R$<>"NO" THEN 00080 00076 PRINT "WE SHALL USE THE SAME VARIABLES AS BEFORE" 00078 GOTO 00126 00080 PRINT "YOU HAVE A CHOICE OF" 00081 PRINT " C OF ANY SURFACE" 00082 PRINT " D,N,V IN SPACE FOLLOWING IT" 00083 PRINT " H,1 U,0 HBAR,1 UBAR,0" 00084 PRINT "TYPE: PARAMETER,SURFACE NO. (E.G. D,7)" 00085 PRINT " END,99 TO FINISH)" 00087 PRINT "DEFAULT INCREMENT IS 1 PERCENT" 00088 V0=1 00090 PRINT "VARIABLE";V0;":="; 00092 INPUT R$,R0 00094 IF R0=99 THEN 00124 00096 GOSUB 00412 00098 IF P0=-1 THEN 00090 00100 V(0,V0)=P(P0,R0) 00102 V(1,V0)=P0 00104 V(2,V0)=R0 00106 PRINT "CURRENT VALUE IS";V(0,V0) 00108 PRINT "INCREMENT="; 00110 INPUT R0 00112 IF R0<>0 THEN 00118 00114 R0=O0+0.01*V(0,V0) 00116 PRINT "DEFAULT INCREMENT IS";R0 00118 V(3,V0)=R0 00120 V0=V0+1 00122 GO TO 00090 00120 V0=V0-1 00126 PRINT "THERE ARE ";V0;" VARIABLES IN ALL" 00127 PRINT 00128 PRINT "DO YOU WISH TO DEFINE THE MERIT FUNCTION:"; 00130 INPUT R$ 00132 IF R$<>"NO" THEN 00138 00134 PRINT "WE SHALL USE THE SAME ONE AS BEFORE" 00136 GOTO 00204 00138 PRINT "THE MERIT FUNCTION IS THE SUM OF" 00140 PRINT "((ABERRATIONS-TARGETS)*WEIGHTS)SQUARED," 00142 PRINT "WHERE ABERRATIONS ARE CHOSEN FROM:" 00144 PRINT " H,HBAR AT ANY SURFACE" 00146 PRINT " U,UBAR,L,LBAR FOLLOWING IT" 00148 PRINT " S,1...S,7 SEIDELS SI...SV CL CT" 00150 PRINT " SBAR,1...SBAR,7 PUPIL SEIDELS" 00152 PRINT " C,D,N,V AS BEFORE" 00154 PRINT "TYPE: PARAMETER,INDEX (E.G. H,8)" 00156 PRINT " (END,99 TO FINISH)" 00158 Q0=1 00160 F0=0 00162 PRINT "ABERRATION";Q0;":="; 00164 INPUT R$,R0 00166 IF R0=99 THEN 00200 00168 GOSUB 00412 00170 IF P0=-1 THEN 00162 00172 Q(0,Q0)=P(P0,R0) 00174 Q(1,Q0)=P0 00176 Q(2,Q0)=R0 00178 PRINT "CURRENT VALUE IS";Q(0,Q0) 00180 PRINT "TARGET="; 00182 INPUT Q(3,Q0) 00184 PRINT "WEIGHT="; 00186 INPUT Q(4,Q0) 00188 F1=(Q(0,Q0)-Q(3,Q0))*Q(4,Q0) 00190 F1=F1*F1 00192 PRINT "MERIT FUNCTION CONTRIBUTION IS";F1 00194 F0=F0+F1 00196 Q0=Q0+1 00198 GOTO 00162 00200 PRINT "TOTAL MERIT FUNCTION IS";F0 00202 Q0=Q0-1 00204 PRINT "THERE ARE";Q0;"ABERRATIONS IN ALL." 00206 PRINT 00208 PRINT "SIMPLEX OPTIMISATION" 00210 PRINT "HOW MANY STEPS="; 00212 INPUT K7 00214 GOSUB 600 00216 PRINT "FINAL VALUE OF MERIT FUNCTION IS";F6 002108 GOTO 00053 00224 REM ** SUBROUTINES ** 00225 REM 00226 REM SYSTEM DATA PRINTOUT 00228 PRINT N$,K9;" SURFACES" 00230 PRINT " CURV SEPN INDEX V" 00232 PRINT " ";P(1,0);P(2,0);P(3,0) 00234 FOR I1=1 TO K9 00236 PRINT I1;P(0,I1) 00238 PRINT " ";P(1,I1);P(2,I1);P(3,I1) 00240 NEXT I1 00242 RETURN 00244 REM CALCULATION OF ALL PARAMETERS 00246 H9=P(4,1) 00248 U9=P(5,0) 00250 N9=P(2,0) 00252 G9=P(7,1) 00254 T9=P(8,0) 00256 V9=P(3,0) 00258 O8=-N9*(G9*U9-H9*T9) 00260 D9=P(1,0) 00262 O9=1.7231584E-20 00264 S8=W8=X8=Y8=Z8=E8=F8=0 00266 S9=W9=X9=Y9=Z9=E9=F9=0 00268 FOR I8=1 TO K9 00270 C8=P(0,I8) 00272 D8=D9 00274 D9=P(1,I8) 00276 N8=N9 00278 N9=P(2,I8) 00280 V8=V9 00282 V9=P(3,I8)+O9 00284 H8=H9 00286 U8=U9 00288 G8=G9 00290 T8=T9 00292 A8=N8*(H8*C8+U8) 00294 B8=N8*(G8*C8+T8) 00296 U9=A8/N9-H8*C8 00298 T9=B8/N9-G8*C8 00300 R8=U9/N9-U8/N8 00302 Q8=T9/N9-T8/N8 00304 P8=(1-1/N9)/V9-(1-1/N8)/V8 00306 S8=S8-A8*A8*H8*R8 00308 W8=W8-A8*B8*H8*R8 00310 X8=X8-B8*B8*H8*R8 00312 Y8=Y8-O8*(A8*Q8-B8*R8) 00314 Z8=Z8-O8*B8*Q8-B8*B8*G8*R8 00316 E8=E8+A8*H8*P8 00318 F8=F8+B8*H8*P8 00320 S9=S9-B8*B8*G8*Q8 00322 W9=W9-B8*A8*G8*Q8 00324 X9=X9-A8*A8*G8*Q8 00326 Z9=Z9+O8*A8*R8-A8*A8*H8*Q8 00328 E9=E9+B8*G8*P8 00330 F9=F9+A8*G8*P8 00332 L8=-H8/(U9+O9) 00334 M8=-G8/(T9+O9) 00336 H9=H8+D9*U9 00338 G9=G8+D9*T9 00340 P(4,I8+1)=H9 00342 P(5,I8)=U9 00344 P(6,I8)=L8 00348 P(7,I8+1)=G9 00350 P(8,I8)=T9 00352 P(9,I8)=M8 00354 NEXT I8 00356 P(10,1)=S8 00358 P(10,2)=W8 00360 P(10,3)=X8 00362 P(10,4)=Y8 00364 P(10,5)=Z8 00366 P(10,6)=E8 00368 P(10,7)=F8 00370 P(11,1)=S9 00372 P(11,2)=W9 00374 P(11,3)=X9 00376 P(11,4)=Y8 00378 P(11,5)=Z9 00380 P(11,6)=E9 00382 P(11,7)=F9 00384 RETURN 00386 REM PARAXIAL TRACE AND SEIDELS OUTPUT 00388 PRINT " H U HBAR UBAR" 00390 PRINT " ";P(5,0);" ";P(8,0) 00392 FOR I1=1 TO K9 00394 PRINT P(4,I1);P(5,I1);P(7,I1);P(8,I1) 00396 NEXT I1 00398 PRINT "SEIDELS" 00400 PRINT P(10,1);P(10,2);P(10,3);P(10,4);P(10,5) 00402 PRINT P(10,6);P(10,7) 00404 PRINT "PUPIL SEIDELS" 00406 PRINT P(11,1);P(11,2);P(11,3);P(11,4);P(11,5) 00408 PRINT P(11,6);P(11,7) 00410 RETURN 00412 REM PARAMETER RECOGNITION 00414 P0=-1 00416 FOR J1=0 TO 11 00418 IF R$=P$(J1) THEN 00426 00420 NEXT J1 00422 PRINT "WHAT?" 00424 RETURN 00426 P0=J1 00428 RETURN 600 REM SIMPLEX OPTIMISATION 620 N6=V0 622 M6=N6+1 624 E6=0.0001 626 FOR I6=1 TO M6 628 FOR J6=1 TO N6 630 X(I6,J6)=V(0,J6) 631 NEXT J6 632 NEXT I6 633 FOR I6=2 TO M6 634 X(I6,I6-1)=X(I6,I6-1)+V(3,I6-1) 636 NEXT I6 656 K6=0 658 REM EVALUATE MERIT FUNCTIONS 660 FOR I6=1 TO M6 668 GO SUB 1000 670 F(I6)=F6 672 NEXT I6 690 REM SELECT L,G,H 692 H6=1 694 L6=1 696 H7=F(1) 698 L7=H7 700 FOR I6=2 TO M6 702 IF F(I6)L7 THEN 714 710 L6=I6 712 L7=F(I6) 714 NEXT I6 716 G6=L6 718 G7=L7 720 FOR I6=1 TO M6 722 IF I6=H6 THEN 730 724 IF F(I6)L7+E6 THEN 748 746 PRINT "CONVERGED" 747 GOTO 950 748 K6=K6+1 750 REM PRINTING 752 PRINT " MERIT =";L7 780 REM FIND REFLECTION AND EXTENDED STEP 782 FOR J6=1 TO N6 784 S6=O 786 FOR I6=1TO M6 788 S6=S6+X(I6,J6) 790 NEXT I6 792 C6 =(S6-X(H6,J6))/N6 794 X(M6+1,J6)=2*C6-X(H6,J6) 796 X(M6+2,J6)=3*C6-2*X(H6,J6) 798 NEXT J6 800 REM EVALUATE FR 802 I6=M6+1 808 GOSUB 1000 810 R7=F6 812 IF R7>(L7+G7)/2 THEN 850 813 REM TRY EXTENDED STEP 814 I6=M6+2 822 GOSUB 1000 824 N7=F6 826 IF N7>R7 THEN 853 827 PRINT "*E"; 828 REM ADOPT EXTENDED OR CONTRACTED STEP 830 FOR J6=1 TO N6 832 X(H6,J6)=X(M6+2,J6) 834 NEXT J6 836 F(H6)=N7 838 GO TO 692 850 IF R7>G7 THEN 866 852 REM ADOPT NORMAL STEP 853 PRINT "*R"; 854 FOR J6=1 TO N6 856 X(H6,J6)=X(M6+1,J6) 858 NEXT J6 860 F(H6)=R7 862 GO TO 692 864 REM CONTRACTED STEP 866 IF R7>H7 THEN 879 868 REM CONTRACTION FROM R 869 PRINT "*RC"; 872 FOR J6=1 TO N6 874 X(M6+2,J6)=0.25*X(H6,J6)+0.75*X(M6+1,J6) 876 NEXT J6 877 GO TO 888 878 REM CONTRACT FROM H 879 PRINT "*HC"; 882 FOR J6=1 TO N6 884 X(M6+2,J6)=0.75*X(H6,J6)+0.25*X(M6+1,J6) 886 NEXT J6 888 REM EVALUATE FN 890 I6=M6+2 896 GOSUB 1000 898 N7=F6 900 IF N7H7 THEN 916 906 REM ADOPT R 908 FOR J6=1 TO N6 910 X(H6,J6)=X(M6+1,J6) 912 NEXT J6 914 H7=R7 916 IF N7>H7 THEN 928 918 REM ADOPT N 920 FOR J6=1 TO N6 922 X(H6,J6)=X(M6+2,J6) 924 NEXT J6 926 REM SHRINK WHOLE SIMPLEX 928 FOR I6=1 TO M6 932 FOR J6=1 TO N6 934 X(I6,J6)=0.5*(X(I6,J6)+X(L6,J6)) 936 NEXT J6 938 NEXT I6 940 GO TO 660 950 REM FINAL VARIABLES AND ABERRATIONS 952 I6=L6 954 FOR J6=1 TO N6 956 V(0,J6)=X(I6,J6) 958 NEXT J6 960 GOSUB 1000 962 FOR J6=1 TO Q0 964 Q(0,J6)=P(Q(1,J6),Q(2,J6)) 966 NEXT J6 968 RETURN 1000 REM MERIT FUNCTION EVALUATION 1002 FOR J6=1 TO V0 1004 P(V(1,J6),V(2,J6))=X(I6,J6) 1006 NEXT J6 1008 GOSUB 00244 1010 F6=0 1012 FOR J6=1 TO Q0 1014 F7=(P(Q(1,J6),Q(2,J6))-Q(3,J6))*Q(4,J6) 1016 F6=F6+F7*F7 1018 NEXT J6 1020 RETURN 9999 END