31:T5ab,function getScatter(size, n) { const {r1,r2} = params; const step = 2 * size / n; const xy = d3.range(n).map(i => step * (i + 0.5) - size); const grid = [ ]; for(const y1 of xy) { for(const x1 of xy) { for(const y2 of xy) { for(const x2 of xy) { if((x1 == x2) && (y1 == y2)) continue; // Calculate initial angles from FVC measurement (x1,y1) assuming correct calib const [T1,P1] = xy2tp(x1, y1, r1, r2, 0., 0.); if(Number.isNaN(T1) || Number.isNaN(P1)) continue; // Calculate target angles assuming correct calib const [T2req,P2req] = xy2tp(x2, y2, r1, r2, 0., 0.); if(Number.isNaN(T2req) || Number.isNaN(P2req)) continue; // Calculate actual angles achieved with these target angles const pos = getPosition(T2req*180/Math.PI, P2req*180/Math.PI); const [x2act,y2act,z2act] = pos.P1; // Calculate final angles from FVC measurement (x2act,y2act) assuming correct calib const [T2act,P2act] = xy2tp(x2act,y2act, r1, r2, 0., 0.); // Calculate requested vs actual changes in angle const dTreq = 180/Math.PI*(T2req - T1), dPreq = 180/Math.PI*(P2req - P1); let dTact = 180/Math.PI*(T2act - T1), dPact = 180/Math.PI*(P2act - P1); if(dTact < dTreq - 350) dTact += 360; grid.push({x1,y1,x2,y2,dTreq,dPreq,dTact,dPact}); } } } } return grid; }32:T52a,function getGrid(size, n) { const {r1,r2} = params; const step = 2 * size / n; const xy = d3.range(n).map(i => step * (i + 0.5) - size); const grid = [ ]; let dxSum = 0, dySum = 0, nSum = 0; for(const y1 of xy) { for(const x1 of xy) { // Convert target (x1,y1) to angles (T,P) assuming calibration arm lengths (r1,r2) are correct and no 3D misalignments const [T, P] = xy2tp(x1, y1, r1, r2, 0., 0.); if(!Number.isNaN(T) && !Number.isNaN(P)) { //const [x2,y2] = tp2xy(T, P, r1 + 1e-3*dr1, r2 + 1e-3*dr2, 1e-3*dx0, 1e-3*dy0); // Calculate true 3D position of fiber tip for angles (T,P) given arm-length errors (dr1,dr2) and any 3D misalignments. const pos = getPosition(T*180/Math.PI, P*180/Math.PI); const [x2,y2,z2] = pos.P1; // Calculate resulting 2D errors (dx,dy) when moving to target (x1,x2) const dx = x2 - x1, dy= y2 - y1; grid.push({x1, y1, T, P, dx, dy}); if(Math.hypot(x1, y1) <= 6) { dxSum += dx; dySum += dy; nSum++; } } } } if(subtractMean) { // subtract mean residual over patrol radius r < 6mm dxSum /= nSum; dySum /= nSum; for(let i=0; i < grid.length; ++i) { grid[i].dx -= dxSum; grid[i].dy -= dySum; } } return grid; }33:T40d,function zPlot(taxis) { const { maxR, maxH } = plotGeom; return Plot.plot({ height: 0.15 * width * maxH / (2 * maxR), width: 0.15 * width, x: { domain:[-maxR,maxR] }, y: { domain:[0,maxH] }, marks: [ Plot.frame(), Plot.line([ [0,0,0], posTP.T0 ], { x:d=>d[taxis], y:d=>d[2], stroke:"black", strokeDasharray: [2,4] }), Plot.line(posGridT, { x:d=>d.T1[taxis], y:d=>d.T1[2], stroke:"blue" }), Plot.line([ posTP.T0, posTP.T1 ], { x:d=>d[taxis], y:d=>d[2], stroke:"blue", strokeWidth: 8, strokeOpacity:0.4 }), Plot.line([ posTP.T1, posTP.P0 ], { x:d=>d[taxis], y:d=>d[2], stroke:"black", strokeDasharray: [2,4] }), Plot.line(posGridP, { x:d=>d.P1[taxis], y:d=>d.P1[2], stroke:"red" }), Plot.line([ posTP.P0, posTP.P1 ], { x:d=>d[taxis], y:d=>d[2], stroke:"red", strokeWidth: 8, strokeOpacity:0.4 }), ] }); }34:T4c6,function getPosition(theta, phi) { // Convert angles from deg to rad theta = theta * Math.PI / 180; phi = phi * Math.PI / 180; const {HT,RT,HP,RP,axisT,axisP} = geometry; // Normalize axis vectors const nT = normalize(axisT); const nP = normalize(axisP); // Offset in frame of P body let P0 = [ 0, 0, HP ], P1 = [ RP, 0, HP ]; // Misalign P body const rotP = buildRotationMatrix(nP); P0 = applyRotation(rotP, P0); P1 = applyRotation(rotP, P1); // Rotate about misaligned P axis (should be no-op for P0) P0 = rodriguesRotate(P0, nP, phi); P1 = rodriguesRotate(P1, nP, phi); // Offset in frame of T body let T0 = [ 0, 0, HT ], T1 = [ RT, 0, HT ]; P0 = offset(P0, T1); P1 = offset(P1, T1); // Apply T body misalignment to T+P assembly const rotT = buildRotationMatrix(nT); T0 = applyRotation(rotT, T0); T1 = applyRotation(rotT, T1); P0 = applyRotation(rotT, P0); P1 = applyRotation(rotT, P1); // Rotate T+P body about misaligned T axis (should be no-op for T0) T0 = rodriguesRotate(T0, nT, theta); T1 = rodriguesRotate(T1, nT, theta); P0 = rodriguesRotate(P0, nT, theta); P1 = rodriguesRotate(P1, nT, theta); return { T0, T1, P0, P1 }; }35:[] 2f:{"title":"Robot Calibration Errors","cells":[{"id":0,"value":"# Robot Calibration Errors\n\nUse the sliders below to introduce errors into the [DESI robot](https://observablehq.com/@dkirkby/desi-robots) calibration (r1,r2,x0,y0) and see the resulting residuals on a grid covering the robot patrol radius:","pinned":false,"hidden":false,"mode":"md"},{"id":214,"value":"viewof ref = Inputs.radio([10, 25, 50, 75, 100], {label: \"Scale [µm]\", value: 25})","pinned":false,"hidden":false,"mode":"ojs"},{"id":451,"value":"viewof subtractMean = Inputs.toggle({label: \"Subtract mean\", value: true})","pinned":false,"hidden":false,"mode":"ojs"},{"id":488,"value":"viewof angles = Inputs.form(\n {\n phi: Inputs.range([0, 180], {label: \"Phi\", value: 0, step: 1}),\n theta: Inputs.range([0, 360], {label: \"Theta\", value: 0, step: 1})\n },\n {\n template: (inputs) => html`

${Object.values(inputs)}`\n }\n)","pinned":false,"hidden":false,"mode":"ojs"},{"id":504,"value":"distPlots = html`\n

\n ${xyPlot(0.5*width)} ${distortionPlot()}\n

`","pinned":false,"hidden":false,"mode":"ojs"},{"id":561,"value":"deltaPlots = html`\n

\n ${dPplot} ${dTplot}\n

`","pinned":false,"hidden":false,"mode":"ojs"},{"id":419,"value":"viewof params = Inputs.form(\n {\n r1: Inputs.range([2.8,3.2], {label:\"r1 (mm)\", step: 0.01, value:3}),\n r2: Inputs.range([2.8,3.2], {label:\"r2 (mm)\", step: 0.01, value:3}),\n dr1: Inputs.range([-20,20], {label:\"dr1 (μm)\", step: 0.1, value:0}),\n dr2: Inputs.range([-20,20], {label:\"dr2 (μm)\", step: 0.1, value:0}),\n Pmis: Inputs.range([0,5], {label: \"Pmis (mrad)\", step: 0.02, value: 0}),\n Paz: Inputs.range([0, 360], {label: \"Paz (deg)\", step: 1, value: 0}),\n Tmis: Inputs.range([0,5], {label: \"Tmis (mrad)\", step: 0.1, value: 0}),\n Taz: Inputs.range([0, 360], {label: \"Taz (deg)\", step: 1, value: 0}),\n },\n {\n template: (inputs) => html`

${Object.values(inputs)}`\n }\n)","pinned":false,"hidden":false,"mode":"ojs"},{"id":266,"value":"plots3D = html`\n

\n ${xyPlot()} ${zPlot(0)} ${zPlot(1)}\n

`","pinned":false,"hidden":false,"mode":"ojs"},{"id":26,"value":"function distortionPlot() {\n const arrowLength = 1.0; // mm\n const zoom = arrowLength / (1e-3 * ref);\n return Plot.plot({\n width: 0.7 * width,\n height: 0.7 * width,\n x: {domain:[-6,6], ticks:[] },\n y: {domain:[-6,6], ticks:[] },\n marks: [\n Plot.arrow(grid, {\n x1: \"x1\",\n y1: \"y1\",\n x2: d => d.x1 + zoom * d.dx,\n y2: d => d.y1 + zoom * d.dy,\n }),\n Plot.arrow([0], { x1: 5 - 0.5 * arrowLength, y1: 5.5, x2: 5 + 0.5 * arrowLength, y2: 5.5 }),\n Plot.text([0], { x: 5.0, y: 5.5, text: d=>`${ref}μm`, dy:12, fontSize:14 }),\n ]\n });\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":15,"value":"grid = getGrid(6, 19)","pinned":false,"hidden":false,"mode":"ojs"},{"id":520,"value":"$31","pinned":false,"hidden":false,"mode":"ojs"},{"id":523,"value":"scatter = getScatter(6, 11);","pinned":false,"hidden":false,"mode":"ojs"},{"id":555,"value":"dPplot = {\n return Plot.plot({\n x: { label:\"Requested dT [deg]\" },\n y: { label:\"Actual dT [deg]\" },\n grid: true,\n marks: [\n Plot.line([{x:-180,y:-180},{x:180,y:180}], {x:\"x\", y:\"y\", stroke:\"blue\"}),\n Plot.dot(scatter, {x:\"dPreq\", y:\"dPact\", r:1}),\n ]\n });\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":546,"value":"dTplot = {\n return Plot.plot({\n x: { label:\"Requested dP [deg]\" },\n y: { label:\"Actual dP [deg]\" },\n grid: true,\n marks: [\n Plot.line([{x:-360,y:-360},{x:360,y:360}], {x:\"x\", y:\"y\", stroke:\"blue\"}),\n Plot.dot(scatter, {x:\"dTreq\", y:\"dTact\", r:1}),\n ]\n });\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":13,"value":"$32","pinned":false,"hidden":false,"mode":"ojs"},{"id":107,"value":"Convert fiber position (x,y) to robot angles (T,P) using the calibration (r1,r2,x0,y0) and assuming no 3D misalignments:","pinned":false,"hidden":false,"mode":"md"},{"id":6,"value":"function xy2tp(x, y, r1, r2, x0, y0) {\n const {sqrt, max, atan2, acos, PI} = Math;\n x -= x0;\n y -= y0;\n if(x == 0 && y == 0) return (r1 == r2) ? [0., PI] : [NaN, NaN]; \n const r12 = r1 + r2, r1sq = r1 * r1, r2sq = r2 * r2;\n const rsq = x * x + y * y, r = sqrt(rsq);\n const phi = atan2(y, x);\n const s = (rsq - r2sq + r1sq) / (2 * r);\n const beta = acos(s / r1);\n //const a = sqrt(r1sq - s * s);\n //const beta = atan2(a, s);\n const T = phi - beta;\n const P = acos((rsq - r1sq - r2sq) / (2 * r1 * r2));\n return [T, P];\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":109,"value":"Convert robot angles (T,P) to fiber position (x,y) using the calibration (r1,r2,x0,y0) and assuming no 3D misalignments:","pinned":false,"hidden":false,"mode":"md"},{"id":37,"value":"function tp2xy(T, P, r1, r2, x0, y0) {\n const {cos, sin} = Math;\n const x = x0 + r1 * cos(T) + r2 * cos(T + P);\n const y = y0 + r1 * sin(T) + r2 * sin(T + P);\n return [x, y];\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":227,"value":"### 3D Misalignments","pinned":false,"hidden":false,"mode":"md"},{"id":251,"value":"function xyPlot(w=0.7*width) {\n const { maxR, maxH } = plotGeom;\n return Plot.plot({\n height: w,\n width: w,\n x: { domain:[-maxR,maxR] },\n y: { domain:[-maxR,maxR] },\n marks: [\n Plot.frame(),\n Plot.line(posGridT, { x:d=>d.T1[0], y:d=>d.T1[1], stroke:\"blue\", strokeDasharray: [2,4] }),\n Plot.line(posGridT, { x:d=>d.P1[0], y:d=>d.P1[1], stroke:\"blue\" }),\n Plot.line([ posTP.T0, posTP.T1 ], { x:d=>d[0], y:d=>d[1], stroke:\"blue\",\n strokeWidth: 8, strokeOpacity:0.4 }),\n Plot.line(posGridP, { x:d=>d.P1[0], y:d=>d.P1[1], stroke:\"red\" }),\n Plot.line([ posTP.P0, posTP.P1 ], { x:d=>d[0], y:d=>d[1], stroke:\"red\",\n strokeWidth: 8, strokeOpacity:0.4 }),\n ]\n });\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":256,"value":"$33","pinned":false,"hidden":false,"mode":"ojs"},{"id":366,"value":"plotGeom = {\n const maxR = d3.max(posGridP, d => Math.hypot(d.P1[0], d.P1[1])) + 0.1;\n const maxH = d3.max(posGridP, d => d.P1[2]) + 0.1;\n return { maxR, maxH };\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":345,"value":"posGridP = phiGrid.map(phi => getPosition(angles.theta, phi))","pinned":false,"hidden":false,"mode":"ojs"},{"id":234,"value":"posGridT = thetaGrid.map(theta => getPosition(theta, angles.phi))","pinned":false,"hidden":false,"mode":"ojs"},{"id":349,"value":"posTP = getPosition(angles.theta, angles.phi)","pinned":false,"hidden":false,"mode":"ojs"},{"id":338,"value":"phiGrid = {\n const n = 91;\n return d3.range(n).map(i => i * 180 / (n-1));\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":239,"value":"thetaGrid = {\n const n = 91;\n return d3.range(n).map(i => i * 360 / (n-1));\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":272,"value":"$34","pinned":false,"hidden":false,"mode":"ojs"},{"id":297,"value":"geometry = {\n const {cos,sin,PI} = Math;\n const {r1,dr1,r2,dr2,Tmis,Taz,Pmis,Paz} = params;\n return {\n HT: 20, // height of theta arm in mm\n RT: r1+1e-3*dr1, // radius of theta arm in mm\n HP: 65, // height of phi arm in mm\n RP: r2+1e-3*dr2, // radius of phi arm in mm\n axisT: unitvec(Tmis,Taz), // un-normalized vector along T axis of rotation in world coords\n axisP: unitvec(Pmis,Paz), // un-normalized vector along P axis of rotation relative to T arm\n };\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":404,"value":"function unitvec(theta, phi) {\n const {cos,sin,PI} = Math;\n theta = theta * 1e-3; // mrad to rad\n phi = phi * PI/180; // deg to rad\n const cth = cos(theta), sth = sin(theta);\n return [ sth * cos(phi), sth * sin(phi), cth ];\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":277,"value":"function rodriguesRotate(v, n, phi) {\n\n const [nx,ny,nz] = n;\n const cosPhi = Math.cos(phi);\n const sinPhi = Math.sin(phi);\n const dot = v[0]*nx + v[1]*ny + v[2]*nz;\n\n return [\n v[0]*cosPhi + (ny*v[2] - nz*v[1])*sinPhi + nx*dot*(1 - cosPhi),\n v[1]*cosPhi + (nz*v[0] - nx*v[2])*sinPhi + ny*dot*(1 - cosPhi),\n v[2]*cosPhi + (nx*v[1] - ny*v[0])*sinPhi + nz*dot*(1 - cosPhi)\n ];\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":294,"value":"function applyRotation(R, v) {\n // Multiply R * v, where R = [col0, col1, col2]\n return [\n R[0][0]*v[0] + R[1][0]*v[1] + R[2][0]*v[2],\n R[0][1]*v[0] + R[1][1]*v[1] + R[2][1]*v[2],\n R[0][2]*v[0] + R[1][2]*v[1] + R[2][2]*v[2]\n ];\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":291,"value":"function buildRotationMatrix(n) {\n // Input: unit vector n = [nx, ny, nz]\n const [nx, ny, nz] = n;\n const eps = 1e-6;\n\n // Step 1: z' = n\n const z = [nx, ny, nz];\n\n // Step 2: choose a reference vector not parallel to n\n let ref = [0, 1, 0]; // ŷ\n if (Math.abs(ny) > 0.999) ref = [1, 0, 0]; // fall back to x̂ if too parallel\n\n // Step 3: x' = ref × z\n let x = [\n ref[1]*z[2] - ref[2]*z[1],\n ref[2]*z[0] - ref[0]*z[2],\n ref[0]*z[1] - ref[1]*z[0]\n ];\n const xnorm = Math.hypot(...x);\n if (xnorm < eps) throw \"buildRotationMatrix: degenerate cross product\";\n x = x.map(v => v / xnorm);\n\n // Step 4: y' = z × x\n const y = [\n z[1]*x[2] - z[2]*x[1],\n z[2]*x[0] - z[0]*x[2],\n z[0]*x[1] - z[1]*x[0]\n ];\n\n // Return 3 columns: [x', y', z']\n return [x, y, z];\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":310,"value":"function offset(a, b) {\n return [ a[0]+b[0], a[1]+b[1], a[2]+b[2] ];\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":303,"value":"function normalize(n) {\n // Ensure n is a unit vector\n const n_norm = Math.hypot(n[0], n[1], n[2]);\n return [ n[0]/n_norm, n[1]/n_norm, n[2]/n_norm ];\n}","pinned":false,"hidden":false,"mode":"ojs"}],"files":"$W35","stdlib":"1","titling":"auto"}