//reflecting spinning mirror calculator
num_lasers=2;
laser_angles=[linspace(0,360*((num_lasers-1)/num_lasers),num_lasers);zeros(1,num_lasers)]/180*%pi; //[rad] theta;phi of laser such that (0;0)=+1 in X, theta rotates towards Y+ in XY plane, phi rotates towards Z+ out of XY plane
dt=laser_angles(1,2)-laser_angles(1,1);
//la2=laser_angles;
//la2(1,:)=la2(1,:);//+dt/2;
//la2(2,:)=-ones(1,num_lasers)*dt;
//laser_angles=[laser_angles,la2];
//num_lasers=num_lasers*2;
//laser_angles=[linspace(0,11*%pi,num_lasers);linspace(-80,-2,num_lasers)/180*%pi];
//laser_angles=[laser_angles,[laser_angles(1,:)+%pi;laser_angles(2,:)]];
//laser_angles=[linspace(0,360*((num_lasers-1)/num_lasers),num_lasers);zeros(1,num_lasers)]/180*%pi;
//laser_angles=[0;0]/180*%pi;
//num_lasers=size(laser_angles,2);
//dxl=laser_angles(2,2)-laser_angles(2,1);
//dxt=(laser_angles(2,1)+%pi/2)/2;
npl=360;
//mirror_angles=[[linspace(0,360,npl),linspace(0,360,npl),linspace(0,360,npl)]/180*%pi+%pi*1.5;dxt-[0*ones(1,npl), dxl/6*ones(1,npl), dxl*2/6*ones(1,npl)]];
mirror_angles=[(180+linspace(0,360,npl))/180*%pi;45*ones(1,npl)/180*%pi];
laser_acceptance_angle_h=88/180*%pi; //+-[rad]
clf();

function rv=reflect(sv,mv)
    //accepts multiple mv, single sv
    //mv is normal of mirror surface
    //sv is pointing from mirror to laser location
    //rv is pointing from mirror to reflected spot location
    //all must be unit length 3D column vectors
    dv=[mv(1,:)-sv(1);mv(2,:)-sv(2);mv(3,:)-sv(3)];
    dp=sum(dv.*mv,1);
    a1=[mv(1,:).*dp;mv(2,:).*dp;mv(3,:).*dp];
    a2=2.*(dv-a1);
    rv=[sv(1)+a2(1,:);sv(2)+a2(2,:);sv(3)+a2(3,:)];
endfunction

function r=ierem(n,d)
    r=n-d*round(n/d); //IEEERemainder() returns +-distance to nearest multiple
endfunction

function [q,w]=contregion(b)
    //returns vector [a1,a2;c1,c2;...] where each row is start and end index (inclusive) of continuous region of T values in b
    //and boolean w is T if the first and last rows should be combined into a single segment
    a=length(b);
    r=b(1);
    r1=1;
    q=[];
    for i=2:a
        if r==b(i)
            continue;
        end
        if r
            r2=i-1;
            q=[q;r1,r2];
        else
            r1=i;
        end
        r=b(i);
    end
    if r then
        r2=a;
        q=[q;r1,r2];
    end
    w=b(1)&b(a);
endfunction

function plottracelines(rvs,lv,cl)
    [q,w]=contregion(lv);
    if(isempty(q))
        return;
    end
    rv=rvs(:,q(1,1):q(1,2));
    nq=size(q,1);
    ni=1;
    while ni<=nq
        param3d(rv(1,:),rv(2,:),rv(3,:));
        gce().foreground = color(cl);
        ni=ni+1;
        if ni<=nq
            rv=rvs(:,q(ni,1):q(ni,2));
        end
    end
endfunction

mv=[cos(mirror_angles(2,:)).*cos(mirror_angles(1,:));cos(mirror_angles(2,:)).*sin(mirror_angles(1,:));sin(mirror_angles(2,:))];

for il=1:num_lasers//[1,2,num_lasers/2,num_lasers/2+1] //1:num_lasers
    tls=mirror_angles(1,:)-laser_angles(1,il);
    tls=ierem(tls,2*%pi);
    tl1=(abs(tls)<=laser_acceptance_angle_h);
    sv=[cos(laser_angles(1,il)).*cos(laser_angles(2,il));sin(laser_angles(1,il)).*cos(laser_angles(2,il));sin(laser_angles(2,il))];
    rv=reflect(sv,mv);
    plottracelines(rv,tl1,'red');
    param3d(sv(1),sv(2),sv(3));
    gce().mark_mode = "on";
    gce().mark_style=9;
end
if(%t)
param3d(0,0,1);
gce().mark_mode = "on";
gce().mark_style=3;
end
isoview on
