OTUP workflow: target specific optimization of the transmit k‐space trajectory for flexible universal parallel transmit RF pulse design

To optimize transmit k‐space trajectories for a wide range of excitation targets and to design “universal pTx RF pulses” based on these trajectories.


| INTRODUCTION
The great potential of MRI scanners operating at ultra-high field (UHF) (i.e., a B 0 field strength of 7 T and above) has been extensively discussed. 1-3 UHF systems provide higher signal-to-noise ratios, 4 reduced acquisition times and improved spatial resolutions. However, the shorter electromagnetic wavelength at UHF results in inhomogeneity of the RF field. The most promising approach to overcome these issues is the parallel transmission (pTx) technique, as it provides improved control over the spatial and temporal RF field. 5,6 PTx RF pulse design is accompanied by an underlying transmit (or excitation) k-space trajectory spanned by appropriate gradients. Usually the pTx pulse is calculated based on this transmit k-space trajectory. Both are played out simultaneously, resulting in excitation of the spins in a scanned object. For a desired target excitation pattern "the choice of a suitable transmit k-space trajectory is crucial." 7 For non-selective (wholebrain) and slice selective excitation, various studies have optimized the positions of "k T -point" 7,8 and "spoke" 9,10 trajectories together with an appropriate RF pulse. Yip et al. 11 optimized EPI trajectories along with the RF pulse.
RF pulses based on more sophisticated trajectories such as "stack of spirals" (SOS), 12 "concentric-shells" 12 and spiral non-selective (SPINS) 13 trajectories have also been investigated in the past, in both simulations 14 and measurements (in rats). 15 The utilized trajectories in these studies were manually adapted to the field of view (FOV) and voxel size of the target excitation pattern, but the optimality of the chosen transmit k-space trajectories with respect to the excitation target was not evaluated. For instance, in the work of Malik et al. 13 "k max was chosen to cover sufficiently high spatial frequencies to correct for typical B þ 1 variations across the adult brain" at 3 T. RF pulses were calculated based on fixed transmit trajectory without any prior k-space trajectory optimization.
Shao et al. 14 introduced the "k-space trajectory container" principle, which utilized echo-planar and spiral k-space trajectories but stretched them in order to cover a k-space segment optimal for a certain target pattern.
Davids et al. 16 presented the first approach to jointly optimize the RF pulse and the transmit trajectories in the form of 3D shells, SOS and crosses. The authors utilized a very versatile shape parametrization by assigning a set of control points for each k-space trajectory. Connecting these control points resulted in the final trajectory, based on which a pTx pulse was designed and evaluated. While this approach is very flexible, the trajectories are potentially over-parametrized and hence authors reported a respective convergence problem and dependence on the convexity of the individual partly ill defined optimization problems. The authors used the MATLAB (MathWorks, Natick, MA, USA) fmincon function for the trajectory optimization, which may be unsuited for such non-convex optimization problems, as the optimization result strongly depends on the initial guess. Furthermore, the optimized trajectories and resulting pulses were tested only in simulations.
Hao et al. 17 suggested an equally versatile approach with respect to the parametrization of k-space trajectories by using a second-order Bspline basis with consideration of realistic gradient strength and slew rate constraints. The optimization was initialized by SOS, SPINS and k T -point trajectories. For optimization they utilized algorithms that are "monotonically decreasing," for instance algorithms that are implemented in the MATLAB fmincon function. Similarly to the work of Davids et al., 16  Finally, Herrler et al. 18 optimized SPINS trajectories jointly with the RF pulse by using the analytical equations (with small adjustments compared with the original SPINS publication 13 ). The trajectory optimization problem was solved with a global search algorithm, an approach that has a higher probability (compared with the approaches of Davids et al. 16 and Hao et al. 17 ) to find the global optimum. However, their procedure was only tested for whole-brain excitation and not for local excitation (LEx) target patterns. Furthermore, their work did not consider types of k-space trajectory other than SPINS.
Complementary to the previous studies on optimal transmit k-space trajectory design for pTx, in this work we optimized four different basis trajectories suitable for a wide range of pTx RF pulse design targets and used respective analytical expressions 13,19 with a low number of parameters. To ensure robust convergence of this non-convex and global optimization problem a particle swarm solver 20 was used. The respective trajectories were a single variable density spiral-in trajectory 19,21,22 (1SOS), a two stack of variable density spiral-in trajectory 12,14 (2SOS), a three stack of variable density spiral-in trajectory (3SOS) and a SPINS 13 trajectory. We considered the analytical equations of spirals 19 and SPINS 13 in order to create the trajectories. The parameters to be optimized were the parameters of these equations. The procedure was demonstrated for three different target excitation pattern (two LEx and one "whole-brain-like" target excitation target pattern).
During the second part of this study, so-called "universal pTx RF pulses" 23 (UPs) based on the optimized k-space trajectories were designed.
After the study by Herrler et al., 18 this work is only the second study combining sophisticated optimized k-space trajectories with UP design. While Herrler et al. only tested UP design for whole-brain excitation based on SPINS, in the current study we optimized SPINS and spiral trajectories and tested it on LEx and whole-brain-like target patterns.
Traditionally the transmit RF field (B 1 + ) and the static magnetic field (B 0 ) distribution is measured for each subject, and based on this information the pTx RF pulse is calculated during the scan session. The concept of UPs was introduced by Gras et al. 23 to circumvent this inherent lengthy pTx calibration procedure. The major idea of the UP concept is to design RF pulses based on a database of B 0 /B 1 + maps from different subjects prior to pTx scan sessions. Due to the similarity of these maps for most subjects, the UPs also perform well on subjects not contained in the pulse design database.
The UP concept was proved 24,25 and developed 8,26,27 for different sequences, applications and RF coils at 7 T. These works have in common that they present UPs for non-selective or slice selective excitation, based on "k T points" 7 or "spokes." 9 Additionally, a recent study 18 performed a combination of UPs and subsequent further optimization to individual subjects' heads for whole-brain excitation using SPINS trajectories.
Most recently, we published a feasibility study for LEx UPs at 9.4 T. 28 We created LEx UPs based on fixed non-optimized spiral transmit kspace trajectories, which locally excited the visual cortex region of the human brain, while the surrounding tissues experienced only minor excitation.
The UPs in the current study were calculated based on the optimized k-space trajectories. As UPs can be calculated offline while no subject is positioned inside the scanner, additional computation time to optimize the trajectories next to the actual RF pulse is not problematic. The entire OTUP workflow was tested on three different target excitation patterns (two LEx and one "whole-brain-like" target patterns). For one target the resulting RF pulses and transmit trajectories were applied in vivo at 9.4 T for proof of principle. The OTUP code was made available as open source.

| Test target patterns
The following three target excitation patterns ( Figure 1) were created in order to test the workflow.
• "targetNuclei." A small central region encompassing the red nuclei should be excited with a homogenous flip angle (FA) of 7 in eight consecutive transversal slices. The surrounding voxels within these slices should not experience excitation. The voxels directly adjacent to the excited region were excluded from the pulse design mask, to avoid a direct transition from excitation to non-excitation regions, thereby reducing the complexity of the optimization problem. The red nuclei are small structures in the human rostral midbrain involved in motor coordination. 29 This structure is of special interest for research on Parkinson's disease 30,31 and may be of interest for practical reduced FOV application.
F I G U R E 1 The three target FA pattern to test the workflow. From left to right: targetNuclei, targetM and targetWB. The patterns are displayed by a representative slice in all three orientations (from top to bottom: transversal, sagittal, coronal). For the pulse design, only the voxels between the two horizontal red lines were taken into account for each target pattern. Blue, desired non-excitation voxels; yellow, voxels with a desired excitation of 7 ; gray, voxels that were not included either in RF pulse design or in evaluation • "targetM." The letter "M" should be excited with an FA of 7 in one central transversal slice. The surrounding voxels within these slices should experience no excitation. Desired excited and non-excited voxels are directly adjacent. This target pattern served as a "proof of concept" pattern.
• "targetWB." All brain voxels located in 16 consecutive transversal slices should be excited with a homogenous FA of 7 (i.e., there are no voxels that were not allowed to experience excitation). Excitation of the whole brain was not possible as the coverage of the utilized RF coil was not sufficient to reach the upper and lower brain (see B 1 + maps in Supplementary Figure S1).
All target patterns had in common that the transversal slices located above or below the slices of interest (see red horizontal lines in Figure 1) are excluded from the pulse design mask; i.e., the performance of the pulse in these excluded slices is neither evaluated nor optimized by the pulse design algorithm. This approach is feasible since the head to feet direction was chosen as the frequency encoding direction for in vivo imaging.
Thereby, folding artifacts from tissues above and below the red lines would not appear.
For each of the three targets pattern the OTUP workflow consisting of "gradient trajectory optimization" and "final UP calculation" was tested. A workflow chart is visible in Figure 2. The workflow was implemented in MATLAB and is available as open source (https://github.com/ ole1965/workflow_OTUP.git).

| Gradient trajectory optimization
In this section we describe how the transmit k-space trajectories were optimized (upper part of the flow chart in Figure 2). For each basis trajectory its own cost function was implemented. Within each of the cost functions k-space trajectories were created, based on the input parameters to be optimized. Each cost function was minimized by the "particle swarm optimization" 20 to determine the optimal k-space trajectory for a given target excitation pattern. For each intermediate trajectory created during the respective optimization process an RF pulse was designed and its excitation performance was evaluated against the target profile. In order to save computation time, we optimized each k-space trajectory for only one head from the design database. Optimizing the trajectory on all database subjects simultaneously would have been possible but time consuming. However, Supplementary Figure S2 depicts a comparison of k-space trajectories that were optimized based on one head and multiple heads, respectively, for targetNuclei. The performance advantages of the trajectory optimized using multiple heads is only very small compared with the trajectory optimized using data from one head. The optimization result for the trajectory strongly depends on the target pattern, while the B 0 /B 1 + F I G U R E 2 OTUP flow chart for the 1SOS trajectory: The chart needs to be read from top to bottom. The first step is to optimize the trajectory for a given target excitation pattern. Therefore, a cost function to be optimized with the particle swarm optimization is introduced. The input parameter of this cost function is the analytical equation parameter of a spiral (see Equation (1)). Within the cost function, a trajectory is calculated by means of the analytical equation. For this trajectory a TP is designed and evaluated. The output parameter of the cost function is the RMSE. In the next step a UP is designed based on the optimized trajectory from before. The UP design consists of MLS optimization and the application of the active-set algorithm maps have only a minor influence. For this reason, the k-space trajectory optimization was only done for one representative subject during this study.
After trajectory optimization the UPs were designed on B 0 /B 1 + maps from the entire design database using the resulting optimal k-space tra- is the radial extent of the spiral, α ℝ > 0 is the density variable and ω ¼ 2πn, with n the number of spiral turns.
For 1SOS Equation (1 was implemented in a cost function having the equation parameters T, λ, α, n and the position z of the spiral on the k z axis as input parameters. Within the cost function, the spiral k was created based on these input parameters. For this trajectory, an RF pulse was designed based on the B 0 /B 1 + map of one subject. For RF pulse design the "spatial domain method" 32 was utilized because of its speed, due to the small tip angle approximation 33 of the Bloch equations. Herein, only the least squares problem was solved using the MATLAB lsqr function. 34 In contrast, once the optimal trajectory was determined, RF pulses will be optimized by solving the "magnitude least squares" (MLS) 35 problem (see the "Final UP calculation" section 2.3). The output parameter of the cost function was the root mean squared error (RMSE) between the FA profile of that pulse (Bloch simulated) and the target excitation pattern.
For 2SOS and 3SOS the procedure was analogous, with the only exception that two or three spirals were designed, respectively. For 2SOS the cost-function input parameters were T 1 ,λ 1 , α 1 , n 1 for the first spiral and T 2 ,λ 2 , α 2 , n 2 for the second spiral, as well as z 1 and z 2 as the positions on the k z axis. For 3SOS a third spiral with corresponding input parameters for the cost function was added.
A SPINS trajectory 13 can be created using Equations (2)-(4): where k r is the radial coordinate and k θ and k φ are the polar and azimuthal angles. The parameter k max ℝ > 0 is the maximum radial extent, u ℝ > 0 and v ℝ > 0 are the respective angular velocities and T ℕ is the length of the trajectory. Additional parameters are α, β ℝ > 0 . Similar to 1SOS, 2SOS and 3SOS a cost function with the SPINS parameters (k max , α, β, u, v, T) as input parameters was implemented in order to create a SPINS trajectory for each set of parameters and to design a pulse based on this trajectory.
Trajectories that exceeded the maximum gradient amplitude or slew rate were punished with a very high cost-function output value, in order to prevent the optimization from considering these trajectories as optimum. The maximum allowed trajectory duration was 10 ms. Reasonable lower and upper bounds for the input parameters were found empirically. Each basis k-space trajectory was optimized separately. From the four optimized trajectories the one that enabled the lowest RMSE performance or a considerably shorter duration with similar RMSE was considered for UP calculation.
To solve the four optimization problems "particle swarm optimization," 20 implemented in the MATLAB particleswarm function, was utilized.
Particle swarm optimization "works by maintaining a swarm of particles that move around in the search-space influenced by the improvements discovered by the other particles." 36 This optimization routine was chosen because it creates a reliable overview of the entire search space and it does not require an initial guess, which would be hard to find in this case. As this algorithm is stochastic and non-deterministic the calculated optimum of the particleswarm function can differ from one call to another. For this reason, the particleswarm function was executed 20 times for each cost function mentioned above. From the 20 optimized parameter sets, the one that provided the lowest RMSE was utilized for further RF pulse calculation. If several parameter sets provided equal or highly similar RMSE values, the one that produced the shortest trajectory duration was used further.
In order to demonstrate the influence of the underlying transmit trajectory on the RF pulse performance, we tested each of the optimized trajectories for each of the test targets. For instance, for targetWB, we calculated pulses based on the optimal trajectory for targetNuclei, targetM and targetWB. The pulses were calculated and applied on one example subject using Bloch simulations.

| Final UP calculation
Based on the optimized transmit trajectories, UPs were designed analogously to our own preliminary work 28 (lower part of the flow chart in Figure 2). The following optimization problem was considered in order to create a pulse p Ã UP that excites the same target pattern on multiple heads (i.e., the UP): where A full,j , with j ¼ 1,…, N DB , is the full system information matrix 32 of subject j and N DB is the size of the design database, i.e., the number of subject heads based on which the pulse will be designed. p is the RF pulse and m tar is the target excitation pattern. This approach is an extension of the spatial domain method. 32 For each test target and optimized trajectory, we solved this optimization problem in Equation (5) with a combination of the MLS optimization 35 and the active-set algorithm implemented in the MATLAB fmincon function.
The MLS optimization relies (similarly to the spatial domain method) on the small tip angle approximation 33 but considers only the resulting magnetization profile's magnitude in the optimization process. The spatial domain method used within the gradient optimization was drastically faster (because it only solves the least squares optimization problem) but includes the profile's phase and magnitude in the optimization. However, the profile's phase was of no interest during this study. Considering only the profile's magnitude and not the phase in the pulse design results in an improved pulse performance for applications where the profile's phase is of no relevance. 35 Similarly to our own preliminary work, 28 the MLS result was utilized as the initial guess for the active-set algorithm, in order to potentially further improve the pulse performance.
Within the cost function minimized by the active-set algorithm the FA profile for a given pulse was simulated with Bloch equations (without relaxation). The output parameter of the cost function was the RMSE between the magnitude of the resulting FA profile and the target FA pattern.
In line with the literature, 23-26 the active-set algorithm was chosen because of its speed and robustness. 8,37 With the fmincon function the solution was constrained to a maximum pulse amplitude of 130 V at plug level (hardware limit). The optimization was stopped if further improvement on 30 consecutive iterations was negligible.
The final UPs were simulated on seven non-database heads using Bloch equations. Their performances were compared with the corresponding tailored pTx RF pulse (TP) performances. The TPs were calculated with the MLS method only.
The UP design database consisted of the B 0 /B 1 + maps from 11 different subjects. The 11 datasets were randomly chosen from the 18 overall acquired B 0 /B 1 + maps at 9.4 T. In the supplementary material (see caption of Figures S3 and S4) we provide an analysis of database size, where we observe that UP performance tends to converge when the size of the database is greater than 11 (see caption of Figures S3 and S4).
All calculations (gradient trajectory optimization, analysis of database size-Supporting Information-and final UP calculation) were performed on a high-performance-compute system node equipped with an Intel Xeon Ice Lake Platinum 8360Y processor (256 GB RAM, 72 cores with 2.4 GHz each) exploiting parallel computing. Tissue masks were created using a masking routine based on a neuronal network. 28 The final UP for targetNuclei was applied in vivo at 9.4 T as a proof of principle example. In order to test the in vivo performance for tar-getNuclei, the calculated pulses were applied in a T 2 *-weighted 3D short T R spoiled gradient echo (GRE) sequence (voxel size 0.8 Â 0.8 Â 0.8 mm 3 , matrix size 224 (left to right) Â 280 (anterior to posterior) Â 280 (head to feet), 3D encoding direction left to right, phase encoding direction anterior to posterior, frequency encoding direction head to feet, T R = 18 ms, T E = 8 ms, BW = 260 Hz/px, GRAPPA 2 Â 2).

| Measurement system and data acquisition
Global and local SAR were supervised using the "virtual observation point" (VOP) method. 41,42 For VOP creation we considered the two virtual family multi-tissue models "Duke" and "Ella" 43 for two slightly different head positions inside the RF coil, respectively. The final VOP file contained 494 VOPs. The maximum local SAR limit over a 10 s window was 20 W/kg; over a 6 min window the limit was 10 W/kg. The corresponding global SAR limits were 6.4 W/kg and 3.2 W/kg.
The dwell time for RF pulses and gradient shapes was 10 μs. The scanner inherent gradient delay of 4 μs was taken into account in the pulse files. To measure the gradient delay, an RF pulse was designed for a simple target excitation pattern (i.e., a vertical line in a transversal slice from anterior to posterior direction) in a phantom experiment. The underlying transmit k-space trajectory was a 1SOS trajectory. For this pulse and trajectory different gradient delay times were tested. For 4 μs the orientation of the target excitation pattern in the phantom was equal to the intended orientation (i.e., the vertical line). For delay times below or above 4 μs, the vertical line was rotated slightly because of the underlying spiral trajectory.
3 | RESULTS Figure 3 shows the optimized transmit k-space trajectories. It should be noted that all plots have the same k x and k y axis limits; however, for the sake of clarity the limits of the k z axes of the SPINS plots differs from the limits of the k z axes of the SOS plots. The RMSE values next to each plot are based on the difference between the target FA profile and the respective Bloch simulated FA profile (i.e., the unit of the RMSE is degrees).

| Gradient trajectory optimization
For targetNuclei the 1SOS optimization resulted in a 9.44 ms long spiral located on À3/m on the k-space z-axis. The spiral had an extent of approximately ±80/m on the k x and k y axes. The first of the two spirals from the 2SOS result (6.5 ms) had the same extent. The extent of the second spirals is approximately half that of the first spiral. One spiral of the 3SOS result (8.37 ms) had an extent of ±100/m in k x and k y directions.
Next to a second, smaller spiral, the optimization created an additional very short third spiral shaped like an arc or a half-turn, rather than a full spiral. The SPINS optimization resulted in a trajectory (9.38 ms) similar to a concentric 3D shell. 15 As the 2SOS trajectory enables the shortest duration (6.5 ms), while having a similar RMSE to the other trajectories (0.26), it was chosen for UP calculations for targetNuclei.
F I G U R E 3 Optimized transmit trajectories. In each column, the same basis trajectory (1SOS, 2SOS, 3SOS or SPINS) was utilized. Each row presents the results for one of the three test target patterns. For each optimized trajectory the duration of the trajectory (ms) and the RMSE ( ) that the trajectories enable is displayed above the images. The green boxes present the trajectories that were set as optimum for each target The trajectories for targetM extended further into the outer k-space compared with the targetNuclei trajectories. The 1SOS and the 2SOS results were very similar, having one large spiral reaching ±150/m on the k x and k y axes (one spiral of the 2SOS result is only an arc). The 3SOS result consisted of a large low density spiral having the same k x and k y extents as 1SOS and 2SOS and a smaller spiral underneath. The third spiral, again, was only an arc. The SPINS result had similar extents on the k x and k y axes, but a large extent on the k z axis (from 0/m to À400/m).
Again, all four results were similarly well suited for RF pulse design; i.e., all four RMSE values were between 0.15 and 0.17. All trajectories exploited the maximum allowed duration of 10 ms. The SPINS result was considered for further pulse calculations for targetM because it showed the best performance.
For targetWB all trajectories exhibited a drastically lower extent in k-space compared with the targetNuclei and targetM trajectories. Additionally, the durations are shorter compared with the durations from the other trajectories. As the 3SOS result constituted the best compromise between duration (3.96 ms) and RMSE (0.36) this trajectory was chosen for further pulse calculations for targetWB.
The calculation time for each of these k-space trajectory optimizations was below 24 h (depending on the target, Figure 1).
As Figure 4 shows, for each target pattern the associated optimal transmit trajectory enabled the best performing RF pulse. However, if the underlying trajectory was not designed for the respective target pattern, the RF pulse performance diminished. This effect is most pronounced if F I G U R E 4 Bloch simulated FA profiles of RF pulses designed for different target patterns and different transmit trajectories. In the first row subject specific TPs aiming at excitation of targetNuclei were applied. Each pulse was designed based on another trajectory (depicted in the columns). In the first column, the underlying trajectory was the optimized trajectory for targetNuclei. In the second (third) column, the underlying trajectory was the optimized trajectory for targetM (targetWB). Analogously, in rows two and three TPs were applied aiming at excitation of targetM and target WB, respectively. Below each profile, the corresponding RMSE and the maximum amplitude of each TP is displayed the optimal trajectory for the targetM pattern is applied to the targetWB pattern and vice versa. The maximum pulse voltages were relatively low if the RF pulse was designed based on the trajectory designed for a certain target. Utilizing a suboptimal trajectory can lead to dramatically higher pulse amplitudes (for instance, 586 V for targetWB using the targetM trajectory, 24 V using the targetWB trajectory).

| Final UP calculation
For targetNuclei, the fmincon function was stopped after 395 iterations ($9 d). Figure 5 presents a comparison between the simulated performances of the UP and the corresponding TPs on seven non-database heads. The TPs slightly outperformed the UP. The computation time of each TP was about 24 s (pulse length 6.5 ms). The FA profiles resulting from the TPs exhibited a uniform excitation of the desired excitation area with mean FAs close to 7 . The desired non-excitation areas were in some voxels only very slightly excited, with mean FAs close to 0 . The mean RMSE performance of the TPs (UP) was 0.25 (0.52). For the UP the mean FAs in the nuclei area were between 6 and 7 (except for Head 14) and the excitation in the desired non-excitation areas was close to 0 .
Analogously, Figure 6 displays the TP and UP performances for targetM. The mean TP (UP) RMSE performance was 0.14 (0.53). Each TP created an excitation perfectly shaped like an "M." For the UP (fmincon function was stopped after 820 iterations, $6 d) the transition between desired excitation and non-excitation voxels was not as sharp as for the TPs (computation time per TP was approximately 5 s, pulse length 9.92 ms). The mean UP excitation within the "M" is mostly between 6 and 7.5 . However, for Head 14 the mean excitation was 5.3 .
In Figure 7 pulse performances for targetWB are depicted. The targetNuclei UP performance was satisfactory on all non-database heads. It provided the lowest difference between the mean RMSE TP value and the mean RMSE UP value (0.27). For this reason, the pulses for targetNuclei were tested in vivo at 9.4 T.
In Figure

| DISCUSSION
This work presents the OTUP workflow to optimize the underlying transmit k-space trajectory for pTx UPs and through this further exploit the advantages of the UP concept.
Because for UP design there is no demand for a fast online pulse calculation procedure while the subject is positioned inside the scanner, the underlying transmit k-space trajectories were optimized to achieve the best possible match with the excitation target before the actual UP design process started. For each of three test target excitation patterns (Figure 1), four basis transmit trajectories were optimized. The optimization time per trajectory was up to 24 h. The optimized trajectories were also utilized for subject specific TPs designed herein for comparison.
As expected, with increasing complexity of the target excitation pattern the extents of the optimized transmit trajectories in k x and k y directions in k-space increased (Figure 3). For the least complex, targetWB (i.e., excitation of 16 entire transversal brain slices), the optimized trajectory ran close to the k-space center, covering only low frequencies. For this reason, the targetWB trajectory is unsuitable for excitation of targetM (i.e., the most complex target pattern, exciting only the voxels in an area shaped like the letter "M") ( Figure 4). Interestingly, none of the optimized trajectories for targetWB exploited the maximum allowed duration of 10 ms. This leads to the conclusion that shorter trajectories have an advantage for "whole-brain-like" targets such as targetWB. RF pulse calculations for longer trajectories possibly lead to over-optimization, which decreases the pulse performance.
For targetM, a very sharp transition between excitation and non-excitation between adjacent voxels needed to be achieved. The optimized targetM trajectory covered high frequencies and expanded till ±150/m. This resulted in a k-space FOV of 300/m in k x and k y directions. Converting this value to the voxel size in the image space (i.e., calculating the reciprocal of the k-space FOV) results in a utilized voxel size of the B 1 + maps of 3.5 mm approximately. As visible in Figure 4, such high frequencies (and possibly such long trajectories) are unsuited for excitation of targetWB.
F I G U R E 8 A, B, Two different subjects. The structure is the same in both boxes. The images in the left-hand column show the pulse design mask with the target excitation area (yellow) and non-excitation areas (blue) with a representative transversal, sagittal and coronal slice for targetNuclei. The black-white images are GRE acquisitions utilizing the UP and the respective TP for targetNuclei on two non-database head.
Only the voxels between the two red horizontal lines on the sagittal and coronal images were taken into account for the pulse design. The numbers above depict the mean signal and standard deviation in desired excitation and non-excitation regions, calculated by means of the masks on the left. The GRE images are expressed in normalized arbitrary units; i.e., the value of each voxel was divided by the maximum signal value in the whole acquisition As some of the spirals in the SOS trajectory optimization reduced to arcs (i.e., the shortest possible spiral) the number of spirals could be included as another parameter to be optimized in future optimizations.
Interestingly, the SPINS trajectories had completely different shapes for the three tested excitation targets. Although originally thought to be well suited for entire brain slices, 13 the proposed optimization method created a SPINS trajectory that was also well suited for the complex LEx targetM. In general, SPINS tend to have a longer duration for similar RMSEs as SOS trajectories.
The transmit k-space trajectory optimization was performed using the small tip angle approximation 32,33 of the Bloch equations within the cost function (see flow chart in Figure 2). In theory, this optimization could also have been performed utilizing the exact Bloch equations (as done for the UP calculation). For that case (i.e., tailored RF pulse design for one head) one function call would need at least 1 d (with parallelization). As the particle swarm algorithm called the cost function at least 6000 times, using the Bloch equations for transmit k-space trajectory optimization would take more than 16 years and was hence not possible during this study.
While in the current work the underlying transmit trajectories were optimized, Shao et al., 14 Malik et al. 13 and numerous other pulse design studies 15,32,35,44,45 have in common that the respective transmit trajectories were only adapted to a given target pattern, FOV and matrix size, rather than performing an actual trajectory optimization.
An actual optimization of the k-space trajectory was performed for spoke 9,10,46 and k T -point trajectories, 8,23 where the locations of the spokes and the k T points were optimized simultaneously with the RF pulse for non-selective or slice selective excitation. Additionally, Davids et al., 16 Hao et al. 17 and Herrler et al. 18 were the first to optimize more sophisticated transmit trajectories such as SOS, SPINS, shells and cross.
The latter works provide normalized RMSE (NRMSE) values between 6% and 7% for whole-brain excitation. In order to create and optimize the k-space trajectories SOS, shells and cross, Davids et al. 16 assigned a set of up to 50 control points for each set of certain shape parameters for the respective basis trajectory. These control points needed to be connected afterwards with additional software. 47 In contrast, in the current work, the trajectories were created by analytical equations for given parameters, which is more straightforward to implement than the approach by Davids et al. due to the lower number of parameters required.
Davids et al. 16 Table T1 in the Supporting Information provides the RMSE values of the 20 initial guesses and the RMSE values of the corresponding active-set results. It can be seen that none of the active-set results provides an RMSE as small as the RMSE delivered after 20 particle swarm executions during this work (0.26 for 2SOS and tar-getNuclei). Furthermore, even if the particle swarm optimum was chosen as the initial guess for the active-set algorithm, the result could not be improved (see last row of Table T1 in the Supporting Information). Both investigations emphasize the benefits of the particle swarm algorithm. It creates a reliable overview of the entire search-space and it does not rely on the initial guess, hence making it more appropriate for k-space trajectory optimization.
Herrler et al. 18 also performed a transmit k-space trajectory optimization based on the analytical equations with a global search optimization routine. However, the authors tested the procedure only for SPINS and whole-brain excitation and fixed the pulse length to 1 ms.
For each target tested in the current study, the optimal trajectory was utilized to calculate UPs. Considering the UP results for targetWB it was conspicuous that the UP performance on Head 14 is drastically lower than the performance on the remaining non-database heads. It should be noted that none of the presented performance values from the literature are directly comparable to the values from the present work, as the utilized RF coils, target pattern and field strengths differed. For instance, all mentioned studies on whole-brain excitation aimed to excite the entire brain homogeneously, while in the current study the very lower and upper parts of the brain needed to be excluded from the design mask ( Figure 1) for RF coil coverage reasons (Supplementary Figure S1).
Considering the flow chart in Figure 2, a possible starting point to modify and potentially improve the OTUP workflow may be the utilization of another output parameter for the k-space trajectory optimization cost function. Instead of calculating the RMSE it should easily be possible to implement other performance measures, such as the summed absolute values of the difference between resulting profile and target profile for each voxel. Another exchangeable module inside OTUP could be the "active-set" algorithm that uses the MLS optimization result as the initial guess. Herein, other optimization methods, such as the "interior-point" method, 48 could be applied. Furthermore, it may be worthwhile to investigate initial guesses other than the MLS solution. Another potential improvement of the OTUP workflow may be consideration of the "gradient impulse response function" (GIRF) 49 or RF amplifier performance 50 or direct consideration of SAR constraints during the design. 51 The GRE acquisitions in Figure 8 demonstrate the power of the UP concept with optimized transmit trajectories: although the depicted subjects were non-database subjects, the UP had a highly similar in vivo performance at 9.4 T to the subject specific TPs. Before the TP could be applied, subject specific B 0 /B 1 + mapping and pulse calculation needed to be performed, while the UP could be applied immediately. The presented UP design software is fully functional.
In all in vivo images (Figure 8), the red nuclei are excited. However, for UP and TPs, there is some excitation in areas where it was not desired and also not predicted by the simulation results. Especially in the TP simulations ( Figure 5), this excitation was not visible at all.
The most likely reason for this effect is experimental inaccuracies within the MRI scanner system. One possibility may be imperfections in the execution of the transmit gradient waveform, caused by spatial non-linearity of the gradient system, 52 suboptimal pre-emphasis calibration, 53 mechanical vibrations after gradient switching 54 or thermal variations in the hardware. 55 To solve this issue, a gradient impulse response function 49 of the individual MRI scanner's gradient system could be incorporated in the k-space trajectory optimization procedure and the RF pulse design in order to match the theoretical and actual gradient waveforms. Only incorporating the scanner inherent gradient delay of 4 μs, as was done during this study, is apparently not sufficient.
Another error source may be inaccuracies in the transmission of the RF pulses. Especially in this study, where the RF pulses mostly have a "noise-like" shape (see Supplementary Figure S5), the RF amplifiers may not be capable of transmitting the exact desired RF shape. Monitoring and correcting 50,56,57 the RF pulses could be one solution for this issue. Another possibility constitutes "Minimum envelope roughness pulse design," 58 an approach "that yields pulses with smoother envelopes." Further reasons for the difference between simulation and in vivo results could be B 0 /B 1 + map inaccuracies 39,59 or motion artifacts. 60 Deviations between simulation and human in vivo results for LEx pulses are a problem that can also be found in the literature. 61 Future investigations need to incorporate the hardware inaccuracies into the trajectory optimization and the RF pulse design, in order to bring simulation and in vivo results in line. As the problems mentioned above could not be addressed during this work, reduced FOV measurements of the red nuclei showed aliasing artefacts (see Supplementary Figure S6). Incorporating the above discussed solutions into the herein presented workflow is a subject of future work.

| CONCLUSION
We introduced and validated a new pTx pulse design method that optimizes the transmit k-space trajectory (SOS and SPINS) to best match the excitation target of interest prior to calculation of respective UPs or subject specific TPs. OTUP was demonstrated for three different target excitation patterns. The importance of utilizing an optimized target dependent k-space trajectory was emphasized.
For one of the target patterns the resulting pulses were tested in vivo at 9.4 T as proof of principle for the integrity of the design method.
TPs and UP produced highly similar GRE images, proving the power of the UP concept.
The workflow code and the B 0 /B 1 + map data from 18 subjects measured at 9.4 T are available online as open source.

SUPPORTING INFORMATION
Additional supporting information may be found in the online version of the article at the publisher's website.