450 likes | 661 Views
Haptic Augmented Reality Implementation Details and Other Issues. Haptics and Virtual Reality Laboratory POSTECH Seokhee Jeon yeager@postech.ac.kr. Covering. Implementation details of the stiffness modulation system Tumor palpation simulator Simple FEM simulation Silicone casting
E N D
Haptic Augmented RealityImplementation Details and Other Issues Haptics and Virtual Reality Laboratory POSTECH SeokheeJeon yeager@postech.ac.kr
Covering • Implementation details of the stiffness modulation system • Tumor palpation simulator • Simple FEM simulation • Silicone casting • From the proposal to the defense
Haptic AR System Description • Library dependencies • Core haptic library: OpenHaptics 3.0 • Force sensor and DAQ: NIDAQmx, ATIDAQ • GUI and graphics: openGL, GLUI • Visual AR, camera capture, and image processing: ARToolkit, PGRFlyCapture, DevIL • Hardware • NANO17 with NIDAQ (FT8178; 35 N maximum) • PHANToM 1.0, 1.5, 1.5 high force • Custom interaction tool
Interaction Tool • The Solidworks files in vrs4 (or Gabjong)
Haptic AR Code Overview Shared memory Position, Velocity Reaction force Reaction force Reaction force, Velocity, position Device controller Sensor manager GUI manager 1 kHz 10 kHz 60 Hz Butterworth LP filter FOAW RLS modules Virtual model Spline curve Units used in the code Length: mm, Force: N, Weight: g
Starting int main(intargc, char* argv[]) { //Graphics manager instance GUIManager GUI; //sensor manager instance SensorManageriSensor; //Device controller instance DeviceControlleriDC; //Creating a window and GUI GUI.initUI(argc, argv,&iSensor,&iDC); //Visual AR init GUI.ARInit(); //force sensor init. iSensor.initATIDAQ(); iSensor.initDAQ(&iDC); //device controller init. iDC.init(); //60 Hz graphics main loop starts. Then, force sensor calibration and 1 kHz haptic loop starts. glutMainLoop(); //System termination GUI.Quit(); hdStopScheduler(); …… return 0; } SensorManager::sensorCalibration() SensorManager::forceSensorCallback() deviceController::phantomCalibration() deviceController::arLoopCallback()
DeviceController //Main callback staticHDCallbackCode HDCALLBACK arLoopCallback(void *data); //self-explainable by names staticvoidcontactDetection(); staticvoidnormalEstimation(); staticvoiddisplacementEstimation(); staticvoidforceControl(); staticvoidcontactDynamicsIdentification(); staticvoidfrictionIdentification(); staticvoidcontactReleaseDetection(); staticvoidtumorForceCalculation(); staticvoidfrictionModulation();
DeviceController::arLoopCallback() • //Main callback function. called in every 1 ms. • HDCallbackCode HDCALLBACK DeviceController::arLoopCallback(void *data) • //force sensor data acquisition • … • //PHANToM ready • … • //encoder data acquisition • … • //structural stiffness compensation • … • //Reaction force calculation • //Force sensor transform, gravity, inertia compensation • … • //FOAW velocity calculation • … • //Contact detection. This function sets 'contactFlag' variable. • contactDetection(); • if(contactFlag == 1) { //if contact occurs • //main AR algorithms start// • //1. Surface normal estimation • normalEstimation();
DeviceController::arLoopCallback() • //2. displacement estimation • displacementEstimation(); • //3. force control • forceControl(); • if(isContactIdent == TRUE) //Entered if 'dynamics ident.' button is pressed. //contact dynamics ident. contactDynamicsIdentification(); • if(isFrictionIdent == TRUE) //Entered if 'friction ident.' button is pressed //friction ident. frictionIdentification(); //4. contact release detection contactReleaseDetection(); if(tumorEnabled) //Breast cancer palpation// tumorForceCalculation(); if(frictionEnabled) //Simple friction modulation example. Hapticvirtuality frictionModulation(); //Force command to the PHANToM hdSetDoublev(HD_CURRENT_FORCE, force); //Data storing for visualization … }
Estimation of Deformation Direction Original surface Algorithm Off-line Identify the frictionmodel for Rendering Deformed surface Tool tip-object Friction
Dahl Friction Model Identification • Dahl friction model Presliding regime Sliding regime Viscous friction coeff. Stiffness coefficient μk frn(t) Coulomb friction coeff. Slope = σ -μk frn(t)
Dahl Friction Model Identification • Dahl friction model • Off-line parameters identification • Acquisition of real geometry info. using dense sampling of object surface in 2D scanning path • Data collection using manual stroking in the 2D path • Linearize the problem by dividing the data into two bins; presliding and sliding bins • Identification of parameters using linear regression technique Takes 20~30 seconds Viscous friction coeff. Stiffness coefficient Coulomb friction coeff.
DeviceController::frictionIdentification voidDeviceController::frictionIdentification() { • hduVector3Dd tanDisplacement = slipPosition-position; trueNormal = -normalSpline.getSlope(position[0]); //get true surface geometry from the spline curve tanVelocity = velocity - (velocity.dotProduct(trueNormal) * trueNormal); tanVelMag = tanVelocity.magnitude(); tanUnit = tanVelocity; forceNormalComponent = sensedForceT.dotProduct(trueNormal); tanUnit.normalize(); forceTanComponent = sensedForceT.dotProduct(tanUnit); if(abs(tanDisplacement) < (SLIP_THRESHOLD)) {//data bin for pre-sliding regime estimatorS.update(forceTanComponent,tanVelMag); slipPosition = position; } • if(abs(tanDisplacement) > SLIP_THRESHOLD2) {//data bin for sliding regimeestimatorF.update(forceTanComponent,forceNormalComponent,tanVelMag); } //verification. comparison between true and estimate doublefc = forceNormalComponent*estimatorF.uk; doubletanFCurrent = (fc) + ((tanFPrev-(fc) )*exp(0.9/fc*tanDisplacement.magnitude() ) ); tanFWithViscous = tanFCurrent + (estimatorF.sigma*tanVelMag);
DeviceController::normalEstimation voidDeviceController::normalEstimation() • //compensate CONTACT_THRESHOLD from the original reaction force • forceNormalComponent = sensedForceT.dotProduct(surfaceNormal) - (CONTACT_THRESHOLD); • doublefc = forceNormalComponent*estimatorF.uk; //kinetic friction hduVector3Dd tanDisplacement = slipPosition-position; • //kinetic friction + Dahl friction tanFCurrent = (fc) + ((tanFPrev-(fc)) * exp(estimatorSlop.sigma2/fc*tanDisplacement.magnitude())); Ff = (tanFCurrent + (estimatorF.sigma*tanVelMag)) * Vtan;// + viscous friction tanFPrev = tanFCurrent; currNormal = sensedForceT - Ff; //surface normal currNormal.normalize();
Estimation of Deformation Magnitude Off-line Identify dynamicsrelation, G Rendering Take the inverse of G to derive x Force Force Displacement Displacement Hunt-Crossley model Stiffness coef. Damping coef. Exponential coef. G Parameter identification Adapted the Haddadi et al.’s method Takes 10~20 seconds for each object
Identification of Hunt-Crossley Model • Off-line parameters identification • Data collection using manual palpation at a fixed point • Linearize the model by taking log function (adapted from Haddadi et al.’s method) • Identification of parameters using linear regression technique
Contact Dynamics Identification • voidDeviceController::contactDynamicsIdentification() hduVector3Dd iDifference = initialContactPoint-position; trueDisplacement = iDifference.dotProduct(surfaceNormal); //linear recurrsive least square estimator parameterVariation= indirectEstimator.update (forceNormalComponent,trueDisplacement,Vnormal)); if((parameterVariation*parameterVariation) < 1.0f) { //convergence test convergeFlag = 1; } doubletempD = forceNormalComponent/ (indirectEstimator.K + (indirectEstimator.B*Vnormal)); displacement = pow(tempD,1.0f/indirectEstimator.n);
DeviceController::displacementEstimation voidDeviceController::displacementEstimation() doubletempD = forceNormalComponent/(indirectEstimator.K + (indirectEstimator.B*filteredVnormal)); displacement = pow(tempD,1.0f/indirectEstimator.n); contactPosition = position+(displacement*surfaceNormal);
GUIManager • User interface • Data visualization • Visual augmentation
GUIManager: Initialization • initUI • GLUI initialization • Object library load • ARinit • setupMarker() • setupUSBCamera() • setupPGCamera() //Load real object library that stores pre-identified dynamics parameters // //Library file format// //<number_of_objects> //<object_name> <K> <B> <n> <uk> <sigma1> <sigma2> 3 FoamBall 0.1105 -0.0011 1.6306 -0.6474 0.0022 0.1173 SpongeBlock 0.0115 -0.0002 2.0779 -0.6363 0.0008 0.0556 SiliconeBreast 0.07 -0.0004 1.2 -1.6145 0.0073 0.0328
GUIManager: Callbacks • //image capturing and marker tracking for visual AR voidGUIManager::myGlutIdle( void ) • Grab a video frame • Detect AR markers in the video frame • Extract the transform matrix of the markers • Normal haptic AR rendering • data visualization • Visuo-haptic AR • Viewing transformation • Virtual object overlay voidGUIManager::myGlutDisplay( void )
GUIManager: Callbacks voidGUIManager::myGlutKeyboard(unsignedchar Key, int x, int y) Keys: System q Quit demo. ? Help n Change real object dynamics parameters Stiffness e Toggle stiffness modulation Modulation [ or ] Increase or decrease the desired stiffness Tumor t Toggle tumor rendering Palpation 1-4 Select tumor rendering method y Place virtual tumor at 15 mm below HIP Friction s Toggle friction modulation Modulation 5-8 Select friction modulation method Visuo-Haptici Toggle AR screen a or z Scale V coord. with regard to H coord. l.,/k; Translate V coord. v or b V-H calibration or reset calibration p Lock tracking u or j Adjust V-H latency d Toggle drawing ! Toggle tool overlay
Visuo-Haptic Calibration:Coordinate Frames PHANToM Coord. Camera Coord. Least square problem K. S. Arun, T. S. Huang, and S. D. Blostein. Least-squares fitting of two 3-d point sets. IEEE Trans. Pattern Anal. Mach. Intell., (5):698– 700, 1987.
Visuo-Haptic Calibration:Synchronization • xc usually lags behind xp
Sensor Manager • initATIDAQ • initDAQ • sensorCalibration • Called at the beginning of the program • Can be manually called by pressing ‘sensor calib.’ button • sensorCallback • Data acquisition • Update shared memory at every 10th tick
Miscellaneous Modules • Butterworth low-pass filter • FOAW filter • Recursive least square estimator • class ButterworthLPFilter //possible curofffrequenct: 250,200,150,100,50,30,20,10 ButterworthLPFilter(intcutoff_frequency); //true data in, filtered data out double update(double); class FOAW //current velocity in, filtered velocity out doublevel(doublecurrPosition); //window doubleyHistory[MAX_WINDOW]; class estimator //internal variables for estimation. see reference. double P[3][3]; double theta[3]; double phi[3]; double L[3]; double sigma; double forget; //Forget factor for recursive process double K, B, n; //stiffness, damping, //exponent coefficient //passing data for estimation. input: normal force, true displacement, normal velocity double update(double force, double displacement, double velocity);
Miscellaneous Modules • Cubic spline curve • Virtual tumor model classCSpline //Spline curve coefficients … //generate spline parameters. Call this after the all the surface points are sampled. voidgenerateSpline( ); //return the slope of the curve at point x longdoublegetSlope( longdouble x ); //return the y position of the curve at point x longdoublegetPosition( longdouble x ); //collecting the surface points voidgetInput(longdoublexPosition, longdoublefxValue); classVirtualModel doubleradious; hduVector3Dd originalPos, transFromOriginal; double stiffness, surroundingStiffness; VirtualModel(doubleir, double is, doubleig, hduVector3Dd it); voidsetFixedPos(hduVector3Dd ip); doubleisCollided(hduVector3Dd hip); voidgetPosition(hduVector3Dd& targetPosition); doublegetRadious(); voidsetTrans(hduVector3Dd it);
Important Constants • DeviceController.h • CONTACT_THRESHOLD: Decision threshold for contact detection • TIP_LENGTH: The length from the gimbal origin to the tool tip • SLIDE_THRESHOLD1,2: Tangential displacement that separates the sliding and presliding regimes in friction • FOAW.h • MAX_WINDOW: Maximum winndow size • T: sampling period • D: position resolution of input data
Structural Stiffness Estimation • Actual displacement • Estimation • Tapping on a virtual wall with stiffness • x(t), xs(t) acquisition • Model fitting From the encoders From LVDT displacement sensor
Stiffness Estimation of PHANToM 1.5 HF : 2.537 N/mm
Ks Compensation Results =1.5 N/mm =1 N/mm =0.2 N/mm
Tumor Palpation Simulator • Force at user’s hand From real breast model From virtual tumor
Preprocessing Tumor Response Real tumor mock-up Force Hunt-Crossley model identification Displacement : No-tumor breast model response Force : Response of tumor alone Displacement Hunt-Crossley model identification again
Rendering • See “void DeviceController::tumorForceCalculation()” Rendering Approximate ft(t) at positions not lying on l0 using the identified response Virtual tumor Tool tip
Simple Linear Elastic FEM Simulation • Tetrahedron elements construction • NGSolve • Constructive solid geometry tetrahedron list • Take a look at an example • ‘Export mesh’ natural mesh format Operations in constructive solid geometry
Simple Linear Elastic FEM Simulation • FEM simulation • Let’s see an example
Silicone Casting • Materials • ECOPLEX 0030 (Smooth-On, Inc.) • Part A and Part B • SILC FLEX: Silicone thinner • Reseller: www.hyup-shin.co.kr (협신물산) • Mixing formula • More thinner softer but more sticky • Always combine equal amount of parts A and B • Preparation • Detergent, Original form, stick, paper cup, rubber gloves
Silicone Casting • Making a mold • Fix the original form in a paper cup (using a string) • Mix part A : part B = 1 : 1 (hardest) without the thinner • Pour the mix into the paper cup • Wait for at least 4 hours to be completely cured • Remove the original form • Casting a copy • See formula and mix A, B, and thinner. • Apply a detergent over the surface of the mold mold separator • Pour the mix into the mold • Wait for 12 hours (more thinner more time to be cured) • Carefully separate the copy out of the mold
Silicone Casting • Note • Mixing part A and B: stir for at least 10 minutes • More thinner, more sticky surface (thinner does not completely dry up) • Making original form: clay, Sculpy • Inhomogeneous silicone • Attach string when making a part inside • Hang the part to be placed inside the surrounding part
From Proposal to Defense • Proposal • D-45: Select the committee members and contact to them with the brief explanation of your topic • D-30: Arrange the proposal schedule with the members, don’t forget to ask if the external member needs accommodation, Hotel reservation • D-15: Finish a proposal paper and send itto the members, 학위논문 외부심사위원 위촉승인서, 출장신청서 제출 • d-5, d-1 phone call or e mail for the reminder • Practice and practice (attack defense!) • D day: 논문연구계획서 싸인, 학과 제출, Dinner with professors!, comment 정리 • D+2: 놀기끝!
From Proposal to Defense • Defense • Need at least one enlisted journal paper (acceptance notification is OK) • D-45-30: Arrange schedule with the members and remind them about your topic (summary of your research) • D-15: Send your thesis (hardcopy) to the members, 논문심사요청서 과사 제출, 학위논문 외부심사위원 위촉승인서, 출장신청서 제출 • D-5, D-1: Reminder, • D day: 논문심사결과보고서 (싸인), 식사 준비, • 12-31, 6-31: 논문심사결과보고서 과사제출 • D+7: Thesis revision, Thesis seal • 1-6, 7-6: Bind the thesis and 1 과사, 4 library, 5 professors • 졸업정산
Finding a Job • Make connections with people! • Start searching for a job as early as possible! • CV and resume management