force plate calibration

%23%20%2F%2F%2F%20script%0A%23%20requires-python%20%3D%20%22%3E%3D3.11%22%0A%23%20dependencies%20%3D%20%5B%22marimo%3E%3D0.23.9%22%2C%20%22numpy%3E%3D2.0%22%2C%20%22matplotlib%3E%3D3.9%22%2C%20%22scipy%3E%3D1.14%22%2C%20%22seaborn%3E%3D0.13%22%2C%20%22urllib3%3E%3D2.2%22%5D%0A%23%20%2F%2F%2F%0A%0Aimport%20marimo%0A%0A__generated_with%20%3D%20%220.25.2%22%0Aapp%20%3D%20marimo.App()%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%20Algorithm%20for%20force%20plate%20calibration%0A%0A%20%20%20%20%3E%20Marcos%20Duarte%2C%20Renato%20Naville%20Watanabe%3Cbr%3E%0A%20%20%20%20%3E%20%5BLaboratory%20of%20Biomechanics%20and%20Motor%20Control%5D(https%3A%2F%2Fbmclab.pesquisa.ufabc.edu.br%2F)%3Cbr%3E%0A%20%20%20%20%3E%20Federal%20University%20of%20ABC%2C%20Brazil%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20How%20to%20use%20this%20notebook%0A%0A%20%20%20%20Follow%20the%20six%20measured%20loads%20through%20the%20matrix%20model%2C%20a%20simulation%20with%20known%20answers%2C%20and%20the%20alternating%20estimates%20of%20calibration%20and%20load-cell%20orientation.%20Keep%20N%2C%20N%20m%20and%20m%20visible%20in%20your%20calculations.%20Try%20each%20exercise%20before%20opening%20its%20derivation%2C%20then%20use%20the%20final%20problems%20as%20an%20independent%20check.%0A%0A%20%20%20%20%23%23%20Opening%20challenge%0A%0A%20%20%20%20A%20plate%20reports%20a%20force%20that%20is%205%25%20too%20large.%20A%20calibration%20can%20correct%20that%20gain.%20But%20what%20if%20a%20sideways%20force%20also%20leaks%20into%20a%20moment%20channel%2C%20and%20the%20reference%20load%20cell%20is%20rotated%3F%20Can%20one%20experiment%20identify%20all%20of%20those%20effects%3F%0A%0A%20%20%20%20%23%23%20Setup%0A%0A%20%20%20%20NumPy%2C%20SciPy%2C%20Matplotlib%20and%20seaborn%20support%20the%20retained%20simulations%20and%20optimization.%20Seeded%20synthetic%20forces%20make%20repeated%20runs%20comparable.%20The%20example%20has%20known%20true%20calibration%20and%20rotation%20angles%3B%20no%20laboratory%20acquisition%20is%20needed.%20Agreement%20with%20this%20synthetic%20case%20verifies%20the%20calculation%20under%20its%20assumptions%2C%20not%20the%20accuracy%20of%20a%20physical%20plate.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20This%20notebook%20demonstrates%20the%20algorithm%20for%20force%20plate%20calibration%20proposed%20by%20Cedraro%20et%20al.%20(2008%2C%202009).%0A%0A%20%20%20%20A%20force%20plate%20(FP)%20is%20an%20electromechanical%20device%20that%20measures%20the%20components%20of%20the%20vectors%20force%20%24(%5Cmathbf%7BF%7D%20%3D%20%5BF_X%2C%5C%2C%20F_Y%2C%5C%2C%20F_Z%5D)%24%20and%20moment%20of%20force%20or%20torque%20%24(%5Cmathbf%7BM%7D%20%3D%20%5BM_X%2C%5C%2C%20M_Y%2C%5C%2C%20M_Z%5D)%24%20applied%20to%20the%20it.%20The%20FP%20is%20composed%20by%20a%20transducer%20that%20transforms%20a%20mechanical%20deformation%20to%20an%20electrical%20signal%20usually%20using%20strain%20gauges%20or%20piezoelectric%20sensors.%20The%20transformation%20from%20electrical%20signals%20(input)%20to%20force%20and%20moment%20of%20force%20(output)%20as%20a%20function%20of%20time%20in%20a%20six-component%20FP%20usually%20is%20given%20by%20the%20following%20linear%20relationship%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BL%7D(t)%20%3D%20%5Cmathbf%7BC%7D%5Cmathbf%7BV%7D(t)%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Where%20%24%5Cmathbf%7BL%7D(t)%24%20is%20the%20force%20plate%20output%20vector%20%24(%5B%5Cmathbf%7BF%7D(t)%2C%20%5Cmathbf%7BM%7D(t)%5D%5ET)%24%2C%20in%20N%20and%20Nm%2C%20%24%5Cmathbf%7BV%7D(t)%24%20is%20the%20vector%20of%20electrical%20signals%20(six%20voltage%20signals%2C%20in%20V)%20and%20%20%24%5Cmathbf%7BC%7D%24%20is%20known%20as%20the%20six-by-six%20(constant)%20calibration%20matrix%20(in%20N%2FV%20or%20Nm%2FV).%20Note%20that%20we%20used%20the%20term%20vector%20here%20to%20refer%20to%20an%20uni-dimensional%20matrix%20(usual%20in%20scientific%20computing)%2C%20which%20is%20different%20from%20vector%2Fscalar%20concept%20in%20Mechanics.%3Cbr%3E%0A%20%20%20%20The%20expansion%20of%20the%20former%20equiation%20at%20a%20given%20instant%20is%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cbegin%7Bbmatrix%7D%0A%20%20%20%20F_x%20%5C%5C%20F_y%20%5C%5C%20F_z%20%5C%5C%20M_x%20%5C%5C%20M_y%20%5C%5C%20M_z%0A%20%20%20%20%5Cend%7Bbmatrix%7D%5C%2C%20%3D%20%5C%2C%0A%20%20%20%20%5Cbegin%7Bbmatrix%7D%0A%20%20%20%20C_%7B11%7D%20%26%26%20C_%7B12%7D%20%26%26%20C_%7B13%7D%20%26%26%20C_%7B14%7D%20%26%26%20C_%7B15%7D%20%26%26%20C_%7B16%7D%20%5C%5C%0A%20%20%20%20C_%7B21%7D%20%26%26%20C_%7B22%7D%20%26%26%20C_%7B23%7D%20%26%26%20C_%7B24%7D%20%26%26%20C_%7B25%7D%20%26%26%20C_%7B26%7D%20%5C%5C%0A%20%20%20%20C_%7B31%7D%20%26%26%20C_%7B32%7D%20%26%26%20C_%7B33%7D%20%26%26%20C_%7B34%7D%20%26%26%20C_%7B35%7D%20%26%26%20C_%7B36%7D%20%5C%5C%0A%20%20%20%20C_%7B41%7D%20%26%26%20C_%7B42%7D%20%26%26%20C_%7B43%7D%20%26%26%20C_%7B44%7D%20%26%26%20C_%7B45%7D%20%26%26%20C_%7B46%7D%20%5C%5C%0A%20%20%20%20C_%7B51%7D%20%26%26%20C_%7B52%7D%20%26%26%20C_%7B53%7D%20%26%26%20C_%7B54%7D%20%26%26%20C_%7B55%7D%20%26%26%20C_%7B56%7D%20%5C%5C%0A%20%20%20%20C_%7B61%7D%20%26%26%20C_%7B62%7D%20%26%26%20C_%7B63%7D%20%26%26%20C_%7B64%7D%20%26%26%20C_%7B65%7D%20%26%26%20C_%7B66%7D%0A%20%20%20%20%5Cend%7Bbmatrix%7D%5C%2C%0A%20%20%20%20%5Cbegin%7Bbmatrix%7D%0A%20%20%20%20V_1%20%5C%5C%20V_2%20%5C%5C%20V_3%20%5C%5C%20V_4%20%5C%5C%20V_5%20%5C%5C%20V_6%0A%20%20%20%20%5Cend%7Bbmatrix%7D%0A%20%20%20%20%24%24%0A%0A%20%20%20%20The%20terms%20off-diagonal%20are%20known%20as%20the%20crosstalk%20terms%20and%20represent%20the%20effect%20of%20a%20load%20applied%20in%20one%20direction%20on%20the%20other%20direction.%20For%20a%20FP%20with%20none%20or%20small%20crosstalk%2C%20the%20off-diagonal%20terms%20are%20zero%20or%20very%20small%20compared%20to%20the%20main-diagonal%20terms.%20Note%20that%20the%20equation%20above%20is%20in%20fact%20a%20system%20of%20six%20linear%20independent%20equations%20with%20six%20unknowns%20each%20(where%20%24V_1%20...%20V_6%24%20are%20the%20measured%20inputs)%3A%0A%0A%20%20%20%20%5Cbegin%7Bcases%7D%0A%20%20%20%20%20%20%20%20F_x%20%26%3D%26%20C_%7B11%7DV_1%20%2B%20C_%7B12%7DV_2%20%2B%20C_%7B13%7DV_3%20%2B%20C_%7B14%7DV_4%20%2B%20C_%7B15%7DV_5%20%2B%20C_%7B16%7DV_6%20%5C%5C%0A%20%20%20%20%20%20%20%20F_y%20%26%3D%26%20C_%7B21%7DV_1%20%2B%20C_%7B22%7DV_2%20%2B%20C_%7B23%7DV_3%20%2B%20C_%7B24%7DV_4%20%2B%20C_%7B25%7DV_5%20%2B%20C_%7B26%7DV_6%20%5C%5C%0A%20%20%20%20%20%20%20%20F_z%20%26%3D%26%20C_%7B31%7DV_1%20%2B%20C_%7B32%7DV_2%20%2B%20C_%7B33%7DV_3%20%2B%20C_%7B34%7DV_4%20%2B%20C_%7B35%7DV_5%20%2B%20C_%7B36%7DV_6%20%5C%5C%0A%20%20%20%20%20%20%20%20M_x%20%26%3D%26%20C_%7B41%7DV_1%20%2B%20C_%7B42%7DV_2%20%2B%20C_%7B43%7DV_3%20%2B%20C_%7B44%7DV_4%20%2B%20C_%7B45%7DV_5%20%2B%20C_%7B46%7DV_6%20%5C%5C%0A%20%20%20%20%20%20%20%20M_y%20%26%3D%26%20C_%7B51%7DV_1%20%2B%20C_%7B52%7DV_2%20%2B%20C_%7B53%7DV_3%20%2B%20C_%7B54%7DV_4%20%2B%20C_%7B55%7DV_5%20%2B%20C_%7B56%7DV_6%20%5C%5C%0A%20%20%20%20%20%20%20%20M_z%20%26%3D%26%20C_%7B61%7DV_1%20%2B%20C_%7B62%7DV_2%20%2B%20C_%7B63%7DV_3%20%2B%20C_%7B64%7DV_4%20%2B%20C_%7B65%7DV_5%20%2B%20C_%7B66%7DV_6%0A%20%20%20%20%5Cend%7Bcases%7D%0A%0A%20%20%20%20Of%20course%2C%20an%20important%20aspect%20of%20the%20FP%20functionning%20is%20that%20it%20should%20be%20calibrated%2C%20i.e.%2C%20the%20calibration%20matrix%20must%20be%20known%20and%20accurate%20(it%20comes%20with%20the%20force%20plate%20when%20you%20buy%20one).%20Cedraro%20et%20al.%20(2008)%20proposed%20a%20method%20for%20in%20situ%20re-calibration%20of%20FP%20and%20their%20algorithm%20is%20presented%20next.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Algorithm%0A%0A%20%20%20%20Consider%20that%20in%20a%20re-calibration%20procedure%20we%20apply%20on%20the%20FP%20known%20forces%2C%20%24%5Cmathbf%7BF%7D_I%20%3D%20%5BF_%7BX_I%7D%2C%5C%2C%20F_%7BY_I%7D%2C%5C%2C%20F_%7BZ_I%7D%5D%5ET%24%2C%20at%20known%20places%2C%20%24%5Cmathbf%7BCOP%7D%20%3D%20%5BX_%7BCOP%7D%2C%5C%2C%20%20Y_%7BCOP%7D%2C%5C%2C%20%20Z_%7BCOP%7D%5D%24%20(the%20center%20of%20pressure%20coordinates%20in%20the%20FP%20reference%20frame).%3Cbr%3E%0A%20%20%20%20The%20moments%20of%20forces%2C%20%24%5Cmathbf%7BM%7D_I%20%3D%20%5BM_%7BX_I%7D%2C%5C%2C%20M_%7BY_I%7D%2C%5C%2C%20M_%7BZ_I%7D%5D%5ET%24%2C%20due%20to%20these%20forces%20can%20be%20found%20using%20the%20equation%20%24%5Cmathbf%7BM%7D_I%20%3D%20%5Cmathbf%7BCOP%7D%20%5Ctimes%20%5Cmathbf%7BF%7D_I%24%2C%20which%20can%20be%20expressed%20in%20matrix%20form%20as%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BM%7D_I%20%3D%0A%20%20%20%20%5Cbegin%7Bbmatrix%7D%0A%20%20%20%200%20%26%26%20-Z_%7BCOP%7D%20%26%26%20Y_%7BCOP%7D%20%5C%5C%0A%20%20%20%20Z_%7BCOP%7D%20%26%26%200%20%26%26%20-X_%7BCOP%7D%20%5C%5C%0A%20%20%20%20-Y_%7BCOP%7D%20%26%26%20X_%7BCOP%7D%20%26%26%200%0A%20%20%20%20%5Cend%7Bbmatrix%7D%5C%2C%20%5Cmathbf%7BF%7D_I%20%5C%2C%20%3D%20%5C%2C%20%5Cmathbf%7BA%7D_%7BCOP%7D%5Cmathbf%7BF%7D_I%0A%20%20%20%20%24%24%0A%0A%20%20%20%20%24%5Cmathbf%7BA%7D_%7BCOP%7D%24%20(a%20%5Bskew-symmetric%20matrix%5D(https%3A%2F%2Fen.wikipedia.org%2Fwiki%2FSkew-symmetric_matrix))%20is%20simply%20the%20COP%20position%20in%20matrix%20form%20in%20order%20to%20calculate%20the%20%5Bcross%20product%20with%20matrix%20multiplication%5D(https%3A%2F%2Fen.wikipedia.org%2Fwiki%2FCross_product).%0A%0A%20%20%20%20These%20known%20loads%20on%20the%20FP%20can%20also%20be%20represented%20as%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BL%7D_I%20%3D%0A%20%20%20%20%5Cbegin%7Bbmatrix%7D%0A%20%20%20%20%5Cmathbf%7BF%7D_I%20%5C%5C%0A%20%20%20%20%5Cmathbf%7BM%7D_I%0A%20%20%20%20%5Cend%7Bbmatrix%7D%0A%20%20%20%20%24%24%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Exercise%201%20%E2%80%94%20known%20moments.**%20A%20force%20%24(0%2C0%2C500)%24%20N%20is%20applied%20at%20%24(0.1%2C-0.2%2C0)%24%20m.%20Calculate%20the%20moment%20about%20the%20origin%20and%20the%20full%20six-component%20known%20load.%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Worked%20solution%20%E2%80%94%20Exercise%201%20%E2%80%94%20known%20moments%0A%0A%20%20%20%20Use%20%24%5Cmathbf%20M%3D%5Cmathbf%20r%5Ctimes%5Cmathbf%20F%3D(-100%2C-50%2C0)%24%20N%20m.%20Thus%20%24%5Cmathbf%20L_I%3D(0%2C0%2C500%2C-100%2C-50%2C0)%5ET%24.%20The%20cross-product%20matrix%20has%20units%20of%20length%2C%20which%20is%20why%20the%20moment%20block%20cannot%20be%20treated%20as%20having%20force%20units.%0A%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%23%20Linear%20re-calibration%0A%0A%20%20%20%20For%20a%20linear%20re-calibration%2C%20the%20relationship%20between%20the%20measured%20FP%20output%2C%20%24%5Cmathbf%7BL%7D%24%2C%20and%20the%20known%20loads%2C%20%24%5Cmathbf%7BL%7D_I%24%2C%20is%20approximated%20by%20a%20linear%20equation%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BL%7D_I%20%3D%20%5Cmathbf%7BC%7D%5Cmathbf%7BL%7D%20%2B%20%5Cmathbf%7BE%7D%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Where%20%24%5Cmathbf%7BC%7D%24%20now%20is%20the%20six-by-six%20re-calibration%20matrix%20(with%20units%20that%20depend%20on%20the%20block%3A%20force%2Fforce%20and%20moment%2Fmoment%20are%20dimensionless%2C%20moment%2Fforce%20has%20units%20m%2C%20and%20force%2Fmoment%20has%20units%201%2Fm)%20and%20%24%5Cmathbf%7BE%7D%24%20is%20a%20gaussian%2C%20uncorrelated%2C%20zero%20mean%0A%20%20%20%20noise%20six-by-one%20matrix.%3Cbr%3E%0A%20%20%20%20The%20re-calibration%20matrix%20can%20be%20found%20by%20solving%20the%20equation%20above%20and%20then%20%24%5Cmathbf%7BC%7D%24%20can%20be%20later%20used%20to%20re-calibrate%20the%20FP%20output%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BL%7D_C%20%3D%20%5Cmathbf%7BC%7D%5Cmathbf%7BL%7D%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Where%20%24%5Cmathbf%7BL%7D_C%24%20is%20the%20re-calibrated%20FP%20output.%20For%20a%20perfectly%20calibrated%20FP%2C%20%24%5Cmathbf%7BL%7D_C%20%3D%20%5Cmathbf%7BL%7D%24%20and%20%24%5Cmathbf%7BC%7D%20%3D%20%5Cmathbf%7BI%7D%24%2C%20the%20six-by-six%20identity%20matrix.%0A%0A%20%20%20%20Cedraro%20et%20al.%20(2008%2C%202009)%20proposed%20to%20use%20a%20calibrated%20three-component%20load%20cell%20(LC)%20to%20measure%20the%20loads%20%24%5Cmathbf%7BF%7D_I(t)%24%20applied%20on%20the%20FP%20at%20%24k%24%20known%20measurements%20sites.%20The%20LC%20measures%20the%20loads%20in%20its%20own%20coordinate%20system%20%24(xyz)%24%3A%20%24%5Cmathbf%7BF%7D_%7BLC%7D(t)%20%3D%20%5BF_x(t)%2C%5C%2C%20F_y(t)%2C%5C%2C%20F_z(t)%5D%5ET%24%2C%20which%20is%20probaly%20rotated%20(by%20an%20unknown%20value%2C%20represented%20by%20rotation%20matrix%20%24%5Cmathbf%7BR%7D%5Ek%24)%20in%20relation%20to%20the%20FP%20coordinate%20system%20(the%20coordinate%20systems%20are%20also%20translated%20to%20each%20other%20but%20the%20translation%20is%20known%20and%20given%20by%20the%20COP%20position).%3Cbr%3E%0A%20%20%20%20For%20each%20measurement%20site%2C%20the%20equation%20for%20the%20determination%20of%20the%20re-calibration%20matrix%20will%20be%20given%20by%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BP%7D%5Ek%5Cmathbf%7BR%7D%5Ek%5Cmathbf%7BF%7D%5Ek_%7BLC%7D(t)%3D%20%5Cmathbf%7BP%7D%5Ek%5Cmathbf%7BF%7D_I%5Ek(t)%20%3D%20%5Cmathbf%7BC%7D%5Cmathbf%7BL%7D%5Ek(t)%20%2B%20%5Cmathbf%7BE%7D%5Ek(t)%20%5Cquad%20k%20%3D%201%2C%20...%2C%20n%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Where%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BP%7D%5Ek%20%3D%0A%20%20%20%20%5Cbegin%7Bbmatrix%7D%0A%20%20%20%20%5Cmathbf%7BI%7D_3%20%5C%5C%0A%20%20%20%20%5Cmathbf%7BA%7D_%7BCOP%7D%0A%20%20%20%20%5Cend%7Bbmatrix%7D%0A%20%20%20%20%24%24%0A%0A%20%20%20%20and%20%24I_3%24%20is%20the%20three-by-three%20identity%20matrix.%0A%0A%20%20%20%20Using%20a%20typical%20load%20cell%2C%20with%20a%20flat%20bottom%2C%20on%20top%20the%20FP%2C%20a%20realistic%20assumption%20is%20to%20consider%20that%20%24z%24%20of%20LC%20is%20aligned%20to%20%24Z%24%20of%20FP%20(the%20vertical%20direction)%3B%20in%20this%20case%20the%20rotation%20matrix%20is%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BR%7D%5Ek%20%3D%0A%20%20%20%20%5Cbegin%7Bbmatrix%7D%0A%20%20%20%20%5Ccos%5Calpha%5Ek%20%26%26%20-%5Csin%5Calpha%5Ek%20%26%26%200%20%5C%5C%0A%20%20%20%20%5Csin%5Calpha%5Ek%20%26%26%20%5Ccos%5Calpha%5Ek%20%26%26%200%20%5C%5C%0A%20%20%20%200%20%26%26%200%20%26%26%201%0A%20%20%20%20%5Cend%7Bbmatrix%7D%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Cedraro%20et%20al.%20(2008)%20propose%20the%20following%20algorithm%20to%20estimate%20%24%5Cmathbf%7BC%7D%24%3A%0A%20%20%20%201.%20The%20misalignments%2C%20%24%5Calpha%5Ek%24%2C%20are%20initialized%3A%20%24%5Cmathbf%7B%5Calpha%7D%20%3D%20%5B%5Calpha%5E1%2C%20%5Ccdots%2C%20%5Calpha%5En%5D%24%3B%0A%20%20%20%202.%20%24%5Cmathbf%7BC%7D%24%20is%20calculated%20by%20a%20least-squares%20approach%3B%0A%20%20%20%203.%20The%20residual%20errors%20are%20estimated%20as%3A%20%24%5Cmathbf%7BE%7D%5Ek(t)%20%3D%20%5Cmathbf%7BP%7D%5Ek%5Cmathbf%7BR%7D%5Ek%5Cmathbf%7BF%7D%5Ek_%7BLC%7D(t)%20-%20%5Cmathbf%7BC%7D%5Cmathbf%7BL%7D%5Ek(t)%24%3B%0A%20%20%20%204.%20The%20increment%20%24%5Cmathbf%7B%5CDelta%5Calpha%7D%24%20is%20calculated%20by%20minimizing%20the%20cost%20function%20%24%5Csum_%7Bk%2Ct%7D%5Cmathbf%7BE%7D%5Ek(t)%5ET%5Cmathbf%7BE%7D%5Ek(t)%24%2C%20assuming%20dimensional%20unitary%20weights%3B%0A%20%20%20%205.%20The%20parameters%20are%20updated%3A%20%24%5Cmathbf%7B%5Calpha%7D%20%3D%20%5Cmathbf%7B%5Calpha%7D%20%2B%20%5Cmathbf%7B%5CDelta%5Calpha%7D%24.%0A%0A%20%20%20%20The%20iteration%20of%20steps%202%E2%80%935%20stops%20when%20each%20%24%7C%5CDelta%5Calpha%5Ek%7C%20%3C%20%5Cvarepsilon_0%24%2C%20where%20%24%5Cvarepsilon_0%3D10%5E%7B-10%7D%24%20is%20the%20chosen%20threshold.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20Before%20fitting%2C%20distinguish%20a%20voltage%20calibration%20from%20a%20re-calibration%20of%20existing%20force%2Fmoment%20channels.%20The%20latter%20still%20mixes%20units%20across%20its%20off-diagonal%20blocks.%20Likewise%2C%20summing%20raw%20squared%20force%20and%20moment%20residuals%20implicitly%20chooses%20a%20numerical%20weighting%20tied%20to%20the%20chosen%20units.%20For%20experimental%20use%2C%20specify%20residual%20scales%20or%20an%20uncertainty-based%20weighting%20matrix%3B%20changing%20metres%20to%20millimetres%20without%20adjusting%20that%20weighting%20changes%20the%20fit.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%23%20Simulation%0A%0A%20%20%20%20Let's%20simulate%20some%20data%20to%20test%20this%20calibration%20procedure.%20Cedraro%20et%20al.%20(2008)%20employed%20sinusoids%2C%20cosenoids%2C%20and%20ramps%20as%20sintetic%20signals%20to%20simulate%20the%20calibration%20process%3A%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20import%20numpy%20as%20np%0A%20%20%20%20from%20numpy.linalg%20import%20inv%0A%20%20%20%20import%20matplotlib.pyplot%20as%20plt%0A%20%20%20%20%23%20'%25matplotlib%20inline'%20command%20supported%20automatically%20in%20marimo%0A%20%20%20%20import%20seaborn%20as%20sns%0A%20%20%20%20sns.set_context(%22notebook%22%2C%20font_scale%3D1.4%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20rc%3D%7B%22lines.linewidth%22%3A%203%2C%20%22lines.markersize%22%3A%208%2C%20%22axes.titlesize%22%3A%20'x-large'%7D)%0A%20%20%20%20return%20inv%2C%20np%2C%20plt%0A%0A%0A%40app.cell%0Adef%20_(np%2C%20plt)%3A%0A%20%20%20%20%23%20simulated%20forces%20measured%20by%20the%20load%20cell%20in%20its%20local%20coordinate%20system%0A%20%20%20%20_rng%20%3D%20np.random.default_rng(0)%0A%20%20%20%20samples%20%3D%20np.linspace(1%2C%206000%2C%206000)%0A%20%20%20%20ns%20%3D%20samples.shape%5B0%5D%0A%20%20%20%20Flc%20%3D%20np.array(%5B100*np.sin(5*2*np.pi*samples%2Fns)%20%2B%202*_rng.standard_normal(6000)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20100*np.cos(5*2*np.pi*samples%2Fns)%20%2B%202*_rng.standard_normal(6000)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20samples%2F15%20%2B%20200%20%2B%205*_rng.standard_normal(6000)%5D)%0A%20%20%20%20%23%20plots%0A%20%20%20%20fig%2C%20axs%20%3D%20plt.subplots(3%2C%201%2C%20figsize%3D(8%2C%205)%2C%20sharex%3D'all')%0A%20%20%20%20axs%5B0%5D.plot(samples%2C%20Flc%5B0%5D)%0A%20%20%20%20axs%5B0%5D.set_ylabel('Fx%20(N)')%0A%20%20%20%20axs%5B0%5D.locator_params(axis%3D'y'%2C%20nbins%3D3)%0A%20%20%20%20axs%5B0%5D.yaxis.set_label_coords(-.08%2C%200.5)%0A%20%20%20%20axs%5B1%5D.plot(samples%2C%20Flc%5B1%5D)%0A%20%20%20%20axs%5B1%5D.set_ylabel('Fy%20(N)')%0A%20%20%20%20axs%5B1%5D.locator_params(axis%3D'y'%2C%20nbins%3D3)%0A%20%20%20%20axs%5B1%5D.yaxis.set_label_coords(-.08%2C%200.5)%0A%20%20%20%20axs%5B2%5D.plot(samples%2C%20Flc%5B2%5D)%0A%20%20%20%20axs%5B2%5D.set_ylabel('Fz%20(N)')%0A%20%20%20%20axs%5B2%5D.set_xlabel('Samples')%0A%20%20%20%20axs%5B2%5D.locator_params(axis%3D'y'%2C%20nbins%3D3)%0A%20%20%20%20axs%5B2%5D.yaxis.set_label_coords(-.08%2C%200.5)%0A%20%20%20%20plt.tight_layout(pad%3D.5%2C%20h_pad%3D.025)%0A%20%20%20%20plt.show()%0A%20%20%20%20return%20Flc%2C%20ns%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20And%20Cedraro%20et%20al.%20(2008)%20also%20proposed%20five%20measurement%20sites%20and%20a%20re-calibration%20matrix%20for%20the%20simulated%20re-calibration%3A%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(Flc%2C%20np)%3A%0A%20%20%20%20C%20%3D%20np.array(%5B%5B1.0354%2C%20-0.0053%2C%20-0.0021%2C%20-0.0289%2C%20-0.0402%2C%200.0081%5D%2C%20%5B0.0064%2C%201.0309%2C%20-0.0031%2C%200.0211%2C%200.0135%2C%20-0.0001%5D%2C%20%5B0.0%2C%20-0.0004%2C%201.0022%2C%20-0.0005%2C%20-0.0182%2C%200.03%5D%2C%20%5B-0.0012%2C%20-0.0385%2C%200.0002%2C%200.9328%2C%200.0007%2C%200.0017%5D%2C%20%5B0.0347%2C%200.0003%2C%200.0008%2C%20-0.0002%2C%200.9325%2C%20-0.0024%5D%2C%20%5B-0.0004%2C%20-0.0013%2C%20-0.0003%2C%20-0.0023%2C%200.0035%2C%201.0592%5D%5D)%0A%20%20%20%20COP%20%3D%20np.array(%5B%5B0%2C%20112%2C%20112%2C%20-112%2C%20-112%5D%2C%20%5B0%2C%20192%2C%20-192%2C%20192%2C%20-192%5D%2C%20%5B-124%2C%20-124%2C%20-124%2C%20-124%2C%20-124%5D%5D)%20%2F%201000%0A%20%20%20%20nk%20%3D%20COP.shape%5B1%5D%0A%20%20%20%20Acop%20%3D%20lambda%20x%2C%20y%2C%20z%3A%20np.array(%5B%5B0.0%2C%20-z%2C%20y%5D%2C%20%5Bz%2C%200.0%2C%20-x%5D%2C%20%5B-y%2C%20x%2C%200.0%5D%5D)%0A%20%20%20%20Flc_1%20%3D%20np.tile(Flc%2C%20nk)%0A%20%20%20%20return%20Acop%2C%20C%2C%20COP%2C%20Flc_1%2C%20nk%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20Let's%20generate%20the%20loads%20measured%20by%20the%20FP%20given%20the%20re-calibration%20matrix%20and%20the%20simulated%20forces%20measured%20by%20the%20load%20cell%20(we%20will%20consider%20no%20rotation%20for%20now).%20For%20that%20we%20will%20have%20to%20solve%20the%20equation%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BL%7D_I%20%3D%20%5Cmathbf%7BC%7D%5Cmathbf%7BL%7D%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Which%20is%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BL%7D%20%3D%20%5Cmathbf%7BC%7D%5E%7B-1%7D%5Cmathbf%7BL%7D_I%0A%20%20%20%20%24%24%0A%0A%20%20%20%20%24%5Cmathbf%7BC%7D%24%20is%20a%20square%20(6-by-6)%20matrix%20and%20the%20computation%20of%20its%20inverse%20is%20straightforward.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(Acop%2C%20C%2C%20COP%2C%20Flc_1%2C%20inv%2C%20nk%2C%20np%2C%20ns)%3A%0A%20%20%20%20Li%20%3D%20np.empty((6%2C%20ns%20*%20nk))%0A%20%20%20%20P%20%3D%20np.empty((6%2C%203%2C%20nk))%0A%20%20%20%20for%20(k%2C%20cop)%20in%20enumerate(COP.T)%3A%0A%20%20%20%20%20%20%20%20P%5B%3A%2C%20%3A%2C%20k%5D%20%3D%20np.vstack((np.eye(3)%2C%20Acop(*cop)))%0A%20%20%20%20%20%20%20%20Li%5B%3A%2C%20k%20*%20ns%3A(k%20%2B%201)%20*%20ns%5D%20%3D%20P%5B%3A%2C%20%3A%2C%20k%5D%20%40%20Flc_1%5B%3A%2C%20k%20*%20ns%3A(k%20%2B%201)%20*%20ns%5D%0A%20%20%20%20L%20%3D%20inv(C)%20%40%20Li%0A%20%20%20%20return%20L%2C%20Li%0A%0A%0A%40app.cell%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20In%20the%20calculations%20above%20we%20took%20advantage%20of%20the%20%5Bnew%20operator%20for%20matrix%20multiplcation%20in%20Python%203%5D(https%3A%2F%2Fwww.python.org%2Fdev%2Fpeps%2Fpep-0465%2F)%3A%20%60%40%60%20(mnemonic%3A%20%60%40%60%20is%20%60*%60%20for%20mATrices).%0A%0A%20%20%20%20We%20can%20now%20simulate%20the%20re-calibration%20procedure%20by%20determining%20the%20re-calibration%20matrix%20using%20these%20loads.%20Of%20course%2C%20the%20re-calibration%20matrix%20to%20be%20determined%20should%20be%20equal%20to%20the%20simulated%20re-calibration%20matrix%20we%20started%20with%2C%20but%20this%20is%20the%20fun%20of%20the%20simulation%20-%20we%20know%20where%20we%20want%20to%20go.%0A%0A%20%20%20%20The%20re-calibration%20matrix%20can%20be%20found%20by%20solving%20the%20following%20equation%20(considering%20the%20angles%20equal%20zero%20for%20now)%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BL%7D_I%20%3D%20%5Cmathbf%7BC%7D%5Cmathbf%7BL%7D%0A%20%20%20%20%24%24%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BL%7D_I%20%5Cmathbf%7BL%7D%5E%7B-1%7D%20%3D%20%5Cmathbf%7BC%7D%5Cmathbf%7BL%7D%20%5Cmathbf%7BL%7D%5E%7B-1%7D%20%3D%20%5Cmathbf%7BC%7D%5Cmathbf%7BI%7D%0A%20%20%20%20%24%24%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BC%7D%20%3D%20%5Cmathbf%7BL%7D_I%5Cmathbf%7BL%7D%5E%7B-1%7D%0A%20%20%20%20%24%24%0A%0A%20%20%20%20The%20problem%20is%20that%20%24%5Cmathbf%7BL%7D%24%20in%20general%20is%20a%20non-square%20matrix%20and%20its%20inverse%20is%20not%20defined%20(unless%20you%20perform%20exactly%20six%20measurements%20and%20then%20%24%5Cmathbf%7BL%7D%24%20would%20be%20a%20six-by-six%20square%20matrix%2C%20but%20this%20is%20too%20restrictive).%20However%2C%20we%20still%20can%20solve%20the%20equation%20with%20some%20extra%20manipulation%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BL%7D_I%20%3D%20%5Cmathbf%7BC%7D%5Cmathbf%7BL%7D%0A%20%20%20%20%24%24%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BL%7D_I%20%5Cmathbf%7BL%7D%5ET%20%3D%20%5Cmathbf%7BC%7D%5Cmathbf%7BL%7D%20%5Cmathbf%7BL%7D%5ET%0A%20%20%20%20%24%24%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BL%7D_I%20%5Cmathbf%7BL%7D%5ET(%5Cmathbf%7BL%7D%5Cmathbf%7BL%7D%5ET)%5E%7B-1%7D%20%3D%20%5Cmathbf%7BC%7D%5Cmathbf%7BL%7D%20%5Cmathbf%7BL%7D%5ET%20(%5Cmathbf%7BL%7D%5Cmathbf%7BL%7D%5ET)%5E%7B-1%7D%20%3D%20%5Cmathbf%7BC%7D%5Cmathbf%7BI%7D%0A%20%20%20%20%24%24%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BC%7D%20%3D%20%5Cmathbf%7BL%7D_I%5Cmathbf%7BL%7D%5ET(%5Cmathbf%7BL%7D%5Cmathbf%7BL%7D%5ET)%5E%7B-1%7D%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Note%20that%20%24%5Cmathbf%7BL%7D%20%5Cmathbf%7BL%7D%5ET%24%20is%20a%20square%20matrix%20and%20is%20invertible%20(also%20%5Bnonsingular%5D(https%3A%2F%2Fen.wikipedia.org%2Fwiki%2FInvertible_matrix))%20if%20%24%5Cmathbf%7BL%7D%24%20is%20L.I.%20(%5Blinearly%20independent%20rows%2Fcolumns%5D(https%3A%2F%2Fen.wikipedia.org%2Fwiki%2FLinear_independence)).%20The%20matrix%20%24%5Cmathbf%7BL%7D%5ET(%5Cmathbf%7BL%7D%5Cmathbf%7BL%7D%5ET)%5E%7B-1%7D%24%20is%20known%20as%20the%20%5Bgeneralized%20inverse%20or%20Moore%E2%80%93Penrose%20pseudoinverse%5D(https%3A%2F%2Fen.wikipedia.org%2Fwiki%2FMoore%25E2%2580%2593Penrose_pseudoinverse)%2C%20a%20generalization%20of%20the%20inverse%20matrix.%20If%20we%20denote%20this%20pseudoinverse%20matrix%20by%20%24%5Cmathbf%7BL%7D%5E%2B%24%2C%20we%20can%20state%20the%20solution%20of%20the%20equation%20simply%20as%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BL%7D_I%20%3D%20%5Cmathbf%7BC%7D%5Cmathbf%7BL%7D%0A%20%20%20%20%24%24%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cmathbf%7BC%7D%20%3D%20%5Cmathbf%7BL%7D_I%20%5Cmathbf%7BL%7D%5E%2B%0A%20%20%20%20%24%24%0A%0A%20%20%20%20To%20compute%20the%20Moore%E2%80%93Penrose%20pseudoinverse%2C%20we%20could%20calculate%20it%20by%20the%20naive%20approach%20in%20Python%3A%0A%20%20%20%20%60%60%60python%0A%20%20%20%20Linv%20%3D%20L.T%20%40%20inv(L%20%40%20L.T)%0A%20%20%20%20%60%60%60%0A%20%20%20%20But%20both%20Numpy%20and%20Scipy%20have%20functions%20to%20calculate%20the%20pseudoinverse%2C%20which%20might%20give%20greater%20numerical%20stability%20(but%20read%20%5BInverses%20and%20pseudoinverses.%20Numerical%20issues%2C%20speed%2C%20symmetry%5D(http%3A%2F%2Fvene.ro%2Fblog%2Finverses-pseudoinverses-numerical-issues-speed-symmetry.html)).%3Cbr%3E%0A%20%20%20%20Of%20note%2C%20%5Bnumpy.linalg.pinv%5D(https%3A%2F%2Fnumpy.org%2Fdoc%2Fstable%2Freference%2Fgenerated%2Fnumpy.linalg.pinv.html)%20calculates%20the%20pseudoinverse%20of%20a%20matrix%20using%20its%20singular-value%20decomposition%20(SVD)%20and%20including%20all%20large%20singular%20values%20(using%20the%20%5BLAPACK%20(Linear%20Algebra%20Package)%5D(https%3A%2F%2Fen.wikipedia.org%2Fwiki%2FLAPACK)%20routine%20gesdd)%2C%20whereas%20%5Bscipy.linalg.pinv%5D(https%3A%2F%2Fdocs.scipy.org%2Fdoc%2Fscipy%2Freference%2Fgenerated%2Fscipy.linalg.pinv.html)%20also%20uses%20the%20SVD%20(it%20used%20a%20least-squares%20solver%20until%20SciPy%201.7%2C%20when%20the%20SVD-based%20%60scipy.linalg.pinv2%60%20was%20folded%20into%20it%3B%20%60pinv2%60%20was%20removed%20in%20SciPy%201.9).%0A%20%20%20%20Let's%20use%20%5Bscipy.linalg.pinv%5D(https%3A%2F%2Fdocs.scipy.org%2Fdoc%2Fscipy%2Freference%2Fgenerated%2Fscipy.linalg.pinv.html)%3A%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(L)%3A%0A%20%20%20%20from%20scipy.linalg%20import%20pinv%20as%20_pinv%0A%20%20%20%20Lpinv%20%3D%20_pinv(L)%0A%20%20%20%20return%20(Lpinv%2C)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20Then%2C%20the%20re-calibration%20matrix%20is%3A%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(Li%2C%20Lpinv)%3A%0A%20%20%20%20C2%20%3D%20Li%20%40%20Lpinv%0A%20%20%20%20return%20(C2%2C)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20Which%20is%20indeed%20the%20same%20as%20the%20initial%20calibration%20matrix%3A%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(C%2C%20C2%2C%20np)%3A%0A%20%20%20%20np.allclose(C%2C%20C2)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20The%20residual%20error%20between%20the%20old%20loads%20and%20new%20loads%20after%20re-calibration%20is%3A%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(C2%2C%20L%2C%20Li%2C%20np)%3A%0A%20%20%20%20E%20%3D%20Li%20-%20C2%20%40%20L%0A%20%20%20%20e%20%3D%20np.sum(E%20*%20E)%0A%20%20%20%20print('Average%20residual%20error%20between%20old%20and%20new%20loads%3A'%2C%20e)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Exercise%202%20%E2%80%94%20rank%20and%20pseudoinverse.**%20A%20calibration%20experiment%20applies%20only%20vertical%20forces%20at%20the%20plate%20centre.%20Can%20it%20identify%20an%20arbitrary%20%246%5Ctimes6%24%20re-calibration%20matrix%3F%20How%20does%20an%20SVD%20pseudoinverse%20change%20the%20answer%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Worked%20solution%20%E2%80%94%20Exercise%202%20%E2%80%94%20rank%20and%20pseudoinverse%0A%0A%20%20%20%20The%20measured%20loads%20span%20only%20one%20direction%20in%20this%20idealized%20experiment%2C%20so%20infinitely%20many%20calibration%20matrices%20fit%20them.%20A%20pseudoinverse%20returns%20one%20minimum-norm%20least-squares%20solution%3B%20it%20does%20not%20create%20missing%20excitation.%20To%20identify%20the%20full%20linear%20map%2C%20the%20six%20load%20rows%20must%20have%20full%20row%20rank%2C%20and%20adequate%20conditioning%20and%20independent%20validation%20are%20also%20needed.%0A%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%23%20Optimization%0A%0A%20%20%20%20Let's%20now%20implement%20the%20full%20algorithm%20considering%20the%20likely%20rotation%20of%20the%20load%20cell%20during%20a%20re-calibration.%0A%0A%20%20%20%20The%20idea%20is%20to%20guess%20initial%20values%20for%20the%20angles%2C%20estmate%20the%20re-calibration%20matrix%2C%20estimate%20new%20values%20for%20the%20angles%20that%20minimize%20the%20equation%20for%20the%20residuals%20and%20then%20estimate%20again%20the%20re-calibration%20matrix%20in%20an%20iterative%20approach%20until%20the%20estimated%20angles%20converge%20to%20the%20actual%20angles%20of%20the%20load%20cell%20in%20the%20different%20sites.%20This%20is%20a%20typical%20problem%20of%20%5Boptimization%5D(https%3A%2F%2Fen.wikipedia.org%2Fwiki%2FMathematical_optimization)%20where%20the%20angles%20are%20the%20design%20variables%20and%20the%20equation%20for%20the%20residuals%20is%20the%20cost%20function%20(see%20this%20%5Bnotebook%20about%20optimization%5D(https%3A%2F%2Fmolab.marimo.io%2Fgithub%2FBMClab%2FBMC%2Fblob%2Fmaster%2Fnotebooks%2Foptimization.py)).%0A%0A%20%20%20%20Let's%20code%20the%20optimization%20in%20a%20complete%20function%20for%20the%20force%20plate%20re-calibration%2C%20named%20%60fpcalibra.py%60%2C%20with%20the%20following%20signature%3A%0A%20%20%20%20%60%60%60python%0A%20%20%20%20%20%20%20C%2C%20ang%20%3D%20fpcalibra(Lfp%2C%20Flc%2C%20COP%2C%20threshold%3D1e-10)%0A%20%20%20%20%60%60%60%0A%0A%20%20%20%20Let's%20run%20the%20embedded%2C%20corrected%20implementation%20with%20the%20original%20simulation%3A%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(fpcalibra%2C%20inv%2C%20np)%3A%0A%20%20%20%20%23%20simulated%20true%20re-calibration%20matrix%0A%20%20%20%20_C%20%3D%20np.array(%5B%5B%201.0354%2C%20-0.0053%2C%20-0.0021%2C%20-0.0289%2C%20-0.0402%2C%20%200.0081%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5B%200.0064%2C%20%201.0309%2C%20-0.0031%2C%20%200.0211%2C%20%200.0135%2C%20-0.0001%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5B%200.0000%2C%20-0.0004%2C%20%201.0022%2C%20-0.0005%2C%20-0.0182%2C%20%200.0300%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5B-0.0012%2C%20-0.0385%2C%20%200.0002%2C%20%200.9328%2C%20%200.0007%2C%20%200.0017%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5B%200.0347%2C%20%200.0003%2C%20%200.0008%2C%20-0.0002%2C%20%200.9325%2C%20-0.0024%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5B-0.0004%2C%20-0.0013%2C%20-0.0003%2C%20-0.0023%2C%20%200.0035%2C%20%201.0592%5D%5D)%0A%20%20%20%20%23%20simulated%205%20measurements%20sites%20(in%20m)%0A%20%20%20%20_COP%20%3D%20np.array(%5B%5B%20%20%200%2C%20%20112%2C%20%20112%2C%20-112%2C%20-112%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5B%20%20%200%2C%20%20192%2C%20-192%2C%20%20192%2C%20-192%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5B-124%2C%20-124%2C%20-124%2C%20-124%2C%20-124%5D%5D)%2F1000%0A%20%20%20%20_nk%20%3D%20_COP.shape%5B1%5D%0A%20%20%20%20%23%20simulated%20forces%20measured%20by%20the%20load%20cell%20(in%20N)%20before%20rotation%0A%20%20%20%20_samples%20%3D%20np.linspace(1%2C%206000%2C%206000)%0A%20%20%20%20_ns%20%3D%20_samples.shape%5B0%5D%0A%20%20%20%20_Flc%20%3D%20np.empty((3%2C%20_nk*_ns))%0A%20%20%20%20_rng%20%3D%20np.random.default_rng(0)%0A%20%20%20%20for%20_k%20in%20range(_nk)%3A%0A%20%20%20%20%20%20%20%20_Flc%5B%3A%2C%20_k*_ns%3A(_k%2B1)*_ns%5D%20%3D%20np.array(%5B100*np.sin(5*2*np.pi*_samples%2F_ns)%20%2B%202*_rng.standard_normal(_ns)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20100*np.cos(5*2*np.pi*_samples%2F_ns)%20%2B%202*_rng.standard_normal(_ns)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20_samples%2F15%20%2B%20200%20%2B%205*_rng.standard_normal(_ns)%5D)%0A%20%20%20%20%23%20function%20for%20the%20COP%20skew-symmetric%20matrix%0A%20%20%20%20_Acop%20%3D%20lambda%20x%2Cy%2Cz%20%3A%20np.array(%5B%5B.0%2C%20-z%2C%20y%5D%2C%20%5Bz%2C%20.0%2C%20-x%5D%2C%20%5B-y%2C%20x%2C%20.0%5D%5D)%0A%20%20%20%20%23%20simulated%20loads%20measured%20by%20the%20force%20plate%0A%20%20%20%20_Li%20%3D%20np.empty((6%2C%20_ns*_nk))%0A%20%20%20%20_P%20%3D%20np.empty((6%2C%203%2C%20_nk))%0A%20%20%20%20for%20_k%2C%20_cop%20in%20enumerate(_COP.T)%3A%0A%20%20%20%20%20%20%20%20_P%5B%3A%2C%20%3A%2C%20_k%5D%20%3D%20np.vstack((np.eye(3)%2C%20_Acop(*_cop)))%0A%20%20%20%20%20%20%20%20_Li%5B%3A%2C%20_k*_ns%3A(_k%2B1)*_ns%5D%20%3D%20_P%5B%3A%2C%20%3A%2C%20_k%5D%20%40%20_Flc%5B%3A%2C%20_k*_ns%3A(_k%2B1)*_ns%5D%0A%20%20%20%20_Lfp%20%3D%20inv(_C)%20%40%20%20_Li%0A%20%20%20%20%23%20simulated%20angles%20of%20rotaton%20of%20the%20measurement%20sites%0A%20%20%20%20_ang%20%3D%20np.array(%5B20%2C%20-10%2C%200%2C%2015%2C%20-5%5D)%2F180*np.pi%0A%20%20%20%20%23%20function%20for%20the%20rotation%20matrix%0A%20%20%20%20_R%20%3D%20lambda%20a%20%3A%20np.array(%5B%5Bnp.cos(a)%2C%20-np.sin(a)%2C%200%5D%2C%20%5Bnp.sin(a)%2C%20np.cos(a)%2C%200%5D%2C%20%5B%200%2C%200%2C%201%5D%5D)%0A%20%20%20%20%23%20simulated%20forces%20measured%20by%20the%20load%20cell%20after%20rotation%0A%20%20%20%20for%20_k%20in%20range(_nk)%3A%0A%20%20%20%20%20%20%20%20_Flc%5B%3A%2C%20_k*_ns%3A(_k%2B1)*_ns%5D%20%3D%20_R(_ang%5B_k%5D).T%20%40%20_Flc%5B%3A%2C%20_k*_ns%3A(_k%2B1)*_ns%5D%0A%0A%20%20%20%20_C2%2C%20_ang2%20%3D%20fpcalibra(_Lfp%2C%20_Flc%2C%20_COP)%0A%0A%20%20%20%20_e%20%3D%20np.sqrt(np.sum((_C2-_C)**2))%0A%20%20%20%20print('Residual%20between%20simulated%20and%20optimal%20re-calibration%20matrices%3A'%2C%20_e)%0A%20%20%20%20_e%20%3D%20np.sqrt(np.sum((_ang2-_ang)**2))%0A%20%20%20%20print('Residual%20between%20simulated%20and%20optimal%20rotation%20angles%3A'%2C%20_e)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20The%20simulation%20works%20as%20expected%20and%20the%20function%20was%20able%20to%20estimate%20accurately%20the%20known%20initial%20re-calibration%20matrix%20and%20angles%20of%20rotation.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Exercise%203%20%E2%80%94%20convergence.**%20Angle%20changes%20at%20five%20sites%20are%20%24(0%2C0.02%2C0.01%2C0.03%2C0.01)%24%20rad%20and%20the%20stopping%20tolerance%20is%20%2410%5E%7B-6%7D%24%20rad.%20Has%20the%20alternating%20fit%20converged%3F%20Which%20condition%20should%20continue%20the%20loop%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Worked%20solution%20%E2%80%94%20Exercise%203%20%E2%80%94%20convergence%0A%0A%20%20%20%20No%3A%20four%20sites%20still%20change.%20Continue%20while%20**any**%20absolute%20change%20exceeds%20tolerance%2C%20%60np.any(delta_ang%20%3E%20threshold)%60.%20Using%20%60all%60%20would%20stop%20as%20soon%20as%20one%20site%20becomes%20stationary.%20Recompute%20the%20calibration%20using%20the%20final%20angles%20so%20both%20returned%20quantities%20belong%20to%20the%20same%20final%20iterate.%0A%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Non-linear%20algorithm%20for%20force%20plate%20calibration%0A%0A%20%20%20%20Cappello%20et%20al.%20(2011)%20extended%20the%20algorithm%20described%20earlier%20and%20proposed%20an%20algorithm%20for%20non-linear%20re-calibration%20of%20FPs.%3Cbr%3E%0A%20%20%20%20The%20idea%20is%20that%20a%20load%20applied%20on%20the%20FP%20produces%20bending%20which%20depends%20on%20the%20point%20of%20force%20application%20and%20in%20turn%20will%20result%20in%20systematic%20errors%20in%20the%20COP%20determination.%20Consequently%2C%20this%20non-linearity%20could%20be%20modeled%20and%20compensated%20with%20a%20re-calibration%20which%20takes%20into%20account%20the%20COP%20coordinates%20measured%20by%20the%20FP%20and%20added%20to%20the%20linear%20re-calibration%20we%20deduced%20above%20(Cappello%20et%20al.%2C%202011).%3Cbr%3E%0A%20%20%20%20The%20re-calibration%20equation%20will%20be%20(Cappello%20et%20al.%2C%202011)%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cbegin%7Barray%7D%7Bl%20l%7D%0A%20%20%20%20%5Cmathbf%7BL%7D_C%20%3D%20%5Cmathbf%7BC%7D_0%5Cmathbf%7BL%7D%20%2B%20%5C%2C%0A%20%20%20%20%5Cbegin%7Bbmatrix%7D%0A%20%20%20%20C_%7Bx_%7B11%7D%7D%20%26%20C_%7Bx_%7B12%7D%7D%20%26%200%20%26%20C_%7Bx_%7B14%7D%7D%20%26%20C_%7Bx_%7B15%7D%7D%20%26%20C_%7Bx_%7B16%7D%7D%20%5C%5C%0A%20%20%20%20C_%7Bx_%7B21%7D%7D%20%26%20C_%7Bx_%7B22%7D%7D%20%26%200%20%26%20C_%7Bx_%7B24%7D%7D%20%26%20C_%7Bx_%7B25%7D%7D%20%26%20C_%7Bx_%7B26%7D%7D%20%5C%5C%0A%20%20%20%20C_%7Bx_%7B31%7D%7D%20%26%20C_%7Bx_%7B32%7D%7D%20%26%200%20%26%20C_%7Bx_%7B34%7D%7D%20%26%20C_%7Bx_%7B35%7D%7D%20%26%20C_%7Bx_%7B36%7D%7D%20%5C%5C%0A%20%20%20%20C_%7Bx_%7B41%7D%7D%20%26%20C_%7Bx_%7B42%7D%7D%20%26%200%20%26%20C_%7Bx_%7B44%7D%7D%20%26%20C_%7Bx_%7B45%7D%7D%20%26%20C_%7Bx_%7B46%7D%7D%20%5C%5C%0A%20%20%20%20C_%7Bx_%7B51%7D%7D%20%26%20C_%7Bx_%7B52%7D%7D%20%26%200%20%26%20C_%7Bx_%7B54%7D%7D%20%26%20C_%7Bx_%7B55%7D%7D%20%26%20C_%7Bx_%7B56%7D%7D%20%5C%5C%0A%20%20%20%20C_%7Bx_%7B61%7D%7D%20%26%20C_%7Bx_%7B62%7D%7D%20%26%200%20%26%20C_%7Bx_%7B64%7D%7D%20%26%20C_%7Bx_%7B65%7D%7D%20%26%20C_%7Bx_%7B66%7D%7D%0A%20%20%20%20%5Cend%7Bbmatrix%7D%5C%2C%0A%20%20%20%20%5Cbegin%7Bbmatrix%7D%0A%20%20%20%20F_x%20%5C%5C%20F_y%20%5C%5C%20F_z%20%5C%5C%20M_x%20%5C%5C%20M_y%20%5C%5C%20M_z%0A%20%20%20%20%5Cend%7Bbmatrix%7D%20COP_x%20%2B%0A%20%20%20%20%5Cbegin%7Bbmatrix%7D%0A%20%20%20%20C_%7By_%7B11%7D%7D%20%26%20C_%7By_%7B12%7D%7D%20%26%200%20%26%200%20%26%20C_%7By_%7B15%7D%7D%20%26%20C_%7By_%7B16%7D%7D%20%5C%5C%0A%20%20%20%20C_%7By_%7B21%7D%7D%20%26%20C_%7By_%7B22%7D%7D%20%26%200%20%26%200%20%26%20C_%7By_%7B25%7D%7D%20%26%20C_%7By_%7B26%7D%7D%20%5C%5C%0A%20%20%20%20C_%7By_%7B31%7D%7D%20%26%20C_%7By_%7B32%7D%7D%20%26%200%20%26%200%20%26%20C_%7By_%7B35%7D%7D%20%26%20C_%7By_%7B36%7D%7D%20%5C%5C%0A%20%20%20%20C_%7By_%7B41%7D%7D%20%26%20C_%7By_%7B42%7D%7D%20%26%200%20%26%200%20%26%20C_%7By_%7B45%7D%7D%20%26%20C_%7By_%7B46%7D%7D%20%5C%5C%0A%20%20%20%20C_%7By_%7B51%7D%7D%20%26%20C_%7By_%7B52%7D%7D%20%26%200%20%26%200%20%26%20C_%7By_%7B55%7D%7D%20%26%20C_%7By_%7B56%7D%7D%20%5C%5C%0A%20%20%20%20C_%7By_%7B61%7D%7D%20%26%20C_%7By_%7B62%7D%7D%20%26%200%20%26%200%20%26%20C_%7By_%7B65%7D%7D%20%26%20C_%7By_%7B66%7D%7D%0A%20%20%20%20%5Cend%7Bbmatrix%7D%5C%2C%0A%20%20%20%20%5Cbegin%7Bbmatrix%7D%0A%20%20%20%20F_x%20%5C%5C%20F_y%20%5C%5C%20F_z%20%5C%5C%20M_x%20%5C%5C%20M_y%20%5C%5C%20M_z%0A%20%20%20%20%5Cend%7Bbmatrix%7D%20COP_y%0A%20%20%20%20%5C%5C%5B6pt%5D%0A%20%20%20%20%5Cmathbf%7BL%7D_C%20%3D%20(%5Cmathbf%7BC%7D_0%20%2B%20%5Cmathbf%7BC%7D_x%20COP_x%20%2B%20%5Cmathbf%7BC%7D_y%20COP_y)%5Cmathbf%7BL%7D%20%3D%20%5Cmathbf%7BC%7D_%7BNL%7D%5Cmathbf%7BL%7D%0A%20%20%20%20%5Cend%7Barray%7D%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Where%20%24%5Cmathbf%7BC%7D_0%24%20is%20the%20linear%20re-calibration%20matrix%2C%20%24%5Cmathbf%7BL%7D%24%20is%20the%20measured%20FP%20output%2C%20%24%5Cmathbf%7BC%7D_x%24%20and%20%20%24%5Cmathbf%7BC%7D_y%24%20are%20the%20non-linear%20re-calibration%20matrices.%0A%0A%20%20%20%20To%20estimate%20%24%5Cmathbf%7BC%7D_%7BNL%7D%24%2C%20Cappello%20et%20al.%20(2011)%20suggest%20to%20employ%20the%20algorithm%20proposed%20by%20Cedraro%20et%20al.%20(2008)%20to%20estimate%20the%20linear%20re-calibration%20described%20earlier.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20The%20position-dependent%20model%20below%20is%20retained%20as%20theory%20from%20Cappello%20et%20al.%20The%20following%20%60_fpcalibra_nl%60%20block%20is%20an%20**unfinished%20historical%20draft**%3A%20its%20body%20repeats%20the%20linear%20algorithm%20and%20does%20not%20fit%20the%20extra%20%24C_x%2CC_y%24%20terms.%20It%20is%20included%20for%20reading%20and%20extension%2C%20not%20called%20as%20a%20completed%20nonlinear%20calibration.%20A%20complete%20implementation%20needs%20the%20expanded%20design%2C%20constraints%2C%20identifiability%20checks%20and%20independent%20validation.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(minimize%2C%20np%2C%20pinv%2C%20time)%3A%0A%20%20%20%20%23%20Draft%20(never%20run%20in%20the%20original%20notebook)%3A%20the%20body%20of%20the%20linear%20algorithm%20of%0A%20%20%20%20%23%20fpcalibra%2C%20as%20a%20starting%20point%20for%20implementing%20the%20non-linear%20re-calibration%0A%20%20%20%20def%20_fpcalibra_nl(Lfp%2C%20Flc%2C%20COP%2C%20threshold%3D1e-10%2C%20method%3D'SVD')%3A%0A%20%20%20%20%20%20%20%20%23%20number%20of%20sites%0A%20%20%20%20%20%20%20%20nk%20%3D%20COP.shape%5B1%5D%0A%20%20%20%20%20%20%20%20%23%20number%20of%20samples%0A%20%20%20%20%20%20%20%20ns%20%3D%20int(Lfp.shape%5B1%5D%2Fnk)%0A%20%20%20%20%20%20%20%20%23%20function%20for%20the%20COP%20skew-symmetric%20matrix%0A%20%20%20%20%20%20%20%20Acop%20%3D%20lambda%20x%2Cy%2Cz%20%3A%20np.array(%5B%5B.0%2C%20-z%2C%20y%5D%2C%20%5Bz%2C%20.0%2C%20-x%5D%2C%20%5B-y%2C%20x%2C%20.0%5D%5D)%0A%20%20%20%20%20%20%20%20P%20%3D%20np.empty((6%2C%203%2C%20nk))%0A%20%20%20%20%20%20%20%20for%20k%2C%20cop%20in%20enumerate(COP.T)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20P%5B%3A%2C%20%3A%2C%20k%5D%20%3D%20np.vstack((np.eye(3)%2C%20Acop(*cop)))%0A%20%20%20%20%20%20%20%20%23%20function%20for%20the%202D%20rotation%20matrix%0A%20%20%20%20%20%20%20%20R%20%3D%20lambda%20a%20%3A%20np.array(%5B%5Bnp.cos(a)%2C%20-np.sin(a)%2C%200%5D%2C%20%5Bnp.sin(a)%2C%20np.cos(a)%2C%200%5D%2C%20%5B%200%2C%200%2C%201%5D%5D)%0A%20%20%20%20%20%20%20%20%23%20Pseudoiverse%20of%20the%20loads%20measured%20by%20the%20force%20plate%0A%20%20%20%20%20%20%20%20if%20method.lower()%20%3D%3D%20'svd'%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20Lpinv%20%3D%20pinv(Lfp)%0A%20%20%20%20%20%20%20%20else%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20Lpinv%20%3D%20pinv(Lfp)%0A%20%20%20%20%20%20%20%20%23%20cost%20function%20for%20the%20optimization%0A%20%20%20%20%20%20%20%20def%20costfun(ang%2C%20P%2C%20R%2C%20Flc%2C%20CLfp%2C%20nk%2C%20ns%2C%20E)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20for%20k%20in%20range(nk)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20E%5B%3A%2Ck*ns%3A(k%2B1)*ns%5D%20%3D%20(P%5B%3A%2C%3A%2Ck%5D%20%40%20R(ang%5Bk%5D))%20%40%20Flc%5B%3A%2Ck*ns%3A(k%2B1)*ns%5D%20-%20CLfp%5B%3A%2Ck*ns%3A(k%2B1)*ns%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20return%20np.sum(E%20*%20E)%0A%20%20%20%20%20%20%20%20%23%20inequality%20constraints%0A%20%20%20%20%20%20%20%20bnds%20%3D%20%5B(-np.pi%2F2%2C%20np.pi%2F2)%20for%20k%20in%20range(nk)%5D%0A%20%20%20%20%20%20%20%20%23%20some%20initialization%0A%20%20%20%20%20%20%20%20ang0%20%3D%20np.zeros(nk)%0A%20%20%20%20%20%20%20%20E%20%3D%20np.empty((6%2C%20ns*nk))%0A%20%20%20%20%20%20%20%20da%20%3D%20%5B%5D%0A%20%20%20%20%20%20%20%20delta_ang%20%3D%2010*threshold%0A%20%20%20%20%20%20%20%20Li%20%3D%20np.empty((6%2C%20ns*nk))%0A%20%20%20%20%20%20%20%20start%20%3D%20time.time()%0A%20%20%20%20%20%20%20%20%23%20the%20optimization%0A%20%20%20%20%20%20%20%20while%20np.any(delta_ang%20%3E%20threshold)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20for%20k%20in%20range(nk)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20Li%5B%3A%2Ck*ns%3A(k%2B1)*ns%5D%20%3D%20(P%5B%3A%2C%3A%2Ck%5D%20%40%20R(ang0%5Bk%5D))%20%40%20Flc%5B%3A%2Ck*ns%3A(k%2B1)*ns%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20C%20%3D%20Li%20%40%20Lpinv%0A%20%20%20%20%20%20%20%20%20%20%20%20CLfp%20%3D%20C%20%40%20Lfp%0A%20%20%20%20%20%20%20%20%20%20%20%20res%20%3D%20minimize(fun%3Dcostfun%2C%20x0%3Dang0%2C%20args%3D(P%2C%20R%2C%20Flc%2C%20CLfp%2C%20nk%2C%20ns%2C%20E)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20bounds%3Dbnds%2C%20method%3D'TNC'%2C%20options%3D%7B'disp'%3A%20False%7D)%0A%20%20%20%20%20%20%20%20%20%20%20%20delta_ang%20%3D%20np.abs(res.x%20-%20ang0)%0A%20%20%20%20%20%20%20%20%20%20%20%20ang0%20%3D%20res.x%0A%20%20%20%20%20%20%20%20%20%20%20%20da.append(delta_ang.sum())%0A%0A%20%20%20%20%20%20%20%20tdelta%20%3D%20time.time()%20-%20start%0A%20%20%20%20%20%20%20%20print('%5CnOptimization%20finished%20in%20%25.1f%20s%20after%20%25d%20steps.%5Cn'%20%25(tdelta%2C%20len(da)))%0A%20%20%20%20%20%20%20%20print('Optimal%20calibration%20matrix%3A%5Cn'%2C%20C)%0A%20%20%20%20%20%20%20%20print('%5CnOptimal%20angles%3A%5Cn'%2C%20res.x*180%2Fnp.pi)%0A%20%20%20%20%20%20%20%20print('%5Cn')%0A%0A%20%20%20%20%20%20%20%20return%20C%2C%20res.x%0A%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Function%20fpcalibra.py%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(np)%3A%0A%20%20%20%20%22%22%22Force%20plate%20calibration%20algorithm.%0A%20%20%20%20%22%22%22%0A%0A%20%20%20%20__author__%20%3D%20'Marcos%20Duarte%2C%20https%3A%2F%2Fgithub.com%2Fdemotu%2FBMC'%0A%20%20%20%20__version__%20%3D%20'fpcalibra.py%20v.1.0.1%202016%2F08%2F19'%0A%20%20%20%20__license__%20%3D%20%22MIT%22%0A%0A%20%20%20%20from%20scipy.linalg%20import%20pinv%0A%20%20%20%20from%20scipy.optimize%20import%20minimize%0A%20%20%20%20import%20time%0A%0A%0A%20%20%20%20def%20fpcalibra(Lfp%2C%20Flc%2C%20COP%2C%20threshold%3D1e-10%2C%20method%3D'SVD')%3A%0A%20%20%20%20%20%20%20%20%22%22%22Force%20plate%20calibration%20algorithm.%0A%0A%20%20%20%20%20%20%20%20For%20a%20force%20plate%20(FP)%20re-calibration%2C%20the%20relationship%20between%20the%0A%20%20%20%20%20%20%20%20measured%20FP%20output%20(L)%20and%20the%20known%20loads%20(Li)%20is%20approximated%20by%3A%0A%20%20%20%20%20%20%20%20Li%20%3D%20C%40L%20%2B%20E%20(%40%20is%20the%20operator%20for%20matrix%20multiplication).%0A%20%20%20%20%20%20%20%20Where%20C%20is%20the%206-by-6%20re-calibration%20matrix%20and%20E%20is%20a%20gaussian%2C%0A%20%20%20%20%20%20%20%20uncorrelated%2C%20zero%20mean%20noise%20six-by-one%20matrix.%0A%0A%20%20%20%20%20%20%20%20The%20re-calibration%20matrix%20can%20be%20found%20by%20solving%20the%20equation%20above%20and%0A%20%20%20%20%20%20%20%20then%20C%20can%20be%20later%20used%20to%20re-calibrate%20the%20FP%20output%3A%20Lc%20%3D%20C%40L.%0A%20%20%20%20%20%20%20%20Where%20Lc%20is%20the%20re-calibrated%20FP%20output.%0A%0A%20%20%20%20%20%20%20%20Cedraro%20et%20al.%20(2008)%20%5B1%5D_%20proposed%20to%20use%20a%20calibrated%20three-component%0A%20%20%20%20%20%20%20%20load%20cell%20to%20measure%20the%20forces%20applied%20on%20the%20FP%20at%20known%20measurement%0A%20%20%20%20%20%20%20%20sites%20and%20an%20algorithm%20for%20the%20re-calibration.%0A%0A%20%20%20%20%20%20%20%20This%20code%20implements%20the%20re-calibration%20algorithm%2C%20see%20%5B2%5D_%0A%0A%20%20%20%20%20%20%20%20Parameters%0A%20%20%20%20%20%20%20%20----------%0A%20%20%20%20%20%20%20%20Lfp%20%3A%20numpy%202-D%20array%20(6%2C%20nsamples*nksites)%0A%20%20%20%20%20%20%20%20%20%20%20%20loads%20%5BFx%2C%20Fy%2C%20Fz%2C%20Mx%2C%20My%2C%20Mz%5D%20(in%20N%20and%20Nm)%20measured%20by%20the%20force%0A%20%20%20%20%20%20%20%20%20%20%20%20plate%20due%20to%20the%20corresponding%20forces%20applied%20at%20the%20measurement%20sites%0A%20%20%20%20%20%20%20%20Flc%20%3A%20numpy%202-D%20array%20(3%2C%20nsamples*nksites)%0A%20%20%20%20%20%20%20%20%20%20%20%20forces%20%5BFx%2C%20Fy%2C%20Fz%5D%20(in%20N)%20measured%20by%20the%20load%20cell%20at%20the%0A%20%20%20%20%20%20%20%20%20%20%20%20measurement%20sites%0A%20%20%20%20%20%20%20%20COP%20%3A%20numpy%202-D%20array%20(3%2C%20nksites)%0A%20%20%20%20%20%20%20%20%20%20%20%20positions%20%5BCOPx%2C%20COPy%2C%20COPz%5D%20(in%20m)%20of%20the%20load%20cell%20at%20the%0A%20%20%20%20%20%20%20%20%20%20%20%20measurement%20sites%0A%20%20%20%20%20%20%20%20threshold%20%20%3A%20float%2C%20optional%0A%20%20%20%20%20%20%20%20%20%20%20%20threshold%20to%20stop%20the%20optimization%20(default%201e-10)%0A%20%20%20%20%20%20%20%20method%20%20%3A%20string%2C%20optional%0A%20%20%20%20%20%20%20%20%20%20%20%20legacy%20pseudoinverse%20selector%2C%20'SVD'%20(default)%20or%20'lstsq'.%20Both%20values%0A%20%20%20%20%20%20%20%20%20%20%20%20now%20use%20scipy.linalg.pinv%2C%20which%20is%20SVD-based.%0A%0A%20%20%20%20%20%20%20%20Returns%0A%20%20%20%20%20%20%20%20-------%0A%20%20%20%20%20%20%20%20C%20%20%20%3A%20numpy%202-D%20(6-by-6)%20array%0A%20%20%20%20%20%20%20%20%20%20%20%20force%20plate%20re-calibration%20matrix.%20Force%2Fforce%20and%20moment%2Fmoment%20blocks%0A%20%20%20%20%20%20%20%20%20%20%20%20are%20dimensionless%3B%20force%2Fmoment%20and%20moment%2Fforce%20blocks%20have%20units%201%2Fm%0A%20%20%20%20%20%20%20%20%20%20%20%20and%20m%2C%20respectively%2C%20when%20loads%20use%20N%20and%20Nm.%0A%20%20%20%20%20%20%20%20ang%20%3A%20numpy%201-D%20array%20%5Bang0%2C%20...%2C%20angk%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20optimal%20angles%20of%20rotation%20(in%20rad)%20of%20the%20load%20cells%20at%20the%0A%20%20%20%20%20%20%20%20%20%20%20%20measurement%20sites%0A%0A%20%20%20%20%20%20%20%20References%0A%20%20%20%20%20%20%20%20----------%0A%20%20%20%20%20%20%20%20..%20%5B1%5D%20Cedraro%20A%2C%20Cappello%20A%2C%20Chiari%20L%20(2008)%20Gait%20%26%20Posture%2C%2028%2C%20488%E2%80%93494.%0A%20%20%20%20%20%20%20%20..%20%5B2%5D%20force_plate_calibration.py%0A%0A%20%20%20%20%20%20%20%20Example%0A%20%20%20%20%20%20%20%20-------%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20from%20fpcalibra%20import%20fpcalibra%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%20%20%20%20%3E%3E%3E%20from%20numpy.linalg%20import%20inv%0A%20%20%20%20%20%20%20%20%3E%3E%3E%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%23%20simulated%20true%20re-calibration%20matrix%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20C%20%3D%20np.array(%5B%5B%201.0354%2C%20-0.0053%2C%20-0.0021%2C%20-0.0289%2C%20-0.0402%2C%20%200.0081%5D%2C%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5B%200.0064%2C%20%201.0309%2C%20-0.0031%2C%20%200.0211%2C%20%200.0135%2C%20-0.0001%5D%2C%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5B%200.0000%2C%20-0.0004%2C%20%201.0022%2C%20-0.0005%2C%20-0.0182%2C%20%200.0300%5D%2C%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5B-0.0012%2C%20-0.0385%2C%20%200.0002%2C%20%200.9328%2C%20%200.0007%2C%20%200.0017%5D%2C%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5B%200.0347%2C%20%200.0003%2C%20%200.0008%2C%20-0.0002%2C%20%200.9325%2C%20-0.0024%5D%2C%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5B-0.0004%2C%20-0.0013%2C%20-0.0003%2C%20-0.0023%2C%20%200.0035%2C%20%201.0592%5D%5D)%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%23%20simulated%205%20measurements%20sites%20(in%20m)%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20COP%20%3D%20np.array(%5B%5B%20%20%200%2C%20%20112%2C%20%20112%2C%20-112%2C%20-112%5D%2C%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5B%20%20%200%2C%20%20192%2C%20-192%2C%20%20192%2C%20-192%5D%2C%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5B-124%2C%20-124%2C%20-124%2C%20-124%2C%20-124%5D%5D)%2F1000%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20nk%20%3D%20COP.shape%5B1%5D%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%23%20simulated%20forces%20measured%20by%20the%20load%20cell%20(in%20N)%20before%20rotation%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20samples%20%3D%20np.linspace(1%2C%206000%2C%206000)%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20ns%20%3D%20samples.shape%5B0%5D%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20Flc%20%3D%20np.empty((3%2C%20nk*ns))%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20for%20k%20in%20range(nk)%3A%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%20%20%20%20Flc%5B%3A%2C%20k*ns%3A(k%2B1)*ns%5D%20%3D%20np.array(%5B100*np.sin(5*2*np.pi*samples%2Fns)%20%2B%202*np.random.randn(ns)%2C%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20100*np.cos(5*2*np.pi*samples%2Fns)%20%2B%202*np.random.randn(ns)%2C%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20samples%2F15%20%2B%20200%20%2B%205*np.random.randn(ns)%5D)%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%23%20function%20for%20the%20COP%20skew-symmetric%20matrix%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20Acop%20%3D%20lambda%20x%2Cy%2Cz%20%3A%20np.array(%5B%5B.0%2C%20-z%2C%20y%5D%2C%20%5Bz%2C%20.0%2C%20-x%5D%2C%20%5B-y%2C%20x%2C%20.0%5D%5D)%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%23%20simulated%20loads%20measured%20by%20the%20force%20plate%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20Li%20%3D%20np.empty((6%2C%20ns*nk))%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20P%20%3D%20np.empty((6%2C%203%2C%20nk))%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20for%20k%2C%20cop%20in%20enumerate(COP.T)%3A%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%20%20%20%20P%5B%3A%2C%20%3A%2C%20k%5D%20%3D%20np.vstack((np.eye(3)%2C%20Acop(*cop)))%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%20%20%20%20Li%5B%3A%2C%20k*ns%3A(k%2B1)*ns%5D%20%3D%20P%5B%3A%2C%20%3A%2C%20k%5D%20%40%20Flc%5B%3A%2C%20k*ns%3A(k%2B1)*ns%5D%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20Lfp%20%3D%20inv(C)%20%40%20%20Li%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%23%20simulated%20angles%20of%20rotaton%20of%20the%20measurement%20sites%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20ang%20%3D%20np.array(%5B20%2C%20-10%2C%200%2C%2015%2C%20-5%5D)%2F180*np.pi%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%23%20function%20for%20the%20rotation%20matrix%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20R%20%3D%20lambda%20a%20%3A%20np.array(%5B%5Bnp.cos(a)%2C%20-np.sin(a)%2C%200%5D%2C%20%5Bnp.sin(a)%2C%20np.cos(a)%2C%200%5D%2C%20%5B%200%2C%200%2C%201%5D%5D)%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%23%20simulated%20forces%20measured%20by%20the%20load%20cell%20after%20rotation%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20for%20k%20in%20range(nk)%3A%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%20%20%20%20Flc%5B%3A%2C%20k*ns%3A(k%2B1)*ns%5D%20%3D%20R(ang%5Bk%5D).T%20%40%20Flc%5B%3A%2C%20k*ns%3A(k%2B1)*ns%5D%0A%20%20%20%20%20%20%20%20%3E%3E%3E%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20C2%2C%20ang2%20%3D%20fpcalibra(Lfp%2C%20Flc%2C%20COP)%0A%20%20%20%20%20%20%20%20%3E%3E%3E%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20e%20%3D%20np.sqrt(np.sum((C2-C)**2))%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20print('Residual%20between%20simulated%20and%20optimal%20re-calibration%20matrices%3A'%2C%20e)%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20e%20%3D%20np.sqrt(np.sum((ang2-ang)**2))%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20print('Residual%20between%20simulated%20and%20optimal%20rotation%20angles%3A'%2C%20e)%0A%20%20%20%20%20%20%20%20%22%22%22%0A%0A%20%20%20%20%20%20%20%20%23%20number%20of%20sites%0A%20%20%20%20%20%20%20%20nk%20%3D%20COP.shape%5B1%5D%0A%20%20%20%20%20%20%20%20%23%20number%20of%20samples%0A%20%20%20%20%20%20%20%20ns%20%3D%20int(Lfp.shape%5B1%5D%2Fnk)%0A%20%20%20%20%20%20%20%20%23%20function%20for%20the%20COP%20skew-symmetric%20matrix%0A%20%20%20%20%20%20%20%20Acop%20%3D%20lambda%20x%2Cy%2Cz%20%3A%20np.array(%5B%5B.0%2C%20-z%2C%20y%5D%2C%20%5Bz%2C%20.0%2C%20-x%5D%2C%20%5B-y%2C%20x%2C%20.0%5D%5D)%0A%20%20%20%20%20%20%20%20P%20%3D%20np.empty((6%2C%203%2C%20nk))%0A%20%20%20%20%20%20%20%20for%20k%2C%20cop%20in%20enumerate(COP.T)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20P%5B%3A%2C%20%3A%2C%20k%5D%20%3D%20np.vstack((np.eye(3)%2C%20Acop(*cop)))%0A%20%20%20%20%20%20%20%20%23%20function%20for%20the%202D%20rotation%20matrix%0A%20%20%20%20%20%20%20%20R%20%3D%20lambda%20a%20%3A%20np.array(%5B%5Bnp.cos(a)%2C%20-np.sin(a)%2C%200%5D%2C%20%5Bnp.sin(a)%2C%20np.cos(a)%2C%200%5D%2C%20%5B%200%2C%200%2C%201%5D%5D)%0A%20%20%20%20%20%20%20%20%23%20Pseudoiverse%20of%20the%20loads%20measured%20by%20the%20force%20plate%0A%20%20%20%20%20%20%20%20if%20method.lower()%20%3D%3D%20'svd'%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20Lpinv%20%3D%20pinv(Lfp)%20%20%23%20pinv2%20was%20removed%20in%20SciPy%201.9%3B%20pinv%20now%20also%20uses%20the%20SVD%0A%20%20%20%20%20%20%20%20else%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20Lpinv%20%3D%20pinv(Lfp)%0A%20%20%20%20%20%20%20%20%23%20cost%20function%20for%20the%20optimization%0A%20%20%20%20%20%20%20%20def%20costfun(ang%2C%20P%2C%20R%2C%20Flc%2C%20CLfp%2C%20nk%2C%20ns%2C%20E)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20for%20k%20in%20range(nk)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20E%5B%3A%2Ck*ns%3A(k%2B1)*ns%5D%20%3D%20(P%5B%3A%2C%3A%2Ck%5D%20%40%20R(ang%5Bk%5D))%20%40%20Flc%5B%3A%2Ck*ns%3A(k%2B1)*ns%5D%20-%20CLfp%5B%3A%2Ck*ns%3A(k%2B1)*ns%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20return%20np.sum(E%20*%20E)%0A%20%20%20%20%20%20%20%20%23%20inequality%20constraints%0A%20%20%20%20%20%20%20%20bnds%20%3D%20%5B(-np.pi%2F2%2C%20np.pi%2F2)%20for%20k%20in%20range(nk)%5D%0A%20%20%20%20%20%20%20%20%23%20some%20initialization%0A%20%20%20%20%20%20%20%20ang0%20%3D%20np.zeros(nk)%0A%20%20%20%20%20%20%20%20E%20%3D%20np.empty((6%2C%20ns*nk))%0A%20%20%20%20%20%20%20%20da%20%3D%20%5B%5D%0A%20%20%20%20%20%20%20%20delta_ang%20%3D%2010*threshold%0A%20%20%20%20%20%20%20%20Li%20%3D%20np.empty((6%2C%20ns*nk))%0A%20%20%20%20%20%20%20%20start%20%3D%20time.time()%0A%20%20%20%20%20%20%20%20%23%20the%20optimization%0A%20%20%20%20%20%20%20%20while%20np.any(delta_ang%20%3E%20threshold)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20for%20k%20in%20range(nk)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20Li%5B%3A%2Ck*ns%3A(k%2B1)*ns%5D%20%3D%20(P%5B%3A%2C%3A%2Ck%5D%20%40%20R(ang0%5Bk%5D))%20%40%20Flc%5B%3A%2Ck*ns%3A(k%2B1)*ns%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20C%20%3D%20Li%20%40%20Lpinv%0A%20%20%20%20%20%20%20%20%20%20%20%20CLfp%20%3D%20C%20%40%20Lfp%0A%20%20%20%20%20%20%20%20%20%20%20%20res%20%3D%20minimize(fun%3Dcostfun%2C%20x0%3Dang0%2C%20args%3D(P%2C%20R%2C%20Flc%2C%20CLfp%2C%20nk%2C%20ns%2C%20E)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20bounds%3Dbnds%2C%20method%3D'TNC'%2C%20options%3D%7B'disp'%3A%20False%7D)%0A%20%20%20%20%20%20%20%20%20%20%20%20delta_ang%20%3D%20np.abs(res.x%20-%20ang0)%0A%20%20%20%20%20%20%20%20%20%20%20%20ang0%20%3D%20res.x%0A%20%20%20%20%20%20%20%20%20%20%20%20da.append(delta_ang.sum())%0A%0A%20%20%20%20%20%20%20%20%23%20Return%20a%20calibration%20consistent%20with%20the%20final%20angle%20estimate.%0A%20%20%20%20%20%20%20%20for%20k%20in%20range(nk)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20Li%5B%3A%2Ck*ns%3A(k%2B1)*ns%5D%20%3D%20(P%5B%3A%2C%3A%2Ck%5D%20%40%20R(ang0%5Bk%5D))%20%40%20Flc%5B%3A%2Ck*ns%3A(k%2B1)*ns%5D%0A%20%20%20%20%20%20%20%20C%20%3D%20Li%20%40%20Lpinv%0A%0A%20%20%20%20%20%20%20%20tdelta%20%3D%20time.time()%20-%20start%0A%20%20%20%20%20%20%20%20print('%5CnOptimization%20finished%20in%20%25.1f%20s%20after%20%25d%20steps.%5Cn'%20%25(tdelta%2C%20len(da)))%0A%20%20%20%20%20%20%20%20print('Optimal%20calibration%20matrix%3A%5Cn'%2C%20C)%0A%20%20%20%20%20%20%20%20print('%5CnOptimal%20angles%3A%5Cn'%2C%20res.x*180%2Fnp.pi)%0A%20%20%20%20%20%20%20%20print('%5Cn')%0A%0A%20%20%20%20%20%20%20%20return%20C%2C%20res.x%0A%0A%20%20%20%20return%20fpcalibra%2C%20minimize%2C%20pinv%2C%20time%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20import%20marimo%20as%20mo%0A%0A%20%20%20%20return%20(mo%2C)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Exercise%204%20%E2%80%94%20a%20position-dependent%20correction.**%20A%20single%20channel%20obeys%20%24F_c%3D(1%2B0.2x)F%24%2C%20with%20%24x%24%20in%20m.%20Evaluate%20it%20at%20%24x%3D0%24%20and%20%24x%3D0.1%24%20m%20for%20%24F%3D500%24%20N.%20What%20units%20must%200.2%20carry%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Worked%20solution%20%E2%80%94%20Exercise%204%0A%0A%20%20%20%20At%20the%20centre%20the%20result%20is%20500%20N.%20At%200.1%20m%20the%20gain%20is%201.02%20and%20the%20result%20is%20510%20N.%20The%20coefficient%20has%20units%20m%E2%81%BB%C2%B9%20so%20that%20its%20product%20with%20position%20is%20dimensionless.%20This%20example%20illustrates%20the%20extra%20position%20dependence%3B%20it%20is%20not%20a%20fit%20of%20the%20full%20Cappello%20model.%0A%0A%20%20%20%20%2F%2F%2F%0A%0A%20%20%20%20%23%23%20Checkpoint%20questions%0A%0A%20%20%20%20-%20Does%20a%20small%20residual%20establish%20that%20all%20calibration%20coefficients%20are%20identifiable%3F%0A%20%20%20%20-%20Why%20are%20load-cell%20rotation%20and%20calibration%20alternated%20in%20the%20example%3F%0A%20%20%20%20-%20What%20further%20evidence%20is%20needed%20before%20a%20more%20flexible%20nonlinear%20model%20is%20preferable%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Checkpoint%20answers%0A%0A%20%20%20%20No%3A%20inspect%20excitation%2C%20rank%20and%20conditioning%2C%20as%20well%20as%20residuals.%20The%20reference%20forces%20must%20be%20expressed%20in%20the%20plate%20frame%3B%20the%20unknown%20orientations%20and%20calibration%20affect%20each%20other's%20estimates.%20A%20nonlinear%20model%20needs%20sufficient%20spatial%2Floading%20coverage%20and%20improved%20performance%20on%20independent%20validation%20data%2C%20not%20just%20a%20smaller%20training%20residual.%0A%0A%20%20%20%20%2F%2F%2F%0A%0A%20%20%20%20%23%23%20Problems%0A%0A%20%20%20%201.%20A%20measured%20force%20is%20525%20N%20when%20the%20true%20force%20is%20500%20N.%20Find%20the%20scalar%20re-calibration%20gain.%0A%20%20%20%202.%20For%20%24%5Cmathbf%20r%3D(0.2%2C0.1%2C0)%24%20m%20and%20%24%5Cmathbf%20F%3D(10%2C20%2C300)%24%20N%2C%20give%20%24%5Cmathbf%20r%5Ctimes%5Cmathbf%20F%24.%0A%20%20%20%203.%20State%20the%20units%20of%20a%20coefficient%20multiplying%20%24F_z%24%20to%20predict%20%24M_x%24%2C%20and%20one%20multiplying%20%24M_y%24%20to%20predict%20%24F_x%24.%0A%20%20%20%204.%20A%20%246%5Ctimes1000%24%20load%20matrix%20has%20rank%204.%20Is%20an%20arbitrary%20six-by-six%20re-calibration%20unique%3F%20What%20can%20the%20pseudoinverse%20provide%3F%0A%20%20%20%205.%20Changes%20are%20%24(10%5E%7B-8%7D%2C10%5E%7B-3%7D)%24%20rad%20with%20tolerance%20%2410%5E%7B-6%7D%24%20rad.%20Stop%20or%20continue%3F%20Does%20the%20retained%20nonlinear%20draft%20estimate%20%24C_x%2CC_y%24%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Final%20answers%20to%20the%20problems%0A%0A%20%20%20%201.%20%2420%2F21%5Capprox0.952381%24.%0A%20%20%20%202.%20%24(30%2C-60%2C3)%24%20N%20m.%0A%20%20%20%203.%20m%3B%20m%E2%81%BB%C2%B9.%0A%20%20%20%204.%20No%3B%20a%20minimum-norm%20least-squares%20solution%2C%20with%20unidentifiable%20directions%20remaining.%0A%20%20%20%205.%20Continue%3B%20the%20draft%20does%20not%20estimate%20those%20matrices.%0A%0A%20%20%20%20%2F%2F%2F%0A%0A%20%20%20%20%23%23%20Go%20deeper%0A%0A%20%20%20%20Compare%20an%20independently%20loaded%20validation%20trial%20before%20and%20after%20calibration.%20Report%20force%20and%20moment%20errors%20separately%20and%20specify%20any%20weighting.%20For%20the%20nonlinear%20extension%2C%20examine%20whether%20your%20calibration%20sites%20cover%20the%20positions%20used%20in%20validation.%20Keep%20the%20original%20Cedraro%20and%20Cappello%20papers%20below%20as%20the%20algorithmic%20sources.%0A%0A%20%20%20%20%23%23%23%20Video%20lectures%20on%20the%20Internet%0A%0A%20%20%20%20-%20%5BMIT%2018.06SC%3A%20projection%20matrices%20and%20least%20squares%5D(https%3A%2F%2Focw.mit.edu%2Fcourses%2F18-06sc-linear-algebra-fall-2011%2Fpages%2Fleast-squares-determinants-and-eigenvalues%2Fprojection-matrices-and-least-squares%2F)%20%E2%80%94%20relate%20the%20projection%20view%20to%20the%20load-matrix%20pseudoinverse.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20References%0A%0A%20%20%20%20-%20%5BCedraro%20A%2C%20Cappello%20A%2C%20Chiari%20L%20(2008)%20A%20portable%20system%20for%20in-situ%20re-calibration%20of%20force%20platforms%3A%20theoretical%20validation.%20Gait%20%26%20Posture%2C%2028%2C%20488%E2%80%93494%5D(https%3A%2F%2Fpubmed.ncbi.nlm.nih.gov%2F18450453%2F).%0A%0A%20%20%20%20-%20%5BCedraro%20A%2C%20Cappello%20A%2C%20Chiari%20L%20(2009)%20A%20portable%20system%20for%20in-situ%20re-calibration%20of%20force%20platforms%3A%20experimental%20validation.%20Gait%20%26%20Posture%2C%2029%2C%20449%E2%80%93453%5D(https%3A%2F%2Fpubmed.ncbi.nlm.nih.gov%2F19111467%2F).%0A%0A%20%20%20%20-%20%5BCappello%20A%2C%20Bagala%20F%2C%20Cedraro%20A%2C%20Chiari%20L%20(2011)%20Non-linear%20re-calibration%20of%20force%20platforms.%20Gait%20%26%20Posture%2C%2033%2C%20724%E2%80%93726%5D(https%3A%2F%2Fpubmed.ncbi.nlm.nih.gov%2F21392999%2F).%0A%0A%20%20%20%20-%20%5BSciPy%20%60pinv%60%5D(https%3A%2F%2Fdocs.scipy.org%2Fdoc%2Fscipy%2Freference%2Fgenerated%2Fscipy.linalg.pinv.html)%20%E2%80%94%20the%20current%20SVD-based%20pseudoinverse%2C%20including%20rank%20tolerance.%0A%0A%20%20%20%20-%20%5BModern%20Robotics%3A%20wrenches%5D(https%3A%2F%2Fmodernrobotics.northwestern.edu%2Fnu-gm-book-resource%2F3-4-wrenches%2F)%20%E2%80%94%20force%2Fmoment%20representation%20and%20change%20of%20reference%20frame.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0Aif%20__name__%20%3D%3D%20%22__main__%22%3A%0A%20%20%20%20app.run()%0A