Controlled-STM: A Two-stage Model to Predict User’s Perceived Intensity for Multi-point Spatiotemporal Modulation in Ultrasonic Mid-air Haptics

Multi-point STM offers a great range of parameters (i.e., drawing frequency, number of points) to produce different tactile sensations. However, existing studies offer limited insight on the effects of these parameters, and ignore their effect on the physical stimuli delivered, limiting effective haptic design. We propose a two-stage model to predict response to multi-point STM. The first stage predicts physical stimulus properties with 7.8% error, while the second stage predicts mean and spread of perceived intensity with 8.0 % and 8.8% error. We report 3 studies conducted to derive this model: one to characterize physical stimuli, another one measuring user perceptual thresholds, and a third one measuring user’s perceptual response to multi-point STM. Besides, we characterize 4 effects that influence device performance, confirm if previous effects reported are due to physical or perceptual effects (or both) and derive recommendations for manufacturers, haptic designers and HCI researchers.


INTRODUCTION
Multi-point Spatiotemporal Modulation (STM) uses one or more points of focused ultrasound (i.e., focal points) quickly moving along the users' skin to create tactile sensations [11,55].Several parameters have been found to afect how users perceive STM stimuli, such as drawing frequency ( ) [11,38], the way the shape is discretised [12] or the number of points ( ) used [55].
However, all prior studies directly map the efect of such highlevel STM parameters (e.g., drawing frequency) to users' perceived intensity ( ), without paying any regard to the efect these highlevel parameters could have on the physical properties of the stimuli being delivered (e.g., pressure/force of each point).
As such, it remains unclear if the diferences found are indeed the result of the users' perceptual response (e.g., tactile receptors) or due to the physical properties of the stimuli (e.g, stronger/weaker focal points).Similarly, the physical properties that better determine users' response remain unclear (e.g., pressure or force?; in absolute terms or relative to users' perceivable thresholds?).These uncertainties impede designers from predicting the efect of the tactile stimuli they craft, and limit HCI researchers' insights into how humans respond to multi-point STM.
We address these limitations by constructing a two-stage model: The frst stage characterizes the physical properties of the focal points (e.g., peak pressure, forces) from the input parameters (i.e., number of points used, drawing speed/frequency, intensity and temperature).The second stage builds on these properties to predict the perceptual response on the user.
This approach allows us to decouple physical and perceptual responses, gaining understanding on the parameters that afect the physical stimuli (i.e., due to device's artefacts/limitations), but also the parameters and physical properties that drive users' response to multi-point STM (i.e., human perception).While the frst part might vary from one device to another, the second one only depends on the users' response to a given physical stimuli and could be reusable across devices/manufacturers.
We conducted 3 experiments in order to construct our model.The frst study measured the physical stimuli delivered by an Open-MPD device [41], as input parameters changed (e.g., command intensity ( ), number of points ( ), drawing frequency ( ) and temperature ( )), and explored their relationships to produce our physical model (i.e., frst stage).To model the physical stimuli relative to users' perceivable thresholds, we conducted a second study measuring participants' minimum perceivable thresholds for STM at several drawing frequencies.Finally, our third study consisted of a magnitude estimation task, with varying number of points, intensities and drawing frequencies.
We derive our model from these studies, built as a simple combination of analytical expressions to: 1) predict the acoustic pressure under diferent , , and conditions with an average relative error of 7.8%; 2) predict users perceived intensity based on , and the physical properties predicted from the physical model, with an average relative error of 8.0%.Our model and study results allow us to derive recommendations for several audiences, such as device manufacturers (i.e., efects hindering devices' response, and their relative relevance), haptic designers (i.e., use our model as a potential tool to predict users' response and inter-subject variability, while they are designing the stimuli) and HCI researchers (i.e., insights on human response to multi-point STM).

RELATED WORK
Ultrasonic mid-air haptics (UMH) renders haptic sensation through focused ultrasound, delivered via Phased Arrays of Transducers (PATs).Each transducer is driven with a phase delay such that the waves of all transducers converge on the same point and at the same time, creating a point of high pressure (i.e., a focal point [6]).
Modulation techniques are required to make such points perceivable [55].Amplitude Modulation (AM) changes the intensity (i.e., pressure) of each focal point over time [21], while Lateral Modulation (LM) retains focal points active at all times, using small lateral oscillations to make them perceptible [47].Spatiotemporal Modulation (STM) leverages focal points moving quickly along a tactile shape, either using one [11,28] or many [39,55] points.
Our review focuses on STM techniques.First we describe studies exploring the perceptual response to STM stimulation.Second, we focus on the efects STM parameters can have on the physical stimuli.Third, we summarize prior studies describing how physical properties of the stimuli afect perceived intensity ( ).

Perceptual response to STM haptics
The acoustic pressure (Pa) of each focal point is arguably the most defning parameter determining the intensity of the stimuli, and algorithms exist controlling the intensity of each point [31,39].
However, the way those points interact with users' skin is similarly relevant.Frier et al. [11] found a quadratic relationship between and the speed of the STM point (i.e., drawing speed or frequency , both terms are interchangeable [39]).Similarly, the way the tactile shape is sampled (e.g., points along the shape, time spent at salient points) have also shown to afect [12,34].
Finally, the number of points ( ) used to render the shape have also been shown to afect [55], along with other emotional responses.Other studies have further looked into the efects of STM on emotional responses (e.g., [3,37]), but this paper will limit its focus to predicting .
In summary, prior literature has shown how several parameters infuence , namely: Command Intensity ( , target pressure for each focal point (in Pa) requested to the algorithm); Drawing frequency ( , number of times each focal point traverses the shape per second); Number of Points used ( ) and sampling (location and timing of focal points along the shape).
However, none of these studies looked into the efects that changing such parameters had on the physical stimuli delivered (e.g., pressure of the focal points, forces).They all attribute the efects observed strictly to users' perceptual response, without considering if the physical stimuli used for their comparisons are, indeed, comparable.As such, it is not currently possible to assess if such efects are due to users' perceptual response, or to the way the physical stimuli is delivered by the PAT device.This direct mapping also triggers generalization problems.For instance, Ablart et al. [2] concluded that drawing speeds of 12.5 m/s introduced higher perceived intensity than 10 m/s, while this result was the opposite in [11].Thus, it further encouraged us to study the underlying physical effect of parameters on perceived intensity to ensure reliable delivery of the haptic experiences.

STM parameters and the physical stimuli
There are a number of efects that can afect the physical stimuli (i.e., the focal points) delivered by PATs , such as non-linear efects, limited power, sub-optimum driving frequencies of the transducers or transducers' temperature.We here discuss such efects and how the STM parameters used can infuence them.Non-linear acoustics is a branch of physics which comes into play when dealing with high amplitude waves [32].Linear models do not apply in these cases and such systems are governed by fuid dynamics instead.
Existing STM algorithms [31,39], however, rely on linear acoustics and analytical models of transducers' directivity (e.g., piston model [41]), to determine the acoustic pressure delivered.While the directivity and intensity of each individual transducer has been experimentally validated by several works [19,35], the joint contribution of all transducers in the PAT could (and does, as show later in the paper) lead to sufciently high pressures for non-linear efects to appear.As such, STM solvers could attempt to produce intensities ( ) beyond the linear regime, and fail to predict the intensity of the focal points they are actually creating.
Secondly, the PAT is limited in terms of the acoustic power that it can deliver, but the underlying algorithm can still attempt to create focal point intensities beyond such limited power (i.e., the actual focal point would be weaker than the intended requested to the solver).This is particularly complex if several points ( ) are being created.In this case, the designer would need to retain a notion of feasible ′ , per number of points used, or resort to tools to assess that the algorithm can indeed achieve such intensities.
Thirdly, the high update rate of STM algorithms can lead to suboptimum driving frequencies of the transducers.Transducers are narrow-band electro-mechanical actuators, designed (in most PATs) to operate at 40 KHz.At the same time, STM algorithms require update rates [11,39] in the same order of magnitude than this frequency (e.g., 10 K updates/s vs 40 KHz transducers' frequency).Each update forces a change in phase every few cycles (i.e., 4 cycles, in our example), forcing transducers to operate outside of their resonant frequency and leading to loss of acoustic power [50].
The larger the phase change required per update, the higher the power loss will be, and several STM parameters can require such high phase changes.Higher speeds ( or ) will require more aggressive phase changes, as the phase change depends on each points' displacement relative to the transducer [19].The use of several points ( ) will also require phase-retrieval methods [39], involving less predictable phase updates.
Finally, these same reasons can afect transducers' temperature.The impedance of transducers can be matched for specifc frequencies [13], leading to minimum power consumption and heat dissipation.Operating outside the nominal frequency, however, would still result in higher consumption and heat.Such temperature changes directly afect transducers' output and resonant frequency [50], reducing the performance and accuracy of the PAT.

Physical stimuli and perceived intensity ( )
There are prior works exploring perceptual responses to vibrotactile stimuli [18,22,52].However, most of them focus on contact-based actuators (e.g., vibrators), and it remains unclear if such fndings are transferable to multi-point STM.
As indicated earlier, most studies focus on the efects that STM parameters have on , not considering the physical properties of the stimuli.As such, the physical properties of the STM stimuli (see Figure 1.centre for a typical pressure distribution of a focal point) that better determine are not widely understood.
Some studies use the peak pressure of the points (i.e., summation) as an estimate of the intensity of the tactile stimuli (i.e., see 0 , in Figure 2.a).In GS-PAT [39], authors used this metric to argue that a small number of points ( ) should provide similar STM quality than a single point, but with no actual validation.A later study on multi-point STM [55] disregarded such hypothesis (i.e., was higher for = 1), and suggested that considering the peak pressure relative to users' minimum perceivable threshold ( ) [18] could be a better predictor (see , in Figure 2.b).
Although not specifc from STM, other approaches [10] focus on the acoustic force (i.e., pressure integrated over unit area) delivered by the stimuli.More specifcally, they used a force scale with a sensing plate of 2.1 in diameter, comparing such direct measurements with their force estimations across such area (see 2.1 , in Figure 2.c), but without validating its actual correlation with .Alternatively, knowledge of the users' for STM stimulation could allow for more refned estimations, integrating force only along the area perceptible to the user, either in absolute or in relative terms (i.e., 0 and respectively, in Figures 2.d&e).
Some works also look into the efects of ultrasound stimuli on the users' skin [10], considering their rheological properties to model other efects, such as indentation [11].We however limit our exploration to characterizing the stimuli, and the properties that better predict (i.e., 0 , , 2.1 , 0 and , in Figure 2; Refer to Supplementary Material to see how all these properties can be derived from the pressure distribution and 0 ).

CONTROLLED-STM: SUMMARY OF APPROACH
Our core contribution is the description of a predictive model for Perceived Intensity ( ) for Multi-Point STM stimulation, addressing all modulation parameters (i.e., , and , approach summarised in Figure 1).As a main distinctive point, we do this by looking at the efects of such parameters both on the physical properties of the STM stimuli and on the users' perceptual response, providing insights on which efects are the result of the way the device delivers the stimulation and which ones are inherent to users' perceptual response to STM stimulation.This approach results is a two-stage model: a Physical Model characterizing the way STM parameters infuence the stimuli, and a Perceptual Model predicting users' from the stimuli's parameters.More specifcally, the Physical Model is composed of two components.The frst one 0 (, , , ) (extracted from Study 1) predicts the peak acoustic pressure of each point (i.e., 0 , in Figure 2) from our STM parameters and the device's temperature ( ), as a combination of several functions: non-linear efects ( ), driving frequencies ( ), transducers' temperature ( ) and limited power (), with each function representing a corresponding physical effect afecting device's performance.Here, (or =6531,, =31 ) represents the pressure measured for a single static focal point at CI=6531 and T=31 • C (i.e., default case): (1) The second component of the Physical Model is ( ), describing users' minimum perceivable threshold to STM stimuli, as a function of (extracted from Study 2).This function is retained as part of the Physical Model as such thresholds are required to extract the properties described in Figure 2 (i.e., particularly and are relative to ), but ( ) indeed plays a pivotal role between the physical and perceptual stages.
The Perceptual Model (i.e., second stage) is composed of two functions: i) , predicting average Perceived Intensity ( ); and, ii) predicting the variance (standard error) of across users.They both depend on , and , as extracted from Study 3: The following sections detail each stage, describing the 3 studies conducted to construct our model, as well as the analysis of results, efects observed and rationale that led to this design.

PHYSICAL MODEL
The physical model aims to characterise the properties of the physical stimuli delivered, when diferent STM parameters are used (i.e., , , , ).The results from Study 1 provide an estimation of    We used the open-source platform OpenMPD [41] for all our tests, both for software (GS-PAT algorithm, at 10 KHz) and hardware (PAT with 256 transducers working at 40 KHz, operated at 18 V), modifed by adding a DS18B20 temperature sensor to the back of the PAT to collect temperature readings.We used this setup to create diferent multi-point STM stimuli, measuring the resulting pressure as we varied all 4 parameters afecting device performance ( , , and ).
More specifcally, we measured intensities ( 's) in a range from 544 Pa to 6531 Pa, with the latter being the highest achievable that the solver can deliver to a single focal point.We explored drawing frequencies = {5, 10, 20, 40, 80, 120, 160, 240, 320 Hz}, as to cover the sensitive range of both Meisnner and Pacinian mechanoreceptors [18].We limit our exploration to = {1,2,3} points, as both physical measurements [39] and prior studies [55] indicate signifcant decreases in performance past this number of points.Finally, we characterised the PAT's performance as temperature increased from 25 • C (room temperature) to 50 • C (safe upper limit).We fxed the remaining STM parameters to avoid confounding factors.The height of the focal point was fxed at 15 cm from the device as to avoid changes in pressure related to location [10,39,40].We use a circular pattern (i.e., as to avoid salient points/corners' infuence [34]) of 6.5 cm (within the typical size of an adult's palm of 7.5-9.5 cm [30]), discretised using 100 positions homogeneously sampled along the circle shape (i.e., avoid timing efects [34]).It must be noted that for higher > 100 Hz, it is impossible to render the shape with 100 positions (i.e., at 320 Hz, one can only use 10 ≈ 31 samples).However, Frier et al. [12] found that is not signifcantly afected by the number of samples for such high frequencies, so we did not consider the diference in sampling rate as a confounding factor in these cases.

Measurement setup and procedure.
The measurement device setup is shown in Figure 3.The PAT was placed above foam squares to absorb environment vibrations during our tests.A motorized stage (i.e., CNC machine) was used to place a Bruel & Kjaer 1/8" Type 4138 Pressure-feld Microphone at one point of our STM circle (i.e., (0mm, 31mm, 150mm)).A small search was performed over an area of 8x8 mm to ensure the microphone was placed at the point of maximum pressure (i.e., right on top of our haptic shape).The microphone was connected to a PicoScope sampling at 10 MHz, and each test recorded enough samples as to capture a full revolution of all focal points around the circle (i.e., 10 / samples), and each measurement was repeated 3 times.
As per the microphone data-sheet [29], a free-feld measurement is approximated from the pressure response using the included graph, and 0.78 dB is subtracted from the measurements, as the sound is measured at a 90 degree angle to the microphone tip.This accounts for the increase in pressure due to difraction of the sound waves around the microphone.
We report 4 experiments, each characterizing one of the 4 physical efects afecting the device's performance described in Related Work.For each efect, we only varied (and measured) the parameters infuencing said efect, retaining other (independent) parameters fxed.More specifcally, we consider the case of a single, fxed focal point as a default case (i.e., =6.5 KPa, =1, =0 Hz, =31 • C), and use these parameters as the default fxed values for the independent parameters, as we explore the infuence of varied parameter on the physical stimuli studied.However, for the LP efect, it is constructed diferently than other efects as the efect of on physical output is also dependent on the CI.In any case, this allows us to approach our physical model as a simple combination (multiplication) of correction factors, each derived from a specifc efect when compared to a default case (i.e., fxed, single focal point).Dependent and independent parameters were identifed during preliminary tests (not reported here), but our evaluation covers arbitrary combinations of all parameters, validating these assumptions and also the general correctness of our physical model.
Our results, shown in Figure 4.a, clearly illustrate the presence of non-linear behaviours.The response is reasonably linear for values < 2 KPa (green dashed line), showing good agreement between the expected pressure ( ) and the measured pressure.However, for higher values the PAT fails to deliver the pressure expected (i.e., actual pressure is lower than demanded ).
This relationship can be ftted using a sigmoid function (orange continuous line in Figure 4.a), and the parameters fully defning function ( ) can be found in Supplementary Material.
Our results are shown in Figure 4.b.As expected, the pressure delivered decreases signifcantly at higher frequencies (i.e., >80 Hz).It is however surprising to see that the pressure delivered can exceed that of a static focal point at some frequencies (i.e., .This could be caused by a mismatch in the real-time clock of the PAT (e.g., its FPGA does not deliver exactly 40 KHz), due to transducer's resonance (i.e., actual frequency slightly higher than 40 KHz), or due to frequency shifts due to operating temperature.We repeated these tests 3 times and, whatever the cause, we can confrm this efect is consistently present every time.
This relationship between and pressure can be ft with a 3rd degree polynomial equation ( ), with its specifc parameters being provided in Supplementary Material.
Our results in Figure 4.c show a clear linear relationship between and 0 .Delivered pressure decreases ∼83 Pa for every degree increase (full parameters of ( ) in Supplementary Material).
-Power Limit efects: This efect comes into play when a device cannot deliver the intensities requested, either because the 's are too high, or because it cannot deliver such for such a high number of points.To characterise this efect we measured combinations of 5 diferent s (range 544 Pa to 6.5 KPa) and 3 (1,2 or 3 points), for a total of 5 x 3 x 3 repetition = 45 measurements.The other parameters remained at =0 Hz and =31 • C.
Our results show that, while the PAT can achieve intended pressures for any number of points at lower (i.e., 544 and 1088 Pa), it cannot deliver =2721 Pa for =3; or any higher s for of 2 or 3 points.Fitting requires a more complex multi-variate cubic function (, ), available in Supplementary Material.

Evaluation of the Physical Model.
We evaluated the ability of our physical model to predict the pressure the device would deliver under various combinations of , , and .More specifcally, we tested our model's predictions against the data gathered for the tests in Section 4.1.3(i.e.training dataset), but also against a validation dataset, exploring a broader range of combinations of parameters.
The latter validation dataset explored factorial combinations of: 5 diferent values of s (544 to 6531 Pa); including 1,2 or 3 points; and including 11 frequencies (i.e., 9 prior frequencies between 5-320 Hz, plus 2 unseen frequencies 60 Hz and 200 Hz).During these tests, device temperature varied between 24.2 and 51.5 degrees, remaining very close to our intended operational range (i.e., 25-50 • C).Thus, the validation dataset measured a total of 5 x 11 x 3 = 165 diferent STM stimuli, each at a variable temperature.
We test the accuracy of our model by looking at the relative error between our prediction and the actual measurement (i.e., absolute diference between predicted and measured, divided by the pressure measured).The relative error was 2.05% for our training dataset and 7.8% (close to 6.49% obtained by [40]) for our validation dataset (i.e., detailed relative errors at each condition are available in Table 1 in Supplementary Material).We repeated this validation on a second OpenMPD device (i.e., see section C.2 in Supplementary Material), obtaining a relative error of 8.3% and confrming the trends observed for each physical efect (e.g., temperature being linear).This illustrates that, even if tuning might be required for each device, the general model could generalize to diferent devices.This error from our model is in the same order of magnitude that the inherent error of the solver used (GS-PAT).That is, GS-PAT produces a mismatch of ∼3% between the and the pressure delivered, even if an ideal PAT (free from any physical artefacts) could be used.Also, our error is smaller than the 12.2% JND pressure [42], suggesting that such diferences would not be perceivable by users.This ofers a stark contrast if no model was to be used (i.e., assume GS-PAT simulated pressure resemble real pressure), which would provide a relative error of ∼102%, as detailed in Figure 3 in Supplementary Material.

STUDY 2: Minimum Perceivable Thresholds
Computing some of the parameters describing the physical stimuli (i.e., , and ), required us to determine the .We assume such thresholds will be independent of other parameters (e.g., and ), and explore ( ) only as a function of , which we derive from this user study.

Experimental Design.
We explored s under nine diferent frequencies ( = 5, 10, 20, 40, 80, 120, 160, 240, and 320 Hz).In order to do this, we used a 1-up-1-down staircase design [4,8,40], interleaving individual steps from 2 staircases.The frst staircase started from =1.5 KPa, decreasing with every step, while the second staircase started at the minimum =500 Pa value.Interweaving both staircases helped reduce predictability between current and the previous response, and provided more data (2 repetitions) for analysis.
In either case, participants were asked to answer if they perceived the stimulus (i.e., 'Yes' or 'No' answer), and we used a change of 3 dB in at each step, reducing it to 1 dB after three reversals [18,40].In order to reduce estimation bias, each staircase ran until at least eight reversals were obtained [53].Please note, the initial staircase values (i.e., 1.5 KPa and 500 Pa) were selected from a pilot study, but this design allows testing of threshold values outside this range (e.g., if = 1.5 KPa is not felt, intensity would increase).Each participant completed nine tests, one for each test frequency [18], counterbalanced using a Latin square design.The full experiment took a total of 90 minutes and each participant was compensated with £15 amazon voucher.Ethical approval was obtained from University College London's internal ethics committee and all participants gave informed consent (approval number: UCLIC_2021_014_ObristPE).

Experimental setup and procedure.
We reused the same platform used for Study 1 (OpenMPD software and hardware, with a temperature sensor).The PAT device was enclosed inside a black foam-made box covered with sound absorbing foam, and with a 10 cm by 10 cm square aperture on the top, so that the haptic stimuli could reach the participants' hands.The box also helped participants to align their hands above the device, keeping it at the right distance of 15 cm.An adjustable chair and elbow support were provided, adjusting them for the participant's comfort, to avoid fatigue and to ensure they kept their hand still for the stimulation period.
Noise-canceling headphones playing pink noise were used to prevent participants from being afected by ambient noises and noise from the ultrasonic mid-air haptic device.Participants were required to use their dominant hand to interact with the mouse to record their responses through a GUI displayed on a 22" monitor in front of them.The overall experiment setup can be seen in Figure 5.
A total of 18 participants were recruited (9 females, mean age ±: 26.94±4.64),6 of them with prior experience with UMH.After welcoming each participant into the room, the setup (i.e., box, elbow support and chair) was adjusted to the participant's needs.A prior temperature-taking process was taken on the participant's palm to ensure the temperature is above 35 • C to reduce confounding factors.Each participant was asked to watch a video explaining the experiment procedure, tasks to complete, and instructions on using the GUI and report whether they could perceive the stimulus by answering 'yes' or 'no'.Moreover, participant responses and instructions were provided via the GUI.This ensured that all participants received the same information, minimizing potential biases introduced by the researcher's instructions.
Although related studies do not include any training sessions [9,16,40,43], we included a staircase with =5 Hz terminating after a single reversal as short training [18,49].This was followed by instructions to help participants become familiar with the testing and rating procedure.After training, participants were exposed to 9 tests (one per frequency), each featuring 2 interleaved staircases.After perceiving one stimulus in a staircase for 5 seconds, the GUI instructed the participants to enter whether they perceived the stimulus or not.This short duration of training, compared to the length of the whole study (i.e., 20 seconds vs an hour and half), should ensure minimum interference of training with the testing stage.
Two-minute breaks were scheduled after every 40 stimuli ratings, including a temperature-checking process to verify the participant's palm temperature remained above 35 • C. Each stimulus was followed by a 5-second interval for users to provide their response, with no stimuli being delivered during this period of time.This aimed to reduce sensory bias [40], and to avoid enhancement efects (increased perceived magnitude of a stimulus due to prior stimuli decays to zero after approximately 500 ms [51]).

Results.
After data collection, values were converted to 0 pressures using our physical model and the average of each staircase's last 8 reversals was used to estimate the 50% [53], with results summarized in Figure 6.We removed one outlier (i.e., response outside the 1.5 * interquartile range for that ).Results were normally distributed (Shapiro-Wilk, > .05)but violated sphericity assumptions (Mauchly's test [36] with < .05).Hence, we used Repeated measures ANOVA with Greenhouse-Geisser sphericity correction [15], and interpreted efect sizes using partial 2 [7] and Holm-Bonferroni corrections [1], fnding signifcant efects of on (f-value = 24.282,< 0.001; partial 2 = 0.603), justifying the use of as a function parameter for our function.Average s, their standard deviation for each , and post-hoc results can be found in Tables 3 and 4 in Supplementary Material.
Participants showed maximum sensitivity (lowest ) when is around 20 -40 Hz, which is to be expected as these frequencies trigger Meissner and Pacinian corpuscles [27].Also, such 's would translate to drawing speeds of 4 to 8 m/s, aligning with the fndings in [11] and [55] that such speeds increase user's .Sensitivity decreased (i.e., higher s) for higher frequencies, contradicting fndings from AM stimulation [40] that 120Hz and 40Hz result in minimum and maximum .This result suggests that humans perceive frequencies rendered by STM and AM diferently.
The standard deviation of the increases for 5, 120, 160, 240, and 320 Hz, frequencies to which participants are less sensitive (i.e., high values).This indicates that users found it hard to provide a congruent rating for those frequencies, encouraging the use of 10 to 80 Hz as an optimum range for multi-point STM.
Finally, this experiment allowed us to derive an analytical expression for our ( ) function (shown in Figure 6), as well as to compute our 5 physical properties, introduced in related work and detailed in Supplementary Material.

PERCEPTUAL MODEL
The Perceptual Model aims to provide an estimate of users' response to multi-point STM, including both average and its standard error (i.e., indicator of spread, or inter-subject variability).Study 3 describes a magnitude estimation study run to retrieve this information, describing both our fndings, and the model produced (and the variables that better predict ).

STUDY 3: Magnitude Estimation of
This section describes a magnitude estimation study looking into the reported by participants, as the input parameters of our multi-point STM technique varied (i.e., , , and ).

Experimental Design.
We created STM stimuli varying all input parameters, according to a full factorial design.More specifcally, we reused the same nine frequencies used for our prior studies.For each frequency, we tested 5 diferent values of , evenly distributed between ( ) (i.e., minimum Pa perceptible to users for that ) and maximum = 6.5 KPa.We also tested of 1, 2 or 3 points.Temperature ( ) was monitored during the study, but is not an independent variable.We fxed the remaining STM modulation parameters in the same way than Study 1 (i.e., sampling rate of 100 points, height of 15 cm, circle diameter of 6.5 cm).Each stimuli was repeated 2 times, resulting in a total of 9 x 5 x 3 x 2 = 270 trials per participant.
We used a within-subjects experimental design and samples were presented in a pseudo-randomized order to each participant (i.e., reduce learning efects).Ethical approval was obtained from University College London's internal ethics committee and all participants gave informed consent (approval number: UCLIC_2021_014_ObristPE).

Experimental setup and procedure.
We re-used the same experimental setup used for Study 2, and recruited a total of 18 participants (9 females, mean age ± : 31.72 ± 10.33), 7 of them with prior experience with UMH.
After welcoming each participant into the room, the setup was adjusted to the participant's needs (i.e., box, elbow support and chair) and an introductory video was shown to the user, to explain the experiment procedure, GUI, tasks to complete and rating procedure.After being instructed, a fxed set of 5 training stimuli was used to help participants familiarize themselves with the stimuli, testing and rating procedures and to consolidate their rating scale (particularly important given the unbounded scale used for , which we explain below).
After training, each participant was presented with each of the 270 test stimuli.Users experienced each stimulus (i.e., a given combination of , and ) for 5 seconds, and reported their estimated value for using the GUI. was rated using the absolute magnitude estimation method [23,40], following a user-defned scale from 0 to infnity, which was then normalized between 0 and 1, allowing us to compare results with previous studies using the same metric [11,12,45,55].Regular 1 minute breaks were provided every 40 trials, and the overall duration of the experiment was 40 minutes.At the end of the experiment, each participant received a £10 amazon gift card for their participation.

Results: Multi-point STM efects on
We explored the efects on of our independent variables (i.e., , , ), but also of the physical properties derived from our Physical Model (i.e., 0 , , 0 , 2.1 , ) on the optimum range.For the frst ones, we used Friedman's test (i.e., data was not normally distributed) and interpret efect sizes using Kendall's Value [48].The latter ones cannot be interpreted in this way (i.e., our physical properties are continuous variables), so we instead computed distance correlation coefcients () [46] between each property and .Such efects are summarized in Table 1.
These results indicate that the intensity of the stimuli ( ) bears most responsibility when determining .However, results from our Physical Model also showed how this is an unreliable metric to characterise the physical stimuli.In order to explore if our physical properties could provide a better predictor for , we investigated the correlation between each physical property and perceived intensity (see values in Table 1).This comparison shows that all the properties derived from our physical model provide very high correlation values, all becoming potentially good predictors for .Pressures seem to provide better prediction than forces.This is reasonable as pressure encodes both physical stimuli and area of contact and provides a more accurate representation of user's tactile experience.The inclusion of s also improves (even if only slightly) prediction results.Among all 5 properties, provided the most accurate results and we chose this property to continue to explore our perceptual model.
We built our model around the 4 variables identifed: , (as a potential replacement to ), and , by looking at the efect that introducing each of them has on our prediction accuracy.Table 2 summarizes the diferent models we tested, describing the variables considered (i.e, , , and ), and reports the relative errors that we could achieve both to predict average (i.e., to model how strongly users' would perceive a STM stimuli) and standard error of (i.e., , to model how diferently users' would perceive it) above the MPT.In all cases, we limited our exploration to quadratic models.The way they were derived and the full analytical expression for each model can be found in Supplementary Material.Looking frst at the accuracy of models including a single variable (i.e., , , or ), we can confrm that indeed provide the highest prediction accuracy.
Figure 7 helps visualize the relative efects of such variable, by comparing users' to , both normalized to a 0-1 range.The fgure shows a very strong correlation between and , but also how provides higher errors for specifc frequencies (e.g., high for =5 Hz, or low for = 80 Hz), pointing at perceptual efects related to other variables (e.g., ).
The fgure also illustrates the low values achieved for frequencies outside the 10-80 Hz range, a result anticipated by our previous studies (i.e., weaker physical stimuli for higher frequencies from Study 1, and higher variance of s from Study 2).
The next models in the table iteratively include additional variables (i.e., and ) on top of , as a way to understand the gains from including additional parameters to the model.Including results reduces relative error to 11.9% and 10.0% for average and , improving accuracy by 7.2% and 0.3% respectively, when compared to using alone.
The inclusion of , however, only results in very marginal gains (< 2% improved accuracy).This seems to indicate that such efect has already been included in the formation of .That is, Study 1 showed how afects the intensity of the stimuli, but this is a factor that is already captured by and it seems alone can mostly model such efects successfully.
Our fnal model, including all variables (i.e., , and ) resulted in relative errors of 8.0% for average and 8.8% for , and its overall ftting is shown in Figure 8 (See Tables 7 and 8 in the supplementary material for models' absolute and relative error at each optimum across 1, 2, 3 for and prediction.).Please notice that our results are comparable to other vibrotactile predictive models.Prior works [44,54] indicated that an increase in amplitude tends to elevate the perceived intensity of vibrations according to Steven's power law.While these models require an expression for each frequency they consider, our model, however, allows predictions across a continuous frequency range (10-80 Hz).Although we applied diferent ftting functions, both of our perceived intensity prediction errors consistently fall within an acceptable range, specifcally lesser than the tactile discrimination at 20% [5] mentioned in [54].
Figure 8 further reinforces the strong connection between the and .The horizontal axis shows the range that the STM stimuli in our study generated, showing how frequencies of 20 Hz and 40 Hz could achieve higher pressures.Unsurprisingly, those frequencies with higher range also achieve higher 's, indicating that the range of 's an STM stimuli can generate is mostly related to the pressure it can deliver.
Supplementary Material explores this efect also for frequencies outside our range (e.g., 5 -320 Hz), further showing how highfrequencies (which participants struggled to perceive) indeed provided stimuli with very low .Hence, it will be hard to truly discern users' response to higher frequencies until devices become capable of accurately delivering these frequencies.
It worth noting that our most perceivable frequencies (20 Hz and 40 Hz) result in equivalent drawing speeds of 4 and 8 m/s, aligning with the results from Frier at al [11].However, their observation that surface waves propagating on the skin were stronger could be the result of measuring a stronger stimuli, rather than simply the optimal drawing speeds that produced the phenomenon of constructive interference.

DISCUSSION
In this paper, we looked at the efects that various modulation parameters ( , , ) have on multi-point STM, looking both at their efects on the physical stimuli (pressures and forces delivered) and on users' perceived intensity ( ).Our exploration helped us characterize 4 efects (i.e., non-linear efects, limited power, suboptimum driving frequencies and temperature) afecting the physical stimuli delivered by the PAT device, and allowed us to derive a 2-stage model predicting users' average and standard error (variability from one users' response to another's).However, there are also limitations in our approach.
First of all, our model exploration was focused on identifying which variables provided a better predictive power to ft the results from our study.Further validation studies would be required to assess the generalizability of our model to real-world situations.Also, other variables such as diferent shapes [17,31], locations on the skin (e.g., studies on forearm [38]), and sampling strategies [12] can be considered for generalization purpose.
While our perceptual model is constructed based on the universally producible physical properties of every UMH device, future works should continually validate and refne the model to ensure its robust applicability.While our use of open-source hardware and software make our results reproducible and directly applicable to the OpenMPD community, our model (i.e., the physical model and Finally, our model uses high-level parameters ( and ) as proxies to assess the efect of phase changes (i.e., frequent phase changes are what cause sub-optimum driving frequencies).While it is useful for our purposes to use and (i.e., they are inherent parameters for multi-point STM), works aiming to characterise PAT response in general should explore the efects of sub-optimum driving frequencies directly in terms of transducers' phase changes.
Even with these limitations, our work provides valuable insight to device manufacturers, haptic designers and HCI researchers.
Results from our physical model are particularly relevant for PAT manufacturers, as a frst exploration of existing efects afecting the performance of PATs.Our model should be understood more like a high-level tool here and not a replacement for more accurate, low-level, physics-based models.However, even if not as accurate as those potential follow up models, and even if our physical model needs to be tuned for specifc devices, the general trends identifed for these physical efects will most likely hold.
This can help guide future eforts to improve devices.That is, improvements in temperature control will provide linear benefts (i.e., see factor), while benefts from improving driving frequencies are likely to be cubic (i.e., factor).Our model sets a preliminary scene to identify aspects requiring further exploration and where engineering solutions will be most beneftial.
Regarding temperature efects, monitoring would not be enough (i.e., compensating for temperature losses by increasing could lead to faster device heating).Accounting for temperature will require either active cooling or tuning techniques.Improving driving frequencies would require algorithms minimizing phase changes [20] , but these would need to be adapted to also cope with the high update rates required by STM [39,41].
Our results can also be useful to haptic designers using multipoint STM in interactive design.Our results suggest constraining the multi-STM stimuli designed to 10-80 Hz, with stronger efects in the 20 -40 Hz range.The ability to predict and the inter-subject variability (i.e., standard error of ) can guide them to design stronger or more consistent stimuli, reducing the need for trial and error.The derived mapping has also provided the possibility of reverse prediction (achieve designed perceived intensity by fnding required STM parameters to achieve essential physical property), providing a powerful tool for designing haptic interactions that are tailored to user expectations and preferences.
The predictivity in perceived intensity can further enhance the design of diferent level of stifness sensation in UMH [33] from a physical perspective, and precisely provide various pulse profles and thus enhance pulse training surgeries [24,25].From an entertainment perspective, our model can signifcantly enrich the movie-watching experiences by synchronizing haptic feedback to match the intensity of visual content [3].This is particularly advantageous for multi-point stimuli, which enables complex shape rendering, thus allowing for a more advanced visual-haptic content synchronization.Similarly, our model opens opportunities to enhance haptic gaming experiences in the gaming industry by ofering a diverse range of adjustable intensity feedback [14,26].
Relevant implications can also be derived for HCI researchers.Our modelling procedure revealed how humans perceive haptic stimuli (i.e., both physical stimuli and MPT matters) and the physical properties that afect the user's perceptions.This encourages researchers in the haptic domain to follow similar procedures in optimizing their models by fnding the most relevant physical properties.Our model enhanced accurate stimuli delivery and perceptual response predictability.We showed that prior investigations in STM parameters' efects were limited, as some of the diferences found could simply be due to uncontrolled factors (e.g., temperature).These factors often resulted in stimuli that were not physically comparable, leading to variations in pressure/force.Consequently, the observed diferences between stimuli may have been infuenced by these variables rather than the intended parameter variations.Our model, however, incorporates these factors, characterizes their efects, and then uses them to predict users' perceptual responses.
Our fndings generally align with results from prior studies, but our look into the physical properties of the underlying stimuli might point towards alternative explanations to those suggested before.For instance, while we can confrm results from Frier at al. [11] that speeds of 5 to 8 m/s maximize user's perceived intensity (i.e., these speeds correspond to range around 20 -40 Hz in our Study 3), this efect could also be justifed by the fact that devices' output was higher within this range (see Figure 7).Revisiting these results might be required, providing accurate characterization of the physical stimuli that are being compared.Researchers in this space should remain aware of PAT's poor performance at higher frequencies, and characterization of such higher frequency range might need to wait until delivery solutions are found.
In general, ratings seem more dependant on the physical properties of the stimuli, than on perceptual responses.Physical parameters can act as replacements of STM parameters (i.e., ), improving overall accuracy, and the gains from considering are relatively smaller.This does not mean perceptual responses are not important.Our best predictor ( ) depended on values related to users' responses ( ( )), and modelling and still provided benefts.However, as a community, we should avoid biases towards attributing observed efects to perceptual responses and be careful to also assess underlying physical efects.

CONCLUSION
We presented an exploration on the efects of various modulation parameters ( , , ) have on stimuli, both at a physical and perceptual level.We reported on 3 diferent studies, which allow us to derive a two stage model predicting (i.e., average values and standard error), as a combination of a physical and a perceptual model.
Through our physical model, we identifed and characterised 4 physical efects (non-linear efects, limited power, sub-optimum driving frequencies and temperature), being able to predict their efects on acoustic pressure with a relative error of 7.8% and thus, 5 diferent physical properties (i.e., 0 , , 0 , 2.1 and ).
Through our perceptual model, we could confrm how the physical intensity of the stimuli is the main efect infuencing users' , above perceptual parameters such as or .We can confrm that any of the physical properties proposed provides strong correlation with (i.e., supporting assumptions made in prior studies [10,39,55]), with standing as a slightly better predictor for our scope (i.e., 1-3 points, 10-80 Hz).We combine all this in a fnal perceptual model combining all variables and resulting in a prediction error of 8.0% and 8.8% for average and standard error of respectively.We fnally derive recommendations for manufacturers, haptic designers and HCI researchers exploring multi-point STM stimulation.

Figure 2 :
Figure 2: Physical properties of a focal point considered (for multi-point, we use summation of individual points' properties): a) Peak pressure, relative to room pressure ( 0 ); b) Peak pressure, relative to users' ( ); c) Force across 2.1 cm area ( 2.1 ); d) Absolute force across the area perceivable to user ( 0 ); e) Force across perceivable area and relative to users' ( ).

Figure 3 :
Figure 3: Measuring setup used to characterize 0 delivered by the device as STM parameters were varied, featuring a PAT device and a measuring microphone on a moving stage. ) .

Figure 4 :
Figure 4: Characterization of physical efects: a) Efect of on acoustic pressure.b) Efect of on acoustic pressure.c) Efect of Temperature on acoustic pressure.d) Efect of on the acoustic pressure of each focal point at 5 s.

Figure 5 :
Figure 5: The overall experiment setup for study 2.

Figure 6 :
Figure 6: Efect of on user's in Pa.Mean and standard deviation of across users shown in blue.The orange line display the ftting of our ( ) function.

Table 1 :
The efect size of each STM parameter on perceived intensity and distance correlation coefcients between physical properties and perceived intensity under optimum range.

Figure 7 :
Figure 7: Relationship between and per frequency range ( ) and number of points ( ).The orange line shows average and standard error (whiskers) reported by users, while the blue bars represent the , both normalized to a 0-1 range.

Figure 8 :
Figure 8: Relationship between and perceived intensity under from 10 to 80 Hz at 3 conditions.The error bars represent the standard error.The orange line represents the ftted perceptual model.

Table 2 :
The Mean relative error of each model developed