* [135] Vincent Duchaine and Clement Gosselin.  Safe, stable and intuitive control for physical human-robot interaction.  In 2009 IEEE International Conference on Robotics and Automation, pages 3383–3388, 2009.
  * [136] Zhijun Li, Bo Huang, Zhifeng Ye, Mingdi Deng, and Chenguang Yang.  Physical human–robot interaction of a robotic exoskeleton by admittance control.  IEEE Transactions on Industrial Electronics, 65(12):9614–9624, 2018.
  * [137] Chiara Talignani Landi, Federica Ferraguti, Lorenzo Sabattini, Cristian Secchi, and Cesare Fantuzzi.  Admittance control parameter adaptation for physical human-robot interaction.  In 2017 IEEE International Conference on Robotics and Automation (ICRA), pages 2911–2916, 2017.
  * [138] Hsieh-Yu Li, Ishara Paranawithana, Liangjing Yang, Terence Sey Kiat Lim, Shaohui Foong, Foo Cheong Ng, and U-Xuan Tan.  Stable and compliant motion of physical human–robot interaction coupled with a moving environment using variable admittance and adaptive control.  IEEE Robotics and Automation Letters, 3(3):2493–2500, 2018.
  * [139] Christopher Yee Wong, Saeid Samadi, Wael Suleiman, Abderrahmane Kheddar, and Christopher Yee Wong.  Touch semantics for intuitive physical manipulation of humanoids.  IEEE transactions on human-machine systems., 52(6), 2022-12.
  * [140] Alexis C. Holgado, Nicola Piga, Tito Pradhono Tomo, Giulia Vezzani, Alexander Schmitz, Lorenzo Natale, and Shigeki Sugano.  Magnetic 3-axis soft and sensitive fingertip sensors integration for the icub humanoid robot.  In 2019 IEEE-RAS 19th International Conference on Humanoid Robots (Humanoids), pages 1–8, 2019.
  * [141] Alexis Carlos Holgado, Tito Pradhono Tomo, Sophon Somlor, and Shigeki Sugano.  A multimodal, adjustable sensitivity, digital 3-axis skin sensor module.  Sensors, 20(11), 2020.
  * [142] Alessandro Albini and Giorgio Cannata.  Pressure distribution classification and segmentation of human hands in contact with the robot body.  The International Journal of Robotics Research, 39(6):668–687, 2020.
  * [143] Mattia Leonori, Juan M. Gandarias, and Arash Ajoudani.  Moca-s: A sensitive mobile collaborative robotic assistant exploiting low-cost capacitive tactile cover and whole-body control.  IEEE Robotics and Automation Letters, 7(3):7920–7927, 2022.
  * [144] P. Mittendorfer, E. Yoshida, and G. Cheng.  Realizing whole-body tactile interactions with a self-organizing, multi-modal artificial skin on a humanoid robot.  Advanced Robotics, 29(1):51–67, 2015.
  * [145] Simon Armleder, Emmanuel Dean-Leon, Florian Bergner, and Gordon Cheng.  Interactive force control based on multimodal robot skin for physical human- robot collaboration.  Advanced Intelligent Systems, 4(2):2100047, 2022.
  * [146] Vincent Duchaine, Nicolas Lauzier, Mathieu Baril, Marc-Antoine Lacasse, and Clement Gosselin.  A flexible robot skin for safe physical human robot interaction.  In 2009 IEEE International Conference on Robotics and Automation, pages 3676–3681, 2009.
  * [147] Marc Teyssier, Brice Parilusyan, Anne Roudaut, and Jürgen Steimle.  Human-like artificial skin sensor for physical human-robot interaction.  In 2021 IEEE International Conference on Robotics and Automation (ICRA), pages 3626–3633, 2021.
  * [148] Shan Wei, Yijian Liu, Lina Yang, Haicheng Wang, Haoran Niu, Chao Zhou, Yanyan Wang, Qiuquan Guo, and Da Chen.  Flexible large e-skin array based on patterned laser-induced graphene for tactile perception.  Sensors and Actuators A: Physical, 334:113308, 2022.
  * [149] Teng Xue, Weiming Wang, Jin Ma, Wenhai Liu, Zhenyu Pan, and Mingshuo Han.  Progress and prospects of multimodal fusion methods in physical human–robot interaction: A review.  IEEE Sensors Journal, 20(18):10355–10370, 2020.
  * [150] Emanuele Magrini, Fabrizio Flacco, and Alessandro De Luca.  Estimation of contact forces using a virtual force sensor.  In 2014 IEEE/RSJ International Conference on Intelligent Robots and Systems, pages 2126–2133, 2014.
  * [151] Arne Wahrburg, Björn Matthias, and Hao Ding.  Cartesian contact force estimation for robotic manipulators-a fault isolation perspective.  IFAC-PapersOnLine, 48(21):1232–1237, 2015.
  * [152] Yazhan Zhang, Guanlan Zhang, Yipai Du, and Michael Yu Wang.  Vtacarm. a vision-based tactile sensing augmented robotic arm with application to human-robot interaction.  In 2020 IEEE 16th International Conference on Automation Science and Engineering (CASE), pages 35–42, 2020.
  * [153] Hang Su, Wen Qi, Zhijun Li, Ziyang Chen, Giancarlo Ferrigno, and Elena De Momi.  Deep neural network approach in emg-based force estimation for human–robot interaction.  IEEE Transactions on Artificial Intelligence, 2(5):404–412, 2021.
  * [154] Stavros Grafakos, Fotios Dimeas, and Nikos Aspragathos.  Variable admittance control in phri using emg-based arm muscles co-activation.  In 2016 IEEE International Conference on Systems, Man, and Cybernetics (SMC), pages 001900–001905, 2016.
  * [155] DSV Bandara, Jumpei Arata, and Kazuo Kiguchi.  A noninvasive brain–computer interface approach for predicting motion intention of activities of daily living tasks for an upper-limb wearable robot.  International Journal of Advanced Robotic Systems, 15(2):1729881418767310, 2018.
  * [156] Luis Roda-Sanchez, Celia Garrido-Hidalgo, Arturo S García, Teresa Olivares, and Antonio Fernández-Caballero.  Comparison of rgb-d and imu-based gesture recognition for human-robot interaction in remanufacturing.  The International Journal of Advanced Manufacturing Technology, pages 1–13, 2021.
  * [157] Joseph Campbell and Katsu Yamane.  Learning whole-body human-robot haptic interaction in social contexts.  In 2020 IEEE International Conference on Robotics and Automation (ICRA), pages 10177–10183. IEEE, 2020.
  * [158] Mustafa Can Bingol and Omur Aydogmus.  Performing predefined tasks using the human–robot interaction on speech recognition for an industrial robot.  Engineering Applications of Artificial Intelligence, 95:103903, 2020.
  * [159] K Ashok, Mohd Ashraf, J Thimmia Raja, Md Zair Hussain, Dinesh Kumar Singh, and Anandakumar Haldorai.  Collaborative analysis of audio-visual speech synthesis with sensor measurements for regulating human–robot interaction.  International Journal of System Assurance Engineering and Management, pages 1–8, 2022.
  * [160] Shuo Gao, Yanning Dai, and Arokia Nathan.  Tactile and vision perception for intelligent humanoids.  Advanced Intelligent Systems, 4(2):2100074, 2022.
  * [161] Emanuele Magrini, Fabrizio Flacco, and Alessandro De Luca.  Control of generalized contact motion and force in physical human-robot interaction.  In 2015 IEEE International Conference on Robotics and Automation (ICRA), pages 2298–2304, 2015.
  * [162] Don Joven Agravante, Andrea Cherubini, Antoine Bussy, Pierre Gergondet, and Abderrahmane Kheddar.  Collaborative human-humanoid carrying using vision and haptic sensing.  In 2014 IEEE international conference on robotics and automation (ICRA), pages 607–612. IEEE, 2014.
  * [163] Don Joven Agravante, Andrea Cherubini, Antoine Bussy, Pierre Gergondet, and Abderrahmane Kheddar.  Collaborative human-humanoid carrying using vision and haptic sensing.  In 2014 IEEE International Conference on Robotics and Automation (ICRA), pages 607–612, 2014.
  * [164] H. Kawamoto, Suwoong Lee, S. Kanbe, and Y. Sankai.  Power assist method for hal-3 using emg-based feedback controller.  In SMC’03 Conference Proceedings. 2003 IEEE International Conference on Systems, Man and Cybernetics. Conference Theme - System Security and Assurance (Cat. No.03CH37483), volume 2, pages 1648–1653 vol.2, 2003.
  * [165] Kai Gui, Honghai Liu, and Dingguo Zhang.  Toward multimodal human–robot interaction to enhance active participation of users in gait rehabilitation.  IEEE Transactions on Neural Systems and Rehabilitation Engineering, 25(11):2054–2066, 2017.
  * [166] Ali Al-Yacoub, YC Zhao, William Eaton, Yee Mey Goh, and Niels Lohse.  Improving human robot collaboration through force/torque based learning for object manipulation.  Robotics and Computer-Integrated Manufacturing, 69:102111, 2021.
  * [167] Lourdes Martínez-Villaseñor and Hiram Ponce.  A concise review on sensor signal acquisition and transformation applied to human activity recognition and human–robot interaction.  International Journal of Distributed Sensor Networks, 15(6):1550147719853987, 2019.
  * [168] Zhe Cao, Tomas Simon, Shih-En Wei, and Yaser Sheikh.  Realtime multi-person 2d pose estimation using part affinity fields.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 7291–7299, 2017.
  * [169] Yu Cheng, Bo Yang, Bo Wang, and Robby T Tan.  3d human pose estimation using spatio-temporal networks with explicit occlusion training.  In Proceedings of the AAAI Conference on Artificial Intelligence, volume 34, pages 10631–10638, 2020.
  * [170] Kyoungoh Lee, Inwoong Lee, and Sanghoon Lee.  Propagating lstm: 3d pose estimation based on joint interdependency.  In Proceedings of the European conference on computer vision (ECCV), pages 119–135, 2018.
  * [171] Ching-Hang Chen, Ambrish Tyagi, Amit Agrawal, Dylan Drover, Rohith Mv, Stefan Stojanov, and James M Rehg.  Unsupervised 3d pose estimation with geometric self-supervision.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 5714–5724, 2019.
  * [172] Z. Cao, G. Hidalgo, T. Simon, S. Wei, and Y. Sheikh.  Openpose: Realtime multi-person 2d pose estimation using part affinity fields.  IEEE Transactions on Pattern Analysis & Machine Intelligence, 43(01):172–186, jan 2021.
  * [173] Jan Docekal, Jakub Rozlivek, Jiri Matas, and Matej Hoffmann.  Human keypoint detection for close proximity human-robot interaction.  In 2022 IEEE-RAS 21st International Conference on Humanoid Robots (Humanoids), pages 450–457, 2022.
  * [174] Kenko Fujii, Gauthier Gras, Antonino Salerno, and Guang-Zhong Yang.  Gaze gesture based human robot interaction for laparoscopic surgery.  Medical image analysis, 44:196–214, 2018.
  * [175] Oriane Dermy, François Charpillet, and Serena Ivaldi.  Multi-modal intention prediction with probabilistic movement primitives.  In Human Friendly Robotics: 10th International Workshop, pages 181–196. Springer, 2019.
  * [176] Osama Mazhar, Sofiane Ramdani, Benjamin Navarro, Robin Passama, and Andrea Cherubini.  Towards real-time physical human-robot interaction using skeleton information and hand gestures.  In 2018 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS), pages 1–6, 2018.
  * [177] Osama Mazhar, Benjamin Navarro, Sofiane Ramdani, Robin Passama, and Andrea Cherubini.  A real-time human-robot interaction framework with robust background invariant hand gesture detection.  Robotics and Computer-Integrated Manufacturing, 60:34–48, 2019.
  * [178] Luis Roda-Sanchez, Teresa Olivares, Celia Garrido-Hidalgo, José Luis de la Vara, and Antonio Fernandez-Caballero.  Human-robot interaction in industry 4.0 based on an internet of things real-time gesture control system.  Integrated Computer-Aided Engineering, 28(2):159–175, 2021.
  * [179] Francesco Romano, Gabriele Nava, Morteza Azad, Jernej Čamernik, Stefano Dafarra, Oriane Dermy, Claudia Latella, Maria Lazzaroni, Ryan Lober, Marta Lorenzini, Daniele Pucci, Olivier Sigaud, Silvio Traversaro, Jan Babič, Serena Ivaldi, Michael Mistry, Vincent Padois, and Francesco Nori.  The codyco project achievements and beyond: Toward human aware whole-body controllers for physical human robot interaction.  IEEE Robotics and Automation Letters, 3(1):516–523, 2018.
  * [180] Lorenzo Vianello, Jean-Baptiste Mouret, Eloise Dalin, Alexis Aubry, and Serena Ivaldi.  Human posture prediction during physical human-robot interaction.  IEEE Robotics and Automation Letters, 6(3):6046–6053, 2021.
  * [181] Jessica Lanini, Hamed Razavi, Julen Urain, and Auke Ijspeert.  Human intention detection as a multiclass classification problem: Application in physical human–robot interaction while walking.  IEEE Robotics and Automation Letters, 3(4):4171–4178, 2018.
  * [182] Zhiguang Liu and Jianhong Hao.  Intention recognition in physical human-robot interaction based on radial basis function neural network.  Journal of Robotics, 2019, 2019.
  * [183] Guangzhu Peng, Chenguang Yang, Wei He, and C. L. Philip Chen.  Force sensorless admittance control with neural learning for robots with actuator saturation.  IEEE Transactions on Industrial Electronics, 67(4):3138–3148, 2020.
  * [184] Lorenzo Vianello, Luigi Penco, Waldez Gomes, Yang You, Salvatore Maria Anzalone, Pauline Maurice, Vincent Thomas, and Serena Ivaldi.  Human-humanoid interaction and cooperation: a review.  Current Robotics Reports, 2(4):441–454, 2021.
  * [185] Martin Lawitzky, José Ramón Medina, Dongheui Lee, and Sandra Hirche.  Feedback motion planning and learning from demonstration in physical robotic assistance: differences and synergies.  In 2012 IEEE/RSJ International Conference on Intelligent Robots and Systems, pages 3646–3652, 2012.
  * [186] Shiqi Li, Ke Han, Xiao Li, Shuai Zhang, Youjun Xiong, and Zheng Xie.  Hybrid trajectory replanning-based dynamic obstacle avoidance for physical human-robot interaction.  Journal of Intelligent & Robotic Systems, 103:1–14, 2021.
  * [187] Marco Faroni, Manuel Beschi, and Nicola Pedrocchi.  An mpc framework for online motion planning in human-robot collaborative tasks.  In 2019 24th IEEE International Conference on Emerging Technologies and Factory Automation (ETFA), pages 1555–1558, 2019.
  * [188] Ajung Moon, Maneezhay Hashmi, HF Machiel Van Der Loos, Elizabeth A Croft, and Aude Billard.  Design of hesitation gestures for nonverbal human-robot negotiation of conflicts.  ACM Transactions on Human-Robot Interaction (THRI), 10(3):1–25, 2021.
  * [189] Sonia Chernova and Andrea L Thomaz.  Robot learning from human teachers.  Synthesis lectures on artificial intelligence and machine learning, 8(3):1–121, 2014.
  * [190] Heni Ben Amor, Gerhard Neumann, Sanket Kamthe, Oliver Kroemer, and Jan Peters.  Interaction primitives for human-robot cooperation tasks.  In 2014 IEEE International Conference on Robotics and Automation (ICRA), pages 2831–2837, 2014.
  * [191] Yujun Lai, Gavin Paul, Yunduan Cui, and Takamitsu Matsubara.  User intent estimation during robot learning using physical human robot interaction primitives.  Autonomous Robots, 46(2):421–436, 2022.
  * [192] Scott Niekum, Sarah Osentoski, George Konidaris, and Andrew G Barto.  Learning and generalization of complex tasks from unstructured demonstrations.  In 2012 IEEE/RSJ International Conference on Intelligent Robots and Systems, pages 5239–5246. IEEE, 2012.
  * [193] Nadia Figueroa, Ana Lucia Pais Ureche, and Aude Billard.  Learning complex sequential tasks from demonstration: A pizza dough rolling case study.  In 2016 11th ACM/IEEE International Conference on Human-Robot Interaction (HRI), pages 611–612, 2016.
  * [194] Rui Huang, Hong Cheng, Jing Qiu, and Jianwei Zhang.  Learning physical human–robot interaction with coupled cooperative primitives for a lower exoskeleton.  IEEE Transactions on Automation Science and Engineering, 16(4):1566–1574, 2019.
  * [195] Dylan P Losey, Andrea Bajcsy, Marcia K O’Malley, and Anca D Dragan.  Physical interaction as communication: Learning robot objectives online from human corrections.  The International Journal of Robotics Research, 41(1):20–44, 2022.
  * [196] Dylan P Losey, Craig G McDonald, Edoardo Battaglia, and Marcia K O’Malley.  A review of intent detection, arbitration, and communication aspects of shared control for physical human–robot interaction.  Applied Mechanics Reviews, 70(1), 2018.
  * [197] Mario Selvaggio, Marco Cognetti, Stefanos Nikolaidis, Serena Ivaldi, and Bruno Siciliano.  Autonomy in physical human-robot interaction: A brief survey.  IEEE Robotics and Automation Letters, 6(4):7989–7996, 2021.
  * [198] Said G Khan, Guido Herrmann, Mubarak Al Grafi, Tony Pipe, and Chris Melhuish.  Compliance control and human–robot interaction: Part 1—survey.  International journal of humanoid robotics, 11(03):1430001, 2014.
  * [199] Hsieh-Yu Li, Audelia G Dharmawan, Ishara Paranawithana, Liangjing Yang, and U-Xuan Tan.  A control scheme for physical human-robot interaction coupled with an environment of unknown stiffness.  Journal of Intelligent & Robotic Systems, 100:165–182, 2020.
  * [200] Mojtaba Sharifi, Saeed Behzadipour, and Gholamreza Vossoughi.  Nonlinear model reference adaptive impedance control for human–robot interactions.  Control Engineering Practice, 32:9–27, 2014.
  * [201] Issac Rhee, Gitae Kang, Seung Jae Moon, Yun Seok Choi, and Hyouk Ryeol Choi.  Hybrid impedance and admittance control of robot manipulator with unknown environment.  Intelligent Service Robotics, 16(1):49–60, 2023.
  * [202] Kevin Haninger, Christian Hegeler, and Luka Peternel.  Model predictive control with gaussian processes for flexible multi-modal physical human robot interaction.  In 2022 International Conference on Robotics and Automation (ICRA), pages 6948–6955. IEEE, 2022.
  * [203] Bryan Whitsell and Panagiotis Artemiadis.  Physical human–robot interaction (phri) in 6 dof with asymmetric cooperation.  IEEE Access, 5:10834–10845, 2017.
  * [204] Brahim Brahmi, Mohamed Hamza Laraki, Maarouf Saad, Cristobal Ochoa-Luna, and Abdelkrim Brahmi.  Compliant adaptive control of human upper-limb exoskeleton robot with unknown dynamics based on a modified function approximation technique (mfat).  Robotics and Autonomous Systems, 117:92–102, 2019.
  * [205] Arvid QL Keemink, Herman van der Kooij, and Arno HA Stienen.  Admittance control for physical human–robot interaction.  The International Journal of Robotics Research, 37(11):1421–1444, 2018.
  * [206] H-P Huang and S-S Chen.  Compliant motion control of robots by using variable impedance.  The International Journal of Advanced Manufacturing Technology, 7(6):322–332, 1992.
  * [207] Mojtaba Sharifi, Amir Zakerimanesh, Javad K. Mehr, Ali Torabi, Vivian K. Mushahwar, and Mahdi Tavakoli.  Impedance variation and learning strategies in human–robot interaction.  IEEE Transactions on Cybernetics, 52(7):6462–6475, 2022.
  * [208] Kyeong Ha Lee, Seung Guk Baek, Hyuk Jin Lee, Seung Ho Lee, and Ja Choon Koo.  Real-time adaptive impedance compensator using simultaneous perturbation stochastic approximation for enhanced physical human–robot interaction transparency.  Robotics and Autonomous Systems, 147:103916, 2022.
  * [209] Wen Yu and Adolfo Perrusquía.  Simplified stable admittance control using end-effector orientations.  International Journal of Social Robotics, 12(5):1061–1073, 2020.
  * [210] Jianwei Dong, Jianming Xu, Qiaoqian Zhou, and Songda Hu.  Physical human–robot interaction force control method based on adaptive variable impedance.  Journal of the Franklin Institute, 357(12):7864–7878, 2020.
  * [211] Francesco Romano, Gabriele Nava, Morteza Azad, Jernej Čamernik, Stefano Dafarra, Oriane Dermy, Claudia Latella, Maria Lazzaroni, Ryan Lober, Marta Lorenzini, Daniele Pucci, Olivier Sigaud, Silvio Traversaro, Jan Babič, Serena Ivaldi, Michael Mistry, Vincent Padois, and Francesco Nori.  The codyco project achievements and beyond: Toward human aware whole-body controllers for physical human robot interaction.  IEEE Robotics and Automation Letters, 3(1):516–523, 2018.
  * [212] Kazuya Otani, Karim Bouyarmane, and Serena Ivaldi.  Generating assistive humanoid motions for co-manipulation tasks with a multi-robot quadratic program controller.  In 2018 IEEE International Conference on Robotics and Automation (ICRA), pages 3107–3113, 2018.
  * [213] Marie Charbonneau, Valerio Modugno, Francesco Nori, Giuseppe Oriolo, Daniele Pucci, and Serena Ivaldi.  Learning robust task priorities of qp-based whole-body torque-controllers.  In 2018 IEEE-RAS 18th International Conference on Humanoid Robots (Humanoids), pages 1–9, 2018.
  * [214] Francesco Tassi, Elena De Momi, and Arash Ajoudani.  An adaptive compliance hierarchical quadratic programming controller for ergonomic human–robot collaboration.  Robotics and Computer-Integrated Manufacturing, 78:102381, 2022.
  * [215] Francesco Tassi and Arash Ajoudani.  Multi-modal and adaptive control of human-robot interaction through hierarchical quadratic programming.  2023\.
  * [216] Enrico Mingo Hoffman, Brice Clement, Chengxu Zhou, Nikos G Tsagarakis, Jean-Baptiste Mouret, and Serena Ivaldi.  Whole-body compliant control of icub: first results with opensot.  In IEEE/RAS ICRA Workshop on Dynamic Legged Locomotion in Realistic Terrains, 2018.
  * [217] Diego Felipe Paez Granados, Jun Kinugawa, Yasuhisa Hirata, and Kazuhiro Kosuge.  Guiding human motions in physical human-robot interaction through com motion control of a dance teaching robot.  In 2016 IEEE-RAS 16th International Conference on Humanoid Robots (Humanoids), pages 279–285, 2016.
  * [218] Ganna Pugach, Artem Melnyk, Olga Tolochko, Alexandre Pitti, and Philippe Gaussier.  Touch-based admittance control of a robotic arm using neural learning of an artificial skin.  In 2016 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS), pages 3374–3380. IEEE, 2016.
  * [219] Sven Cremer, Sumit Kumar Das, Indika B. Wijayasinghe, Dan O. Popa, and Frank L. Lewis.  Model-free online neuroadaptive controller with intent estimation for physical human–robot interaction.  IEEE Transactions on Robotics, 36(1):240–253, 2020.
  * [220] Nguyen-Van Toan, Phan-Bui Khoi, and Soo-Yeong Yi.  A mlp-hedge-algebras admittance controller for physical human–robot interaction.  Applied Sciences, 11(12), 2021.
  * [221] Xinbo Yu, Bin Li, Wei He, Yanghe Feng, Long Cheng, and Carlos Silvestre.  Adaptive-constrained impedance control for human–robot co-transportation.  IEEE Transactions on Cybernetics, 52(12):13237–13249, 2022.
  * [222] Xinbo Yu, Sisi Liu, Shuang Zhang, Wei He, and Haifeng Huang.  Adaptive neural network force tracking control of flexible joint robot with an uncertain environment.  IEEE Transactions on Industrial Electronics, 71(6):5941–5949, 2024.
  * [223] Dong Wei, Lipeng Chen, Longfei Zhao, Hua Zhou, and Bidan Huang.  A vision-based measure of environmental effects on inferring human intention during human robot interaction.  IEEE Sensors Journal, 22(5):4246–4256, 2022.
  * [224] Loris Roveda, Jeyhoon Maskani, Paolo Franceschi, Arash Abdi, Francesco Braghin, Lorenzo Molinari Tosatti, and Nicola Pedrocchi.  Model-based reinforcement learning variable impedance control for human-robot collaboration.  Journal of Intelligent & Robotic Systems, 100(2):417–433, 2020.
  * [225] Ruofan Wu, Minhan Li, Zhikai Yao, Wentao Liu, Jennie Si, and He Huang.  Reinforcement learning impedance control of a robotic prosthesis to coordinate with human intact knee motion.  IEEE Robotics and Automation Letters, 7(3):7014–7020, 2022.
  * [226] Fotios Dimeas and Nikos Aspragathos.  Reinforcement learning of variable admittance control for human-robot co-manipulation.  In 2015 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS), pages 1011–1016, 2015.
  * [227] Timothy P Lillicrap, Jonathan J Hunt, Alexander Pritzel, Nicolas Heess, Tom Erez, Yuval Tassa, David Silver, and Daan Wierstra.  Continuous control with deep reinforcement learning.  arXiv preprint arXiv:1509.02971, 2015.
  * [228] Quan Liu, Zhihao Liu, Bo Xiong, Wenjun Xu, and Yang Liu.  Deep reinforcement learning-based safe interaction for industrial human-robot collaboration using intrinsic reward function.  Advanced Engineering Informatics, 49:101360, 2021.
  * [229] Jong In Han, Jeong-Hoon Lee, Ho Seon Choi, Jung-Hoon Kim, and Jongeun Choi.  Policy design for an ankle-foot orthosis using simulated physical human–robot interaction via deep reinforcement learning.  IEEE Transactions on Neural Systems and Rehabilitation Engineering, 30:2186–2197, 2022.
  * [230] Mojtaba Sharifi, Vahid Azimi, Vivian K. Mushahwar, and Mahdi Tavakoli.  Impedance learning-based adaptive control for human–robot interaction.  IEEE Transactions on Control Systems Technology, 30(4):1345–1358, 2022.
  * [231] Klas Kronander and Aude Billard.  Online learning of varying stiffness through physical human-robot interaction.  In 2012 IEEE International Conference on Robotics and Automation, pages 1842–1849, 2012.
  * [232] Klas Kronander and Aude Billard.  Learning compliant manipulation through kinesthetic and tactile human-robot interaction.  IEEE transactions on haptics, 7(3):367–380, 2013.
  * [233] Wentao Liu, Ruofan Wu, Jennie Si, and He Huang.  A new robotic knee impedance control parameter optimization method facilitated by inverse reinforcement learning.  IEEE Robotics and Automation Letters, 7(4):10882–10889, 2022.
  * [234] Neil C Thompson, Kristjan Greenewald, Keeheon Lee, and Gabriel F Manso.  The computational limits of deep learning.  arXiv preprint arXiv:2007.05558, 2020.
  * [235] Dmitry Baskakov and Dmitry Arseniev.  On the computational complexity of deep learning algorithms.  In Proceedings of International Scientific Conference on Telecommunications, Computing and Control: TELECCON 2019, pages 343–356. Springer, 2021.
  * [236] William J Dally, Stephen W Keckler, and David B Kirk.  Evolution of the graphics processing unit (gpu).  IEEE Micro, 41(6):42–51, 2021.
  * [237] Rubens Luiz Rech and Paolo Rech.  Reliability of google’s tensor processing units for embedded applications.  In 2022 Design, Automation & Test in Europe Conference & Exhibition (DATE), pages 376–381, 2022.
  * [238] P.H.W. Leong and K.H. Tsoi.  Field programmable gate array technology for robotics applications.  In 2005 IEEE International Conference on Robotics and Biomimetics - ROBIO, pages 295–298, 2005.
  * [239] Brian Plancher, Sabrina M. Neuman, Thomas Bourgeat, Scott Kuindersma, Srinivas Devadas, and Vijay Janapa Reddi.  Accelerating robot dynamics gradients on a cpu, gpu, and fpga.  IEEE Robotics and Automation Letters, 6(2):2335–2342, 2021.
  * [240] Yuexuan Tu, Saad Sadiq, Yudong Tao, Mei-Ling Shyu, and Shu-Ching Chen.  A power efficient neural network implementation on heterogeneous fpga and gpu devices.  In 2019 IEEE 20th international conference on information reuse and integration for data science (IRI), pages 193–199. IEEE, 2019.
  * [241] David Kortenkamp, Reid Simmons, and Davide Brugali.  Robotic systems architectures and programming.  Springer handbook of robotics, pages 283–306, 2016.
  * [242] Félix Ingrand and Malik Ghallab.  Deliberation for autonomous robots: A survey.  Artificial Intelligence, 247:10–44, 2017.
  * [243] A. Meystel.  Planning in a hierarchical nested controller for autonomous robots.  In 1986 25th IEEE Conference on Decision and Control, pages 1237–1249, 1986.
  * [244] James S Albus et al.  The nist real-time control system (rcs): An application survey.  In Proc. of the AAAI 1995 Spring Symposium Series, Stanford University, Menlo Park, CA, 1995.
  * [245] Ronald C Arkin.  Behavior-based robotics.  MIT press, 1998.
  * [246] Faisal Qureshi, Demetri Terzopoulos, and Ross Gillett.  The cognitive controller: a hybrid, deliberative/reactive control architecture for autonomous robots.  In International Conference on Industrial, Engineering and Other Applications of Applied Intelligent Systems, pages 1102–1111. Springer, 2004.
  * [247] R. Brooks.  A robust layered control system for a mobile robot.  IEEE Journal on Robotics and Automation, 2(1):14–23, 1986.
  * [248] Ian Horswill.  Polly: A vision-based artificial agent.
  * [249] Maja J Matarić.  Integration of representation into goal-driven behavior-based robots.  In The artificial life route to artificial intelligence, pages 165–186. Routledge, 2018.
  * [250] Daniel Toal, Colin Flanagan, Caimin Jones, and Bob Strunz.  Subsumption architecture for the control of robots.  IMC-13, Limerick, 1996.
  * [251] Michele Amoretti and Monica Reggiani.  Architectural paradigms for robotics applications.  Advanced Engineering Informatics, 24(1):4–13, 2010.
  * [252] PU Chavan, M Murugan, and PP Chavan.  A review on software architecture styles with layered robotic software architecture.  In 2015 International Conference on Computing Communication Control and Automation, pages 827–831. IEEE, 2015.
  * [253] Morgan Quigley, Ken Conley, Brian Gerkey, Josh Faust, Tully Foote, Jeremy Leibs, Rob Wheeler, Andrew Y Ng, et al.  Ros: an open-source robot operating system.  In ICRA workshop on open source software, volume 3, page 5. Kobe, Japan, 2009.
  * [254] Alberto Martinetti, Peter K Chemweno, Kostas Nizamis, and Eduard Fosch-Villaronga.  Redefining safety in light of human-robot interaction: A critical review of current standards and regulations.  Frontiers in chemical engineering, 3:32, 2021.
  * [255] Laurel Riek and Don Howard.  A code of ethics for the human-robot interaction profession.  Proceedings of we robot, 2014.
  * [256] Yusuf Aydin, Ozan Tokatli, Volkan Patoglu, and Cagatay Basdogan.  A computational multicriteria optimization approach to controller design for physical human-robot interaction.  IEEE Transactions on Robotics, 36(6):1791–1804, 2020.
  * [257] Rossitza Setchi, Maryam Banitalebi Dehkordi, and Juwairiya Siraj Khan.  Explainable robotics in human-robot interactions.  Procedia Computer Science, 176:3057–3066, 2020.
  * [258] Francisco J Rodríguez Lera, Camino Fernández Llamas, Ángel Manuel Guerrero, and Vicente Matellán Olivera.  Cybersecurity of robotics and autonomous systems: Privacy and safety.  Robotics-legal, ethical and socioeconomic impacts, 2017.
  * [259] Adriana Hamacher, Nadia Bianchi-Berthouze, Anthony G. Pipe, and Kerstin Eder.  Believing in bert: Using expressive communication to enhance trust and counteract operational error in physical human-robot interaction.  In 2016 25th IEEE International Symposium on Robot and Human Interactive Communication (RO-MAN), pages 493–500, 2016.
  * [260] Yue Hu, Naoko Abe, Mehdi Benallegue, Natsuki Yamanobe, Gentiane Venture, and Eiichi Yoshida.  Toward active physical human–robot interaction: Quantifying the human state during interactions.  IEEE Transactions on Human-Machine Systems, 52(3):367–378, 2022.
  * [261] Nora Fronemann, Kathrin Pollmann, and Wulf Loh.  Should my robot know what’s best for me? human–robot interaction between user experience and ethical design.  AI & SOCIETY, 37(2):517–533, 2022.
  * [262] Ayanna Howard and Jason Borenstein.  The ugly truth about ourselves and our robot creations: the problem of bias and social inequity.  Science and engineering ethics, 24:1521–1536, 2018.
  * [263] Tom Williams.  The eye of the robot beholder: Ethical risks of representation, recognition, and reasoning over identity characteristics in human-robot interaction.  In Companion of the 2023 ACM/IEEE International Conference on Human-Robot Interaction, HRI ’23, page 1–10, New York, NY, USA, 2023. Association for Computing Machinery.
  * [264] Antonia Meissner, Angelika Trübswetter, Antonia S Conti-Kufner, and Jonas Schmidtler.  Friend or foe? understanding assembly workers’ acceptance of human-robot collaboration.  ACM Transactions on Human-Robot Interaction (THRI), 10(1):1–30, 2020.
  * [265] Aimee van Wynsberghe, Madelaine Ley, and Sabine Roeser.  Ethical aspects of human–robot collaboration in industrial work settings.  The 21st Century Industrial Robot: When Tools Become Collaborators, pages 255–266, 2022.
  * [266] Reza Etemad-Sajadi, Antonin Soussan, and Théo Schöpfer.  How ethical issues raised by human–robot interaction can impact the intention to use the robot?  International journal of social robotics, 14(4):1103–1115, 2022.
  * [267] AJung Moon, Shalaleh Rismani, and HF Machiel Van der Loos.  Ethics of corporeal, co-present robots as agents of influence: a review.  Current Robotics Reports, 2:223–229, 2021.
  * [268] Nora Fronemann, Kathrin Pollmann, and Wulf Loh.  Should my robot know what’s best for me? human–robot interaction between user experience and ethical design.  AI & SOCIETY, 37(2):517–533, 2022.


# Quantum-assisted distortion-free audio Signal Sensing

Chen Zhang<EMAIL_ADDRESS>3rd Institute of Physics, University
of Stuttgart, Pfaffenwaldring 57, Stuttgart 70569, Germany  Durga Dasari
<EMAIL_ADDRESS>3rd Institute of Physics, University of
Stuttgart, Pfaffenwaldring 57, Stuttgart 70569, Germany  Matthias Widmann  3rd
Institute of Physics, University of Stuttgart, Pfaffenwaldring 57, Stuttgart
70569, Germany  Jonas Meinel  3rd Institute of Physics, University of
Stuttgart, Pfaffenwaldring 57, Stuttgart 70569, Germany  Vadim Vorobyov  3rd
Institute of Physics, University of Stuttgart, Pfaffenwaldring 57, Stuttgart
70569, Germany  Polina Kapitanova  Department of Physics and Engineering, ITMO
University, Saint Petersburg 197101, Russia  Elizaveta Nenasheva  Giricond
Research Institute, Ceramics Co. Ltd., Saint Petersburg 194223, Russia  Kazuo
Nakamura  Hydrogen and Carbon Management Technology Section, Hydrogen and
Carbon Management Technology Strategy Department, Tokyo Gas Co., Ltd.,
Yokohama 230-0045, Japan  Hitoshi Sumiya  Advanced Materials Labotratory,
Sumitomo Electric Industries, Ltd., Itami 664-0016, Japan  Shinobu Onoda
Takasaki Advanced Radiation Research Institute, National Institutes for
Quantum Science and Technology, Takasaki 370-1292, Japan  Junichi Isoya
Faculty of Pure and Applied Sciences, University of Tsukuba, Tsukuba 305-8573,
Japan  Jörg Wrachtrup  3rd Institute of Physics, University of Stuttgart,
Pfaffenwaldring 57, Stuttgart 70569, Germany

###### Abstract

Quantum sensors are keeping the cutting-edge sensitivities in metrology.
However, for high-sensitive measurements of arbitrary signals, limitations in
linear dynamic range could introduce distortions when sensing the frequency,
magnitude and phase of unknown signals. Here, we overcome these limitations
with advanced sensing protocol that combines quantum phase-sensitive detection
with the heterodyne readout. We present theoretical and experimental
investigations using nitrogen-vacancy centers in diamond, showing the ability
to sense radio signals with a 98 dB linear dynamic range, a 31 pT/Hz1/2
sensitivity, and arbitrary frequency resolution. Further, we perform the
quantum-assisted distortion-free audio signal (melody, speech) sensing with
high fidelity. The methods developed here could broaden the horizon for
quantum sensors towards applications in telecommunication, where high-fidelity
and low-distortion at multiple frequency bands within small sensing volumes
are required.

##  I Introduction

Quantum sensors are setting new frontiers of sensing techniques with their
extraordinary performances in sensitivity and stability [1, 2, 3, 4, 5]. These
techniques rely on either measuring the line-shift of spin or atomic
transition frequencies or reading out the relative populations of the occupied
energy levels using interferometric methods [6, 7]. In most cases, there are
trade-off relations between the sensitivity and other features in metrology
[8]. For example, a high-sensitive measurement acquired by detecting the
transition line shift requires a narrow linewidth, which, on the other hand,
will limit the dynamic range. Interferometric measurements detect a sinusoidal
response, and linearity is only achieved when the phase signal is in a small
dynamic range. It sets a massive limitation on the sensitivity when sensing an
unknown signal that gets measured beyond this linear regime, for example, when
the working point of the sensor is at the maxima or minima of the sinusoidal
signal response. Thus, it becomes a bottleneck for high sensitivity
measurements that are required in many cutting-edge applications. Operating
within the linear dynamic range (LDR) can be crucial for reconstructing
unknown signals. One way to directly extract the phase factor, which is linear
to the physical quantity to be detected, is to use phase-sensitive detection
known as the classical lock-in technique. In this work, using a modified
sensing scheme that introduces an external readout phase modulation, we
acquire the target quantum phase signal after demodulation. Therefore, the LDR
is no longer limited to the small-angle approximation. Hence our method
combines large dynamic range with maximum sensitivity.

Nitrogen-vacancy (NV) centers in diamond have been at the forefront in
performing high-sensitive measurements of various physical quantities, viz.,
magnetic and electric field, temperature, and strain distributions internal
and external to diamond [9, 10, 11, 12, 13]. The NV magnetometry has been
performed under bias fields ranging from zero-field to a few Tesla, and for
sensing signals with frequencies ranging from DC to a few GHz [14, 15, 16,
17]. While dynamical-decoupling techniques are usually employed for high
sensitivity [9, 18, 19], arbitrary frequency resolution can be achieved with
the quantum heterodyne (Q-dyne) detection [20, 21]. However, both methods
suffer from a limited LDR when they are applied to measure arbitrary signals.

For high dynamic range measurements, a closed-loop frequency-locking scheme
together with optically detected magnetic resonance (ODMR) can be used to
track resonance frequency shifts [22]. However, this scheme cannot be used for
ac field measurements in combination with interferometric methods, if the
signal frequency is higher than the readout sampling frequency. Phase-
estimation algorithms can effectively improve the LDR in Ramsey measurements
by varying the sequence with adaptive feedback schemes [23, 24]. However, for
the case of ac sensing schemes e.g. Hahn-echo, as varying the sequence itself
will change the sensor response to the ac signals, such methods become less
applicable. Therefore, a technique is still missing, that addresses the LDR
while maintaining high sensitivity and frequency resolution, for example, in
sensing arbitrary radio-frequency fields within a broad bandwidth.

Sensing radio-frequency signals by electric-field sensors, either conventional
electronic receivers or the Rydberg atom sensors, need antennas to collect and
guide the electric signals towards the sensors [25, 26, 27, 28]. Although the
receivers can be highly integrated, the dimension of antennas can scale to
meters due to the signal wavelength. The size becomes critical when there is
limited space for the sensor, for example, in a satellite. In this regard,
quantum magnetometers can be very attractive due to their small sensing volume
and high sensitivity [29]. A flux concentrator can be used as a substitute to
conventional antennas for obtaining high signal gain. Independent of the
signal wavelength, the dimensions of such flux concentrators can be as small
as a few centimeters [30, 31].

In this paper, we demonstrate the quantum-assisted distortion-free audio
signal sensing with NV center ensembles in diamond using a quantum-phase-
sensitive detection (QPSD) technique combined with heterodyne readout.
Firstly, we introduce the QPSD technique, which can provide an extended LDR in
interferometry measurements by using two synchronized driving fields. Then, we
present the heterodyne readout, which can interpret e.g. radio signals to get
frequency information. Taking advantage of the bandwidth of the Hahn-echo
sequence and the frequency comb induced by the continuous sampling, we
demonstrate measurements of audio signals around 10 kHz, beyond the coherence
limit without losing sensitivity. Finally, we present arbitrary radio signal
measurements with a LDR of 98 dB at a sensitivity of 31 pT/Hz1/2. A piece of
melody and a speech are encoded as magnetic field signals and measured by the
NV magnetometer. By using the sensor as a quantum radio receiver, we
demonstrate the application potentials for areas such as quantum-assisted
telecommunication and unknown signal exploration.

##  II Results

###  II.1 Quantum Phase Sensitive Detection

Figure 1:  Phase-sensitive NV magnetometry (A) Continuous sampling induced
phase reviving signals, known as the quantum heterodyne (Q-dyne) detection.
The phase reviving frequency changes with the external field and cannot be
used for phase-sensitive detection. The signal responses $G(s_{\phi},k)$
should be small to ensure measurement linearity. $R$ is the detected photon
rate, and $C$ is the detected signal contrast. (B) Unlike the Q-dyne
detection, the quantum phase-sensitive detection (QPSD) is based on the
rotating frame modulation induced by the evolving phase difference of the two
driving MW fields. Two frequency-offset MWs acquire a phase difference of
$\alpha-\beta=2\pi\delta f\Delta t$ after the sampling time interval $\Delta
t$. In the Bloch sphere picture, it can be understand as the MW2 defined
rotating frame $x_{2}y_{2}z$ rotates with rate of $\delta f$ referring to the
MW1 defined rotating frame $x_{1}y_{1}z$. The acquired quantum phase is
$\theta=2\pi\delta f\Delta t-\phi$ at sample 2. Through the quantum phase
modulation, the acquired readout representing the Bloch vector projections is
as shown in (C), where we present the measurements of the quantum phase
$\phi=0$ and $\phi\neq 0$. By demodulating the acquired signal with a lock-in
amplifier, we can get the phase values. The dashed box shows the measurement
sequences we applied in experiments. Except for the last $\pi/2$ pulse applied
with MW2, all the other driving pulses are generated by MW1. The fluorescence
signal is demodulated at the frequency of $1/(2T_{seq})$ to get a fluorescence
intensity readout for a sample. The QPSD readout is acquired with the
demodulation at $\delta f$. (D) Schematic of the experiment. NV centers
ensemble in diamond is used to perform the QPSD readout.

In interferometric measurements, the quantum phase to be detected is usually
converted to a quantum state population difference, resulting in a sinusoid
readout [9]. A way to extract the phase factor from the sinusoidal readout is
to modulate the phase with a specific frequency and perform phase-sensitive
modulation. Such a quantum phase modulation can be introduced by using the
difference between the quantum phase of the sensor to an external oscillator.
The Q-dyne method uses such a strategy for resolving frequency of signals
better than the relaxation time of the sensor, as shown in Fig. 1a [20, 21].
However, it cannot be used for phase-sensitive detection because the Q-dyne
frequency is also what to be resolved and an extra modulation is still needed
[21]. Another way to introduce such a phase modulation is to use the frequency
offset between two different driving fields of the sensor [32, 33]. These
driving fields define two rotating frames, and the evolution of the spin as
seen from one rotating frame will lead to a quantum phase modulation due to
the relative rotation of the two frames, i.e., rotating frame modulation. The
modulation frequency only depends on the frequency difference of the two
driving fields, as shown in Fig. 1b and c. By performing multiple measurements
within a modulation cycle and upon using lock-in detection, we will achieve
phase-sensitive detection. Below we mathematically describe this relative
evolution of the sensor under such interferometric measurement with two-
frequency driving fields.

Aligning an external field $B_{0}$ along the NV axis, we use the two-level
subspace of the NV ground triplet in the derivation. Thus, the Hamiltonian of
the system can be written as:

$\mathcal{H}=\omega_{0}S_{z}+\gamma_{e}B_{1}\cos\left(2\pi
f_{1}t+\alpha\right)S_{x},$ (1)

where $\omega_{0}$ is the transition frequency of the two-level subspace,
$B_{1}$ is the oscillating magnetic field perpendicular to the NV axis,
$f_{1}$ and $\alpha$ are the frequency and phase of the driving field, and
$\gamma_{e}$ is the electron gyromagnetic ratio. In the rotating frame defined
by the resonance frequency, the time-dependent Hamiltonian is

$\mathcal{H}_{1}^{\prime}=\Omega_{1}\cos\left(\delta\omega_{1}+\alpha\right)S_{x}+\Omega_{1}\sin\left(\delta\omega_{1}t+\alpha\right)S_{y},$
(2)

where $\delta\omega_{1}=2\pi(f_{0}-f_{1})$, and $\Omega_{1}=\gamma_{e}B_{1}/2$
is the Rabi frequency introduced by MW1. In interferometry measurements, a
$\pi/2$ pulse prepares the spin state from the polarized state to an equalized
population, and another $\pi/2$ pulse projects the quantum phase as a
population difference after the sensing procedure. We use the second driving
field, MW2, to offset the frequency of the second $\pi/2$ pulse.
$\delta\omega_{2}$,$\beta$ and $\Omega_{2}$ are used to denote the frequency
offset, Rabi frequency, and phase of MW2. After this, the measured spin-
expectation value is

$\braket{S_{z}}=\sin\left[\phi+\frac{\pi}{2}\left(\frac{\delta\omega_{1}}{\Omega_{1}}-\frac{\delta\omega_{2}}{\Omega_{2}}\right)+\alpha-\beta\right],$
(3)

where $\phi$ is the acquired quantum phase which contains the information we
want to measure, both of the MWs are near-resonant with
$\delta\omega_{1}\ll\Omega_{1}$, $\delta\omega_{2}\ll\Omega_{2}$. Therefore,
the off-resonant term can be neglected, and the phase difference term
$\alpha-\beta$ will evolve with time so that there is

$\braket{S_{z}}\approx\sin\left(\phi+2\pi\delta f\cdot t\right),$ (4)

where $\delta f$ is the frequency difference of the two MWs.

The above result can be seen as a modulation of the rotating frame itself. As
schematically shown in Fig. 1b (left Bloch sphere), we can assume that the two
driving fields have the same phase at the duration of the second $\pi/2$
pulse, and this defines an instantaneous rotating frame with coordinates
$x_{1}y_{1}z$. Thus, the readout is similar to that of the regular Ramsey
interferometry using a single driving field. After an interval of $\Delta t$,
MW2 develops a phase difference of $2\pi\delta f\Delta t$. Since the quantum
phase is finally measured by MW2, the Bloch vector rotates in the new
instantaneous rotating frame with coordinates $x_{2}y_{2}z$, as shown in Fig.
1b (the right Bloch sphere). The rotating frame defined by MW2 rotates
continuously around the $z$-axis with the frequency of $\delta f$. Due to
this, the fluorescence readout also modulates in a sinusoidal fashion, as
shown in Fig. 1c. While the readout signal frequency depends on $\delta f$ and
its amplitude depends on the signal contrast, the initial phase, $\varphi$, is
linear to the field to be measured. Through the external modulation induced by
the MWs, the working point of the sensor evolves in the entire phase range,
which ensures the LDR of the initial phase measurement covering $[-\pi,\pi]$.
By fitting or demodulating the fluorescence signal, we can resolve the
changing of the phase factor $\varphi$ between each modulation cycle and find
measurement linearity for the external field.

The measurement sequence we applied in the experiment is depicted in Fig. 1c,
in which $T_{\phi}$ is the field sensing time, $T_{seq}$ is the sequence
length of one measurement, and we use a second measurement with the final
pulse changed to $(\pi/2)_{-x}$. As a result, the fluorescence signal is
modulated with a frequency of $f_{s}=1/(2T_{seq})$, which is also the sampling
frequency of the fluorescence readout. The demodulation of the fluorescence
signal, denoted by Demod. 1 in Fig. 1d, has a readout bandwidth $f_{s}/2$ set
by the Shannon sampling theorem. The readout is further demodulated by another
demodulator of the lock-in amplifier (LIA), denoted as Demod. 2. Upon
measuring $N$ samples of the fluorescence readout, the bandwidth of the phase
readout is narrowed down to $f_{s}/(2N)$. These measurements are schematically
shown in Fig. 1d.

The sensitivity of such measurements can be derived based on the fitting of
the $N$ samples in the measurement time of $N\cdot 2T_{seq}$. The minimum
detectable phase is derived as

$\delta\phi=\frac{2}{\sqrt{N}}\frac{1}{C\sqrt{\mathcal{N}}},$ (5)

where $C$ is the fluorescence signal contrast, $\mathcal{N}$ is the detected
photon counts in each measurement. The sensitivity to external magnetic field,
however, is still subject to the applied MW sequence, can be derived as

$\eta=\frac{2}{\gamma_{e}|G(\omega)|C}\sqrt{\frac{2T_{seq}}{{\mathcal{N}}}},$
(6)

where $|G(\omega)|$ is the MW filter function which is usually used to
describe the transfer function of such a sensor from magnetic field to quantum
phase. In comparison to the conventional fluorescence readout, the sensitivity
of QPSD readout deteriorates by a factor of $\sqrt{2}$. Details about the
sensitivity derivation can be seen in the Supplementary Materials.

Figure 2:  Sensing performance of the QPSD. (A) Spectra and linearity
comparison of the normal Ramsey readout and the QPSD readout. We apply
$T_{\phi}=6.25\mu s$ in both measurements. The applied peak-to-peak field and
the readout are plotted showing the linearity of the measurements. (B) Spectra
and linearity comparison of normal Hahn-echo readout and the QPSD readout. The
applied phase accumulation time $T_{\phi}=12.5\mu s$. Thus, the Hahn-echo
measurements performs a higher sensitivity but smaller dynamic range than the
Ramsey measurements. (C) Signal response to different sampling frequencies.
The measurements use the same calibration field, and the readouts are
normalized to be plotted in the same vertical axis. (D) Measurement bandwidth.
Ramsey and Hahn-echo sequences are applied to measure test fields at different
frequencies with the same magnitude. The heterodyne frequency responses are
limited in bandwidth by the cut-off frequency of the LIA.

In Fig. 2a and b, we compare the regular interferometry (single driving field)
and with the measurements obtained from the QPSD readout described above. The
strength of the applied external ac fields ranges from 0 to 3
$\mathop{}\\!\mathrm{\mu T}$. For Ramsey measurements, the applied fields are
at a frequency of 46 Hz, and we use a field sensing time
$T_{\phi,Ramsey}=6.25\mathop{}\\!\mathrm{\mu s}$. For Hahn-echo measurements,
we use external fields at 80 kHz+46 Hz and the field sensing time
$T_{\phi,Hahn}=2T_{\phi,Ramsey}$. The test fields are sent to the diamond by a
calibrated loop antenna. The signal readout of the regular interferometry
measurements is proportional to $\sin(\phi)$, where $\phi\propto\gamma_{e}B$
is the accumulated quantum phase. The response is linear only when $\phi$ is
small, limiting the dynamic range. Thus, the regular Ramsey and Hahn-echo
readout quickly saturate due to this limited LDR. We plot the Fourier
transform of the readout signals also in the figures. The harmonics of the 46
Hz signal rise significantly in the fluorescence readout spectral due to the
saturation induced by the limited LDR, compared to the QPSD readout which
shows the linearity over the entire field range. The high-order harmonics of
the signal detected by the QPSD readout are small and mainly arise from the
function generator. In the measurements, one could see the linewidth
broadening induced by the increasing signal power. The peak at 100 Hz, which
is consistently seen in both the Ramsey and Hahn-echo measurements, comes from
the electronics instrumentation. Other side peaks seen near the original
signal frequency in the QPSD readout spectra are due to the mixing of the 100
Hz power line harmonics and the 92 Hz signal harmonics in the LIA.

Besides the LDR, the method also demonstrates measurement robustness to
changing of $T_{seq}$. The motivation of using different $T_{seq}$ is to get
different sampling frequencies as well as measurement bandwidth. Signal
responses to different sampling frequencies, i.e., different $1/2T_{seq}$, are
plotted in Fig. 2c. Characterized by the same test field, fluorescence readout
shows varying signal responses over the sampling frequency range, while the
QPSD readout almost stays at the same level because the measured phase factor
only changes with the external field and the sensing time $T_{\phi}$. It also
indicates that the QPSD readout does not change with varying of fluorescence
signal contrast, which is affected by the low spin polarization rate when
$T_{seq}$ is short in the regime of low excitation laser power.

The measurement bandwidth of the QPSD readout is shown in Fig. 2d, where the
signal responses to different test field frequencies are plotted. The plotted
values are the magnitudes at the corresponding frequencies in the Fourier
transform of the QPSD readout. For the measurements based on the Hahn-echo
sequence, we detected the heterodyne signal for the ac fields. The applied
sequence length, $T_{seq}=100\mathop{}\\!\mathrm{\mu s}$, gives the
referencing frequency $f_{s}=$5 kHz for Demod. 1. We apply the second driving
field with the frequency offset at $\delta f=500$ Hz to have $N=10$ samples in
a modulation cycle. Due to this, there will be flexibility in deciding the
single measurement bandwidth by setting the time constant of Demod. 2. We
choose different settings corresponding to the cut-off frequency of the filter
at 100 Hz and 200 Hz. Finally, one can conclude that the rotating frame
modulation provides QPSD readout magnetometry that has enhanced LDR and
robustness in a flexible bandwidth. As we show below, this can be used for
measurement of arbitrary fields with low distortion.

###  II.2 Frequency Offset Heterodyne Readout

Heterodyne readout has been used to improve the frequency resolution
remarkably in nuclear magnetic resonance spectroscopy. It is also a way to
achieve high precision microwave sensing [34, 35, 36]. High-order dynamical
decoupling sequences are used to narrow the spectral linewidth by decoupling
the sensor response from unwanted signal frequencies [20, 21]. Here arises a
trade-off between the measurable signal bandwidth and fidelity. High-order
dynamical decoupling can ensure a high sensitivity but only allows to measure
signals within the narrow bandwidth defined by the sequence. On the other
hand, the lower limit on the detectable signal frequency is set by the
decoherence time of the sensor. Here, we will use the Hahn-echo sequence in
combination with the QPSD readout to measure signals at frequencies that go
beyond the coherence time of the sensor.

In Qdyne, the sampling time usually satisfies $T_{seq}\neq mT_{\phi}$ so as to
get the heterodyne signal [21]. The frequency of this heterodyne signal
depends on the timing offset. Here, we choose the measurement sampling time
$T_{seq}=mT_{\phi}$ to obtain the heterodyne readout depending on the signal
frequency offset from $1/T_{\phi}$. As a result, the detected phase of signals
at frequencies of $n/T_{\phi}$ is locked by the sequence, where $n$ can be a
random integer. On the other hand, the frequency offset of signals can also
introduce phase revivals, i.e. frequency offset heterodyne signal, as shown in
Fig. 3a. The detected heterodyne frequency would be the exact offset of the
signal frequency to $1/T_{\phi}$.

Figure 3:  Frequency Offset Heterodyne readout. (a) Hahn-echo sequence is used
for this demonstration. The detected phase of the ac signal is locked by the
sequence when the frequency is $f_{ac}=f_{\phi}$. Otherwise, a heterodyne
signal of $\left|f_{ac}-f_{\phi}\right|$ can be measured. The colored regions
mark where the quantum phase is accumulated, while phase accumulations at the
other areas are canceled in the spin evolution. The figure shows identical
heterodyne signals due to $f_{\phi}-f_{1}=f_{2}-f_{\phi}$. (b) We apply ac
fields at different frequencies with an offset of 5 Hz to the sensor so that 5
Hz peaks can be detected as the signal response. The signal frequency response
of the Hahn-echo sequence and CPMG-2 sequence are plotted after normalization,
respectively. In both measurements $T_{\phi}=50\mathop{}\\!\mathrm{\mu s}$,
and the sampling frequency is 2 kHz according to the applied sequence length.
The red lines are the filter functions in theory. (c) Signal frequency
response of Hahn-echo measurements with $1/T_{seq}=10$ kHz. The dash line
indicates the filter introduced by the lock-in amplifier. (d) Sequence
dependency of the frequency resolution. In meas. 1, a 20.005 kHz filed is
applied and measured by sequences with $T_{\phi}=50\mathop{}\\!\mathrm{\mu
s}\pm 4\mathop{}\\!\mathrm{ns}$, and $T_{seq}=20T_{\phi}$. In meas. 2, we keep
$T_{\phi}$ unchanged, and offset $T_{seq}$ with $\pm$4 ns. In meas. 3, the
frequency of the applied field is changed to 16.005 kHz while the other
parameters are the same as meas. 2.

The frequency offset heterodyne readout is modeled based on the MW sequence
filter [37, 38]. Sampling happens in each time interval of
$\left[NmT_{\phi},(Nm+1)T_{\phi}\right]$, where $N\in\mathbb{Z}$. For a random
ac signal component $B_{ac}(t)=B(\omega)e^{-i[\omega t+\varphi(\omega)]}$ and
a measurement with the MW $\pi$-pulse number of $n$, we can get the
accumulated quantum phase as (derivation see in the Supplementary Materials)

$\phi_{r}(N)=\left|G_{n}(\omega)\right|e^{i\left(-\frac{\omega
T_{\phi}}{2}-\frac{P}{2}\pi\right)}\gamma_{e}B(\omega)e^{-i\varphi(\omega)}e^{-i\omega
NmT_{\phi}},$ (7)

where $N$ denotes the sampling timestamp,
$G_{n}(\omega)=\left|G_{n}(\omega)\right|e^{i\left(-\frac{\omega
T_{\phi}}{2}-\frac{P}{2}\pi\right)}$ is the MW filter function, $P=1$ when the
$\pi$-pulse number $n$ is odd and $P=2$ when $n$ is even. Particularly, when
$n=1$ i.e. Hahn-echo sequence is applied, the filter function satisfies
$|G_{1}(2\pi/T_{\phi})|=|G_{1}(\pi/T_{\phi})|$. In principle, measurements of
signals at a wide frequency range is feasible by choosing the appropriate
$T_{\phi}$ in Hahn-echo measurements. For example, by using
$T_{\phi}<1\mathop{}\\!\mathrm{\mu s}$, one can achieve detection of signals
at frequencies higher than 1 MHz. It is more challenging to measure a signal
at a lower frequency, such as a signal at 10 kHz, for the reason that a longer
$T_{2}$ is required. With the property described above, it is feasible to use
$T_{\phi}=50\mathop{}\\!\mathrm{\mu s}$ rather than
$T_{\phi}=100\mathop{}\\!\mathrm{\mu s}$ to achieve the measurement with a
better sensitivity due to the higher signal contrast when $T_{\phi}$ is
smaller. For diamonds which have NV center ensembles with
$T_{2}<100\mathop{}\\!\mathrm{\mu s}$, the property makes it feasible to
measure signals at the frequencies lower than $1/T_{2}$ beyond the coherence
limit.

Given a reference frequency $\omega_{ref}=k\omega_{s},k\in\mathbb{N}$, where
$\omega_{s}=2\pi/(mT_{\phi})$ and
$\omega\in\left(\omega_{ref}-\omega_{s}/2,\omega_{ref}+\omega_{s}/2\right)$,
the evolving phase factor can be rewritten as $e^{-i\omega
NmT_{\phi}}=e^{-i\omega_{H}t}\delta(t-NT_{s})$, where
$\omega_{H}=\omega-\omega_{ref}$ is the heterodyne frequency, $\delta(t)$ is
the Dirac function, and $T_{s}=mT_{\phi}$ is the sampling period. Thus, the
readout signal turns to be

$\phi_{r}(t)=G(\omega)\sum_{N=-\infty}^{\infty}\gamma_{e}B_{H}(t)\delta(t-NT_{s}),$
(8)

where $B_{H}(t)=B(\omega)e^{-i(\omega_{H}t+\varphi)}$ contains all the
information from the origin signal to be measured. As discussed in previous
section, the quantum phase readout bandwidth is limited by the cut-off
frequency $f_{c}$ of the filter of LIA. For any signal with a frequency range
of $[(k-1)f_{s}+f_{c},(k+1)f_{s}-f_{c}]$, aliasing can be filtered. Although a
smaller $f_{c}$ makes the measurement bandwidth narrower, it ensures signals
that in a larger frequency range can be detected without aliasing. By changing
$T_{\phi}$ together with $T_{seq}$, we can resolve a spectrum in multiple
frequency bands with a series of sequences.

We present two specific examples of the measured frequency responses by using
the Hahn-echo and CPMG-2 sequence. We plot both the theoretical MW filter
function and the experimentally measured signal responses together in Fig. 3b.
The field sensing time for both experiment and theory calculations is set to
be $T_{\phi}=50\mathop{}\\!\mathrm{\mu s}$. In the experiments, we measured
the amplitudes of the frequency offset heterodyne signals with
$T_{seq}=250\mathop{}\\!\mathrm{\mu s}$, i.e., the magnetic field sampling
rate is 4 kHz. Due to this reason, the measured MW filters are combed with a
frequency distance of 4 kHz. Aliasing signals exist between the main lobes at
a distance of 2 kHz, because the readout sampling frequency is $f_{s}=2$ kHz.

In order to measure signals that distribute in larger bandwidth, we can
increase the sampling frequency, for example, to $f_{s}=5$ kHz. The spectrum
is plotted in Fig. 3c in decibel, from which one can see that magnitudes are
the same at 10 kHz and 20 kHz, i.e., $1/(2T_{\phi})$ and $1/T_{\phi}$ as
discussed in the derivation. The insets of Fig. 3c depict the signals that the
quantum sensor detects during $T_{\phi}$ at the two frequencies. In this
measurement, the bandwidth limited by the filter of the LIA is at 200 Hz,
i.e., the single measurement bandwidth is 400 Hz, and the detectable signal
frequency range is 9600 Hz.

We notice that a single measurement cannot tell if the ac field frequency
offset is positive or negative from the heterodyne readout. Additional
measurement is needed to distinguish the direction of the frequency offset. By
adding a difference to the phase accumulation time $T_{\phi}$ as well as the
sequence time, we can change the reference frequency $\omega_{ref}$ to get a
different heterodyne frequency. By seeing if the heterodyne frequency
increases or decreases, we can determine if the signal frequency is larger or
smaller than the reference frequency. As the measurements presented in Fig. 3d
that $T_{\phi}=50\mathop{}\\!\mathrm{\mu s}$ is offset by a difference of 4 ns
and $T_{seq}=10T_{\phi}$ changes accordingly, the detected heterodyne
frequency of the signal shift in two different directions. We further
investigated the dependency of the heterodyne frequency on the parameters by
performing measurements that vary (i) $T_{seq}$, (ii) $T_{seq}$ and
$\omega_{ref}$. When $T_{\phi}$ keeps unchanged, the heterodyne frequency
shifts by

$\Delta\omega_{H}=\omega_{ref}\Delta T_{seq}/T_{seq}.$ (9)

Using the equation, we can estimate the frequency fidelity of the given
sequence. For example, with a timing error $<3$ ps, the frequency fidelity of
a signal around 10 kHz could be only 0.06 mHz. The frequency resolution can be
arbitrarily high with a long $T_{seq}$ at the cost of bandwidth.

###  II.3 Sensing of Arbitrary Audio Signals

Figure 4:  Detection of arbitrary audio signals. (A) Phase response of the
QPSD measurement. A 20.08 kHz signal with sequential phase changing is applied
to the sensor. The bars show the phases at different time. Stars mark the
readout of the sensor, and the curve is the simulated readout. (B) An ac field
is applied with the frequency, amplitude and phase switched every 100 ms. The
light blue areas corresponding to the right $y$-axis shows the applied field
of around 10 kHz, and the red curve shows the QPSD readout. (C) Spectral
comparison of the applied signal and the detected magnetic field in a narrow
bandwidth. The applied signal is a sum of 20 different sine signals within 400
Hz bandwidth. (D) A signal with wide bandwidth between 10 to 15 kHz is applied
and detected by varying the sequence. We set an 800 Hz bandwidth for the
measurement of each sequence and use 6 measurements to cover the entire
bandwidth. The red dash line shows the spectrum of the output of the AWG, and
the solid black line is the spectrum of the detected magnetic field signal.
The inset figure is the phase-noise power spectrum density plotted within 1 Hz
bandwidth cut.

We demonstrate measurements of arbitrary audio signals by combining the QPSD
readout and the frequency offset heterodyne detection. We first generate a
signal at 20.08 kHz with its phase varying with time (see Fig. 4a). The MW
filter is set by the Hahn-echo sequence with
$T_{\phi}=50\mathop{}\\!\mathrm{\mu s}$. With the reference frequency is at 20
kHz, the heterodyne readout is at 80 Hz, as seen from the simulated curve. The
phase of the external field is switched with a cycle of 80 Hz and 40 Hz so
that the experimental readout displays the phase change, as shown in Fig. 4a.

Next, we apply a field with its frequency, amplitude and phase all arbitrarily
changing. The signal frequency is around 10 kHz and the signal bandwidth is
within 400 Hz. Using $1/T_{seq}=10\mathop{}\\!\mathrm{kHz}$, we can measure
the signals close to 10 kHz with the same sensitivity as the 20 kHz signal.
The signal length is one second and consists of ten 100 ms parts. In Fig. 4b,
both the applied field waveform and the QPSD readout are plotted. The
heterodyne frequencies well resolve the frequency differences in the original
waveform. The amplitudes of the readout also correspond to the applied field
strength.

As discussed previously, the measurement bandwidth used in the experiment is
400 Hz. For this, we perform a spectrum analysis as shown in Fig. 4c. The
signal to be measured is a sum of 20 tones with random frequencies, amplitudes
and phases. In order to distinguish the sign of frequency offsets for each
component, we measure the signal using an alternative sequence with
$T_{\phi}=50\mathop{}\\!\mathrm{\mu s}+2\mathop{}\\!\mathrm{ns}$. The sharp
peaks observed in the spectrum should shift according to the changes of the
measurement sequence, else we exclude them as noise signals generated from our
electronics. As shown in Fig. 4c, the applied frequencies are properly
resolved. Additionally, we find a 9.93 kHz noise spike from the environment.
The errors in magnitude of the measured signal could be induced by the LIA
filter, as shown earlier in Fig. 2d. The errors could also be caused by an
insufficient sampling number for demodulating the rotating frame modulation.
In the measurements, we apply sequences with their lengths corresponding to a
sampling frequency $f_{s}=5\mathop{}\\!\mathrm{kHz}$. The frequency difference
of the two MWs is $\delta f=500$ Hz and $N=10$ for reading out a phase sample.
To increase the measurement precision, if we use a smaller $\delta f$, it will
decrease the bandwidth.

Though the bandwidth of each measurement sequence is limited, we can still
measure a signal with a wider bandwidth by merging several measurements. The
condition is that the signal bandwidth should not be larger than the sampling
frequency to avoid frequency aliasing. In Fig. 4d, we perform a spectrum
analysis for a signal within a bandwidth from 10 kHz to 15 kHz. The signal to
be detected is a sum of 10 components with their frequencies randomly
distributed in the bandwidth. The signal is generated by an arbitrary signal
generator (AWG) and sent to the test field coil. The dotted curve in the
figure displays the spectrum of the electrical signal from the AWG. There are
some harmonics near each main component due to the limited AWG internal clock
and signal length. The components at different frequencies are measured by
varying $T_{\phi}$ to get different referencing frequencies for heterodyne
detection. The inset figure shows the power spectrum of the QPSD readout noise
within 1 Hz bandwidth, from which we calculate the square root of the standard
deviation $\sigma_{phase}=0.0022^{\circ}$. The magnetic field sensitivity
depends on the applied sequence and the corresponding frequency response. In
the case of the Hahn-echo sequence with $T_{\phi}=50\mathop{}\\!\mathrm{\mu
s}$, we have a calibrated scalar factor of
$k_{sf}=0.071^{\circ}/\mathop{}\\!\mathrm{nT}$, and the magnetic field
sensitivity is $\eta=31\mathop{}\\!\mathrm{pT/\sqrt{Hz}}$. Taking the scalar
factor into the calculation of LDR $[-\pi,\pi]$, we can get the dynamic range
in decibels as $20\log[\pi/(k_{sf}\eta)]=98\mathop{}\\!\mathrm{dB}$. The
sensitivity can be further optimized by, e.g., using higher laser power,
applying high-order dynamical decoupling sequences, and implementing flux
concentration.

Finally, we demonstrate the detection and demodulation of audio signals.
Although 10 kHz is within the audio frequency band, most of the daily audio
sounds have frequencies ranging from hundreds of Hz to a few kHz. Therefore,
signals need to be modulated to a detectable frequency range. For this, we
have used (i) a melody piece composed of 3 tones and (ii) a small part from
Dr. Martin Luther King Jr.’s famous speech ”I have a dream”, to test the
waveform reconstruction by the diamond quantum sensor. The tones of the melody
have frequencies distributed between 500 Hz and 700 Hz. Therefore, we mix it
with a 9.5 kHz reference to get the signal modulated around 10 kHz and
broadcast the mixed signal to the diamond. For case (ii), we have to compress
the signal bandwidth into 200 Hz and then modulate it with a 10 kHz reference.
The audio reconstructed from the diamond sensor can be heard and compared with
the original audio. We have provided additional audio files in the
Supplementary Materials confirming the same.

##  III Discussions

In this work, we overcome the LDR limitation of the conventional
interferometric readout through a new technique that includes the QPSD scheme
and the frequency offset heterodyne readout. The technique allows one to
measure unknown signals with maximal sensitivity independent of their dynamic
range. It improves the feasibility for quantum sensors to perform high-
sensitive measurements of different physical quantities using interferometric
methods, beyond magnetometry.

Theoretically, the extended LDR comes from the multiple measurements that have
the quantum phase evolving through the entire phase range $[-\pi,\pi]$ so that
the initial phase factor that contains the external field information can be
resolved. Such an extended phase range affects the measurement bandwidth as
well as the sensitivity. In theory, the sensitivity does not deteriorate a lot
from the conventional fluorescence readout except for a factor of $\sqrt{2}$.
While in the experiment, we suffer from a low contrast $C=0.19\%$ due to the
low excitation laser power (80 mW). The contrast and the fluorescence photon
count can significantly increase when the laser reaches saturation power [39].
Different dynamical decoupling sequences can also improve the magnetic field
sensitivity through the filter function $G(\omega)$. Flux concentration could
further improve the signal-noise ratio [31, 40]. The flux concentrator can be
very small compared to conventional dipole antennas because the gain no longer
depends on the signal wavelength. With the millimeter size diamond dimension,
the flux concentration factor can easily reach a factor of hundreds when using
a concentrator in centimeters.

The QPSD readout can also enhance the capability of vector magnetometry.
Conventionally, fluorescence emitted from NV centers in multiple orientations
is measured sequentially to acquire the vector components. Similar to the
methods developed here, one could also modulate the signal on each orientation
with different modulation frequencies [41]. Performing measurements on
different NV orientations with appropriate synchronization can suppress the
phase errors in vector reconstruction.

In conclusion, we demonstrated high sensitive detections of arbitrary magnetic
field signals, including audio-range frequency signals like melody and speech,
using the QPSD scheme in combination with the heterodyne readout. A further
improvement in sensitivity can be achieved by using flux concentrators. One
could also generalize the current methods to achieve vector magnetometry with
extended LDR. We envisage that the techniques developed here will have the
potential to develop low-distortion, small-volume quantum sensors for various
applications in science and technology.

##  IV Methods

###  IV.1 Experimental setup

The diamond used in the experiment is a (111)-oriented (0.5 mm)3 cube obtained
from a single crystal grown by the temperature gradient method at high-
pressure high-temperature (HPHT) conditions. The diamond is $99.97\%$ 12C
enriched, and has an initial nitrogen concentration of 1.4 ppm. The final NV
concentration is 0.4 ppm after electron irradiation and annealing. Dephasing
time of the NV ensemble is obtained as $T_{2}^{*}=8.5\mathop{}\\!\mathrm{\mu
s}$ by Ramsey sequence, and a decoherence time
$T_{2}=200\mathop{}\\!\mathrm{\mu s}$ is measured by Hahn-echo sequence. The
diamond is positioned at center of a home-built three dimension coils system,
and is illuminated by a 532 nm laser (Lighthous Sprout-G) at around 80 mW.
Microwave signals are generated from two sources (Rohde&Schwarz, SMIQ03B) and
are individually cut by two switches before the combination. Measurement
sequences are generated by a data timing generator (Tekreonix, DTG5274). After
the combination and amplification of the MW signals, MW pulses are feed to the
diamond through a dielectric resonator antenna [42]. The detected fluorescence
signal is demodulated by a LIA (Zurich Instruments, HF2LI) which has two
independent differential input channels and demodulators. To generate
arbitrary magnetic fields, we write signals to an AWG (Tektronix, AWG520) with
$10^{5}$ samples per second output sampling rate. The test signals are
continuously repeated and sent to a copper loop near the diamond.

###  IV.2 Spectrum analyzing

The spectrum to be analyzed is divided into several sections with the
bandwidth set by the LIA for data acquisition. In each section, the center
frequency determines $T_{\phi}$ of the measurement sequence. Usually, the
center frequency satisfies $f_{c}=1/T_{\phi}+\varepsilon/T_{seq}$, where
$\varepsilon=0,\pm 1$. A time trace is recorded after running the sequence,
and a spectrum is acquired from the Fourier transform of the time trace.
However, the spectrum is a fold of the two sidebands with repect to the center
frequency. The sequence with $T_{\phi}^{\prime}=T_{\phi}+t_{clk}$ and
$T_{seq}^{\prime}=mT_{\phi}^{\prime}$ is applied to get an alternate spectrum
with analyzed frequencies shift by $\Delta
f=\pm\left|1/T_{\phi}-1/T_{\phi}^{\prime}\right|$. The direction of the
frequency shift shows which side the signal component belongs to. In the
algorithm, we set a threshold to separate signal spikes from noise, and use
the known sequences induced spectrum frequency shift to distinguish the signs
of the signal offset frequency to the center frequency. The signal spikes that
do not shift according $\Delta f$ are recognized as systematic noise spikes.
Then, the spectrum of the selected section can be replotted as the example
shown in Fig. 4c. After measuring the spectra of all the sections, we can get
the final spectrum by merging them together.

##  V Acknowledgments

We acknowledge financial support by European Union’s Horizon 2020 research and
innovation program ASTERIQS under grant No. 820394, European Research Council
advanced grant No. 742610, SMel, Federal Ministry of Education and Research
(BMBF) project MiLiQuant and Quamapolis, German Research Foundation grant GRK
2198 and 2642, and Japan Society for the Promotion of Science (JSPS) KAKENHI
No. 17H02751.

## References

  * Gustavson  _et al._ [1997] T. L. Gustavson, P. Bouyer, and M. A. Kasevich, Precision rotation measurements with an atom interferometer gyroscope, Phys. Rev. Lett. 78, 2046 (1997).
  * Allred  _et al._ [2002] J. C. Allred, R. N. Lyman, T. W. Kornack, and M. V. Romalis, High-sensitivity atomic magnetometer unaffected by spin-exchange relaxation, Phys. Rev. Lett. 89, 130801 (2002).
  * Hinkley  _et al._ [2013] N. Hinkley, J. A. Sherman, N. B. Phillips, M. Schioppo, N. D. Lemke, K. Beloy, M. Pizzocaro, C. W. Oates, and A. D. Ludlow, An atomic clock with 10–18 instability, Science 341, 1215 (2013).
  * Kominis  _et al._ [2003] I. K. Kominis, T. W. Kornack, J. C. Allred, and M. V. Romalis, A subfemtotesla multichannel atomic magnetometer, Nature 422, 596 (2003).
  * Kornack  _et al._ [2005] T. W. Kornack, R. K. Ghosh, and M. V. Romalis, Nuclear spin gyroscope based on an atomic comagnetometer, Phys. Rev. Lett. 95, 230801 (2005).
  * Kitching  _et al._ [2011] J. Kitching, S. Knappe, and E. Donley, Atomic sensors–a review, IEEE Sens. J 11, 1749 (2011).
  * Barry  _et al._ [2020] J. F. Barry, J. M. Schloss, E. Bauch, M. J. Turner, C. A. Hart, L. M. Pham, and R. L. Walsworth, Sensitivity optimization for nv-diamond magnetometry, Rev. Mod. Phys. 92, 015004 (2020).
  * Shah  _et al._ [2010] V. Shah, G. Vasilakis, and M. V. Romalis, High bandwidth atomic magnetometery with continuous quantum nondemolition measurements, Phys. Rev. Lett. 104 (2010).
  * Taylor  _et al._ [2008] J. M. Taylor, P. Cappellaro, L. Childress, L. Jiang, D. Budker, P. R. Hemmer, A. Yacoby, R. Walsworth, and M. D. Lukin, High-sensitivity diamond magnetometer with nanoscale resolution, Nat. Phys. 4, 810 (2008).
  * Neumann  _et al._ [2013] P. Neumann, I. Jakobi, F. Dolde, C. Burk, R. Reuter, G. Waldherr, J. Honert, T. Wolf, A. Brunner, and J. H. Shim, High-precision nanoscale temperature sensing using single defects in diamond, Nano Lett. 13, 2738 (2013).
  * Michl  _et al._ [2019] J. Michl, J. Steiner, A. Denisenko, A. Bulau, A. Zimmermann, K. Nakamura, H. Sumiya, S. Onoda, P. Neumann, J. Isoya, and J. Wrachtrup, Robust and accurate electric field sensing with solid state spin ensembles, Nano Lett. 19, 4904 (2019).
  * Chen  _et al._ [2017] E. H. Chen, H. A. Clevenson, K. A. Johnson, L. M. Pham, D. R. Englund, P. R. Hemmer, and D. A. Braje, High-sensitivity spin-based electrometry with an ensemble of nitrogen-vacancy centers in diamond, Physical Review A 95, 053417 (2017).
  * Tisler  _et al._ [2013] J. Tisler, T. Oeckinghaus, R. J. Stöhr, R. Kolesov, R. Reuter, F. Reinhard, and J. Wrachtrup, Single defect center scanning near-field optical microscopy on graphene, Nano Lett. 13, 3152 (2013).
  * Zheng  _et al._ [2019] H. J. Zheng, J. Y. Xu, G. Z. Iwata, T. Lenz, J. Michl, B. Yavkin, K. Nakamura, H. Sumiya, T. Ohshima, J. Isoya, J. Wrachtrup, A. Wickenbrock, and D. Budker, Zero-field magnetometry based on nitrogen-vacancy ensembles in diamond, Phys. Rev. Appl. 11 (2019).
  * Wickenbrock  _et al._ [2016] A. Wickenbrock, H. J. Zheng, L. Bougas, N. Leefer, S. Afach, A. Jarmola, V. M. Acosta, and D. Budker, Microwave-free magnetometry with nitrogen-vacancy centers in diamond, Appl. Phys. Lett. 109 (2016).
  * Barry  _et al._ [2016] J. F. Barry, M. J. Turner, J. M. Schloss, D. R. Glenn, Y. Song, M. D. Lukin, H. Park, and R. L. Walsworth, Optical magnetic detection of single-neuron action potentials using quantum defects in diamond, Proc. Nat. Acad. Sci. , 201601513 (2016).
  * Wang  _et al._ [2015] P. Wang, Z. Yuan, P. Huang, X. Rong, M. Wang, X. Xu, C. Duan, C. Ju, F. Shi, and J. Du, High-resolution vector microwave magnetometry based on solid-state spins in diamond, Nat. Commun. 6, 6631 (2015).
  * Farfurnik  _et al._ [2015] D. Farfurnik, A. Jarmola, L. M. Pham, Z. H. Wang, V. V. Dobrovitski, R. L. Walsworth, D. Budker, and N. Bar-Gill, Optimizing a dynamical decoupling protocol for solid-state electronic spin ensembles in diamond, Phys. Rev. B 92 (2015).
  * de Lange  _et al._ [2011] G. de Lange, D. Riste, V. V. Dobrovitski, and R. Hanson, Single-spin magnetometry with multipulse sensing sequences, Phys. Rev. Lett. 106 (2011).
  * Schmitt  _et al._ [2017] S. Schmitt, T. Gefen, F. M. Sturner, T. Unden, G. Wolff, C. Muller, J. Scheuer, B. Naydenov, M. Markham, S. Pezzagna, J. Meijer, I. Schwarz, M. B. Plenio, A. Retzker, L. P. McGuinness, and F. Jelezko, Submillihertz magnetic spectroscopy performed with a nanoscale quantum sensor, Science 356, 832 (2017).
  * Boss  _et al._ [2017] J. M. Boss, K. S. Cujia, J. Zopes, and C. L. Degen, Quantum sensing with arbitrary frequency resolution, Science 356, 837 (2017).
  * [22] H. Clevenson, L. M. Pham, C. Teale, K. Johnson, D. Englund, and D. Braje, Robust high-dynamic-range vector magnetometry with nitrogen-vacancy centers in diamond, Appl. Phys. Lett. 112.
  * Nusran  _et al._ [2012] N. M. Nusran, M. U. Momeen, and M. V. G. Dutt, High-dynamic-range magnetometry with a single electronic spin in diamond, Nat. Nano. 7, 109 (2012).
  * Joas  _et al._ [2021] T. Joas, S. Schmitt, R. Santagati, A. A. Gentile, C. Bonato, A. Laing, L. P. McGuinness, and F. Jelezko, Online adaptive quantum characterization of a nuclear spin, Npj Quantum Inf. 7, 1 (2021).
  * Gurses  _et al._ [2021] B. V. Gurses, K. T. Whitmore, and M. B. Cohen, Ultra-sensitive broadband ”awesome” electric field receiver for nanovolt low-frequency signals, Rev. Sci. Instrum. 92 (2021).
  * Cohen  _et al._ [2018] M. B. Cohen, R. K. Said, E. W. Paschal, J. C. McCormick, N. C. Gross, L. Thompson, M. Higginson-Rollins, U. S. Inan, and J. Chang, Broadband longwave radio remote sensing instrumentation, Rev. Sci. Instrum. 89 (2018).
  * Sedlacek  _et al._ [2012] J. A. Sedlacek, A. Schwettmann, H. Kubler, R. Low, T. Pfau, and J. P. Shaffer, Microwave electrometry with rydberg atoms in a vapour cell using bright atomic resonances, Nat. Phys. 8, 819 (2012).
  * Meyer  _et al._ [2021] D. H. Meyer, P. D. Kunz, and K. C. Cox, Waveguide-coupled rydberg spectrum analyzer from 0 to 20 ghz, Phys. Rev. Appl. 15 (2021).
  * Boto  _et al._ [2018] E. Boto, N. Holmes, J. Leggett, G. Roberts, V. Shah, S. S. Meyer, L. D. Muñoz, K. J. Mullinger, T. M. Tierney, S. Bestmann, G. R. Barnes, R. Bowtell, and M. J. Brookes, Moving magnetoencephalography towards real-world applications with a wearable system, Nature 555, 657 (2018).
  * Gerginov  _et al._ [2017] V. Gerginov, F. C. S. da Silva, and D. Howe, Prospects for magnetic field communications and location using quantum sensors, Rev. Sci. Instrum. 88 (2017).
  * Fescenko  _et al._ [2020] I. Fescenko, A. Jarmola, I. Savukov, P. Kehayias, J. Smits, J. Damron, N. Ristoff, N. Mosavian, and V. M. Acosta, Diamond magnetometer enhanced by ferrite flux concentrators, Phys. Rev. Res. 2 (2020).
  * Seidel and Muga [2007] D. Seidel and J. G. Muga, Two-frequency ramsey interferometry, Phys. Rev. A 75 (2007).
  * Vutha and Hessels [2015] A. C. Vutha and E. A. Hessels, Frequency-offset separated oscillatory fields, Phys. Rev. A 92 (2015).
  * Meinel  _et al._ [2021] J. Meinel, V. Vorobyov, B. Yavkin, D. Dasari, H. Sumiya, S. Onoda, J. Isoya, and J. Wrachtrup, Heterodyne sensing of microwaves with a quantum sensor, Nat. Commun. 12, 2737 (2021).
  * Chu  _et al._ [2021] Y. M. Chu, P. C. Yang, M. S. Gong, M. Yu, B. Y. Yu, M. B. Plenio, A. Retzker, and J. M. Cai, Precise spectroscopy of high-frequency oscillating fields with a single-qubit sensor, Phys. Rev. Appl. 15 (2021).
  * Staudenmaier  _et al._ [2021] N. Staudenmaier, S. Schmitt, L. P. McGuinness, and F. Jelezko, Phase-sensitive quantum spectroscopy with high-frequency resolution, Phys. Rev. A 104 (2021).
  * Staudacher  _et al._ [2013] T. Staudacher, F. Shi, S. Pezzagna, J. Meijer, J. Du, C. A. Meriles, F. Reinhard, and J. Wrachtrup, Nuclear magnetic resonance spectroscopy on a (5-nanometer)(3) sample volume, Science 339, 561 (2013).
  * Bar-Gill  _et al._ [2013] N. Bar-Gill, L. M. Pham, A. Jarmola, D. Budker, and R. L. Walsworth, Solid-state electronic spin coherence time approaching one second, Nat. Commun. 4, 1743 (2013).
  * Zhang  _et al._ [2021] C. Zhang, F. Shagieva, M. Widmann, M. Kubler, V. Vorobyov, P. Kapitanova, E. Nenasheva, R. Corkill, O. Röhrle, K. Nakamura, H. Sumiya, S. Onoda, J. Isoya, and J. Wrachtrup, Diamond magnetometry and gradiometry towards subpicotesla dc field measurement, Phys. Rev. Appl. 15 (2021).
  * Griffith  _et al._ [2009] W. C. Griffith, R. Jimenez-Martinez, V. Shah, S. Knappe, and J. Kitching, Miniature atomic magnetometer integrated with flux concentrators, Appl. Phys. Lett. 94, 023502 (2009).
  * Schloss  _et al._ [2018] J. M. Schloss, J. F. Barry, M. J. Turner, and R. L. Walsworth, Simultaneous broadband vector magnetometry using solid-state spins, Phys. Rev. Appl. 10 (2018).
  * Kapitanova  _et al._ [2018] P. Kapitanova, V. V. Soshenko, V. V. Vorobyov, D. Dobrykh, S. V. Bolshedvorskii, V. N. Sorokin, and A. V. Akimov, 3d uniform manipulation of nv centers in diamond using a dielectric resonator antenna, Jetp Lett. 108, 588 (2018).

## Supplemental information: Low Distortion Radio Signal Sensing with a
Quantum Sensor

###  V.1 Derivation of the rotating frame modulation

Based on Eq. (2) in the main text, the evolution operator with $\pi/2$ pulse
duration $\tau_{1}=\pi/(2\Omega_{1})$ is

$e^{-i\frac{\mathcal{H}_{1}^{\prime}}{\hbar}t}=\frac{1}{\sqrt{2}}\begin{pmatrix}1&-ie^{-i\left(\delta\omega_{1}\tau_{1}+\alpha\right)}\\\
-ie^{i\left(\delta\omega_{1}\tau_{1}+\alpha\right)}&1\end{pmatrix}.$ (S1)

In the sensing duration, the spin state acquires an extra phase factor due to
the external field and can be expressed as

$\ket{\psi(T_{\phi})}=\frac{1}{\sqrt{2}}\begin{pmatrix}-ie^{-i(\delta\omega_{1}\tau_{1}+\alpha+\phi)}\\\
1\end{pmatrix}.$ (S2)

The evolution operator during the MW2 pulse is similar to (S1), and the spin
state after the second $\pi/2$ pulse is

$\ket{\psi_{meas}}=\frac{1}{\sqrt{2}}\begin{pmatrix}e^{-i\Phi}-ie^{-i\delta\omega_{2}\tau_{2}+\beta)}\\\
-ie^{-i\Phi}e^{i(\delta\omega_{2}\tau_{2}+\beta)}+1\end{pmatrix}.$ (S3)

where $\Phi=\delta\omega_{1}\tau_{1}+\alpha+\phi+\pi/2$,
$\tau_{2}=\pi/(2\Omega_{2})$ is the $\pi/2$ pulse duration of MW2. Then, we
can calculate the measured expectation value as

$\braket{S_{z}}=\frac{\hbar}{2}\bra{\Psi_{meas}}\sigma_{z}\ket{\Psi_{meas}}=\sin\left[\phi+\frac{\pi}{2}\left(\frac{\delta\omega_{1}}{\Omega_{1}}-\frac{\delta\omega_{2}}{\Omega_{2}}\right)+\alpha-\beta\right],$
(S4)

where $\alpha-\beta=\delta\omega\cdot t$,
$\delta\omega=\delta\omega_{1}-\delta\omega_{2}=2\pi\delta f$

###  V.2 Microwave filter function

Here we give a general switching function based on CPMG-n sequence. Typically,
a microwave sequence consisting of $\pi$-pulses can be described by a
switching function

$g(t)=\left\\{\begin{matrix}1,&t\in\left[0,\displaystyle{\frac{T_{\phi}}{2n}}\right)\\\
(-1)^{k},&t\in\left[\displaystyle{\frac{2k-1}{2n}}T_{\phi},\displaystyle{\frac{2k+1}{2n}}T_{\phi}\right)\\\
(-1)^{n},&t\in\left[\displaystyle{\frac{2n-1}{2n}}T_{\phi},T_{\phi}\right)\\\
0,&t\in\left[T_{\phi},mT_{\phi}\right]\\\ \end{matrix}\right.,$ (S5)

where $k\leq n-1$ and $k,n,m\in\mathbb{N}$.The switching function is repeated
in each time interval of $\left[NmT_{\phi},(N+1)mT_{\phi}\right]$, where
$N\in\mathbb{Z}$. Given a random sinusoid signal $B_{ac}(t)=Be^{-i(\omega
t+\varphi)}$, the readout of each cycle is given by
$\phi_{r}(t)=\int_{-\infty}^{\infty}s(t)g(t)dt$. For each sampling point, the
accumulated phase factor is

$\phi_{r}(N)=\left\\{\begin{matrix}\displaystyle{\frac{4\left(\sin\frac{\omega
T_{\phi}}{4n}\right)^{2}\left(\cos\frac{\omega
T_{\phi}}{2}\right)}{\omega\cos\frac{\omega T_{\phi}}{2n}}}e^{i(-\frac{\omega
T_{\phi}}{2}-\frac{\pi}{2})}Be^{-i\varphi}e^{-i\omega NmT_{\phi}},&n\ is\
odd;\\\ \\\ \displaystyle{\frac{4\left(\sin\frac{\omega
T_{\phi}}{4n}\right)^{2}\left(\sin\frac{\omega
T_{\phi}}{2}\right)}{\omega\cos\frac{\omega T_{\phi}}{2n}}}e^{i(-\frac{\omega
T_{\phi}}{2}-\pi)}Be^{-i\varphi}e^{-i\omega NmT_{\phi}},&n\ is\ even.\\\
\end{matrix}\right.$ (S6)

The equation is simplified as

$\phi_{r}(N)=|G_{n}(\omega)|e^{i\left(-\frac{\omega
T_{\phi}}{2}-\frac{P}{2}\pi\right)}\gamma_{e}B(\omega)e^{-i\varphi(\omega)}e^{-i\omega
NmT_{\phi}},$ (S7)

where $G_{n}(\omega)=|G_{n}(\omega)|e^{i\left(-\frac{\omega
T_{\phi}}{2}-\frac{P}{2}\pi\right)}$ is the MW filter function of the
sequence. When $n=1$ and 2, i.e., Hahn-echo and CPMG-2 sequence, the sensor
responses are

$G_{1}(\omega)=\displaystyle{\frac{4}{\omega}\left(\sin{\frac{\omega
T_{\phi}}{4}}\right)^{2}}e^{-i\left(\frac{\omega
T_{\phi}}{2}-\frac{\pi}{2}\right)}$ (S8)
$G_{2}(\omega)=\displaystyle{\frac{8}{\omega}\left(\sin{\frac{\omega
T_{\phi}}{8}}\right)^{2}\left(\sin{\frac{\omega
T_{\phi}}{4}}\right)e^{-i\left(\frac{\omega T_{\phi}}{2}-\pi\right)}}.$ (S9)

When $\omega=2\pi/T_{\phi}$, the two sequences have the equal response as
$|G_{1}|=|G_{2}|=2T_{\phi}/\pi$. The phase response difference is $\pi/2$.
Therefore, the Hahn-echo sequence can be used to measure the $X(\omega,t)$
component for the quantum-LIA, in correspondence to a classical LIA. The
CPMG-2 sequence can be used to measure the $Y(\omega,t)$ component.
Furthermore, we can get $|G_{1}(2\pi/T_{\phi})|=|G1(\pi/T_{\phi})|$ from (S8).

###  V.3 Derivation of shot-noise limited sensitivity

According to the schematic described in the main text, we try to extract the
phase factor from the N samples in a modulation cycle. The shot-noise
contribution to each sample is $\delta\mathcal{F}=\sqrt{2\mathcal{N}}$, where
$\mathcal{N}$ is the detected photon number in each measurement, and there is
a factor of $\sqrt{2}$ because each sample include two measurements. The
output from the photodetector can be described as

$s_{k}=\mathcal{N}C\sin\left(\frac{2k\pi}{N}+\phi\right),$ (S10)

where $C$ is the signal contrast and $k$ denotes the number in the samples.
Since the photodetector output include the shot-noise, the measured
$\overline{s}_{k}$ include the noise
$\delta\overline{s}_{k}\approx\delta\mathcal{F}$. We can define a minimum
error square as

$\chi^{2}=\sum_{k=1}^{N}(\overline{s}_{k}-s_{k})^{2}.$ (S11)

Additionally, we have $\sum_{k=1}^{N}\sin^{2}(2k\pi/N+\phi)=N/2$, and
$\sum_{k=1}^{N}\sin(2k\pi/N+\phi)\cos(2k\pi/N+\phi)=0$ when the each
modulation cycle is equally sampled with an even number. Let
$\partial\chi^{2}/\partial\phi=0$ to minimize $\chi^{2}$,

$\sum_{k=1}^{N}\overline{s}_{k}\cos\left(\frac{2k\pi}{N}+\phi\right)=0$ (S12)

Take a differential,

$\sum_{k=1}^{N}\left[\delta\overline{s}_{k}\cos\left(\frac{2k\pi}{N}+\phi\right)-\overline{s}_{k}\sin\left(\frac{2k\pi}{N}+\phi\right)\delta\phi\right]=0$
(S13)

Sample $\overline{s}_{k}$ can be replaced by $s_{k}+\delta\overline{s}_{k}$.
We can neglect the high order error $\delta\overline{s}_{k}\delta\phi$ and get

$\sum_{k=1}^{N}\delta\overline{s}_{k}\cos\left(\frac{2k\pi}{N}+\phi\right)=\delta\phi\sum_{k=1}^{N}\mathcal{N}C\sin^{2}\left(\frac{2k\pi}{N}+\phi\right)=\mathcal{N}C\frac{N}{2}\delta\phi$
(S14)

The statistical average of the phase variances square is

$\displaystyle\braket{(\delta\phi)^{2}}=\frac{4}{N^{2}}\frac{1}{(\mathcal{N}C)^{2}}\sum_{k=1}^{N}\sum_{l=1}^{N}\braket{\delta\overline{s}_{k}\delta\overline{s}_{l}}\cos\left(\frac{2k\pi}{N}+\phi\right)\cos\left(\frac{2l\pi}{N}+\phi\right)$
(S15)

The measurement errors should be uncorrelated, satisfying
$\braket{\delta\overline{s}_{k}\delta\overline{s}_{l}}=0$. Therefore,

$\braket{(\delta\phi)^{2}}=\frac{4}{N^{2}}\frac{1}{(\mathcal{N}C)^{2}}\sum_{k=1}^{N}\braket{(\delta\overline{s}_{k})^{2}}\cos^{2}\left(\frac{2k\pi}{N}+\phi\right)=\frac{2}{N}\frac{1}{(\mathcal{N}C)^{2}}\sigma_{F}^{2}.$
(S16)

The short-noise determined phase noise in the modulation cycle $2NT_{seq}$ is

$\delta\phi=\sqrt{\frac{2}{N}}\frac{\delta\mathcal{F}}{\mathcal{N}C}=\frac{2}{\sqrt{N}}\frac{1}{C\sqrt{\mathcal{N}}}.$
(S17)

Taking $\delta\phi=\gamma_{e}G(\omega)\delta B$ into the equation, we can get
a sensitivity

$\eta_{Hahn}=\frac{\sqrt{2}\pi}{\gamma_{e}C\sqrt{\mathcal{N}}}\frac{\sqrt{T_{seq}}}{T_{\phi}}.$
(S18)

In order to compare with the sensitivity of the direct fluorescence readout,
we can get the minimum detectable phase when the signal satisfy the small-
angle approximation as

$\delta\phi^{\prime}=\frac{\delta\mathcal{F}}{C\mathcal{N}}=\frac{\sqrt{2}}{C\sqrt{\mathcal{N}}}$
(S19)

Therefore, the sensitivity of Hahn-echo measurement with fluorescence readout
will be

$\eta_{Hahn}^{\prime}=\frac{\pi}{\gamma_{e}C\sqrt{\mathcal{N}}}\frac{\sqrt{T_{seq}}}{T_{\phi}}.$
(S20)

The shot-noise limited sensitivity of the QPSD readout deteriorates from the
value of the fluorescence readout by a factor of $\sqrt{2}$.

Taking the values in experiments, that $\mathcal{N}=4.6\times 10^{11}$,
$C=0.19\%$, $T_{\phi}=50\mathop{}\\!\mathrm{\mu s}$,
$T_{seq}=100\mathop{}\\!\mathrm{\mu s}$, and with
$\gamma_{e}=2\pi\times(28\mathop{}\\!\mathrm{Hz/nT})$, we calculate that the
shot-noise limited sensitivity using QPSD readout is $\eta\approx
4\mathop{}\\!\mathrm{pT/Hz^{1/2}}$.

###  V.4 Frequency resolution

The following is the detailed description of the measurements of Fig. 3d that
discusses the frequency resolution of the experiment. Firstly, we apply
sequences with $T_{\phi}=50\mathop{}\\!\mathrm{\mu s}$,
$50\mathop{}\\!\mathrm{\mu s}\pm 4\mathop{}\\!\mathrm{ns}$ and
$T_{seq}=10T_{\phi}$ to measure a signal at 20.005 kHz, denoted as Meas. 1 in
the figure. The reference frequency shifts 1.6 Hz due to the changing of
$1/T_{\phi}$. As a result, a heterodyne frequency of 6.6 Hz can be detected by
measuring a 20.005 kHz field with $T_{\phi}=50\mathop{}\\!\mathrm{\mu
s}-4\mathop{}\\!\mathrm{ns}$. Correspondingly, a 3.4 Hz signal will be
detected when measuring a 19.995 kHz field so that it becomes distinguishable
to the 20.005 kHz field. For fields at frequencies near the other references
$\omega_{ref}$, the shifts can also be calculated accordingly by
$k/(mT_{\phi})$. To further understand how sensitive is $\omega_{ref}$ to the
sequence parameters, we only differ $\Delta T_{seq}=\pm 4$ns to measure
external fields at 20.005 kHz and 16.005 kHz while keeping
$T_{\phi}=50\mathop{}\\!\mathrm{\mu s}$. In this case, the heterodyne signal
is introduced only by the accumulated phase error in each cycle. The reference
frequency shift is calculated by Eq. (7) in the main text. For example, by
measuring a signal around 16 kHz with $\Delta T_{seq}=4$ ns, we detect a
heterodyne frequency shift of 128 mHz as shown in Fig. 3d. Therefore, the
frequency fidelity is majorly limited by the jitter performance of the pulse
generator. We calculate an example that has a timing error $<3$ ps, according
to the datasheet of the pulse generator we use (Tektronics, DTG5274). The
frequency fidelity is 0.06 mHz when the measured signal is around 10 kHz.

###  V.5 Audio signal measurement

The first file is a piece of melody that we composed with tones C5, D5, and
E5. The corresponding frequencies are $f_{1}=523$Hz, $f_{2}=587$Hz,
$f_{3}=659$Hz. The signal waveform plotted in the time domain can be seen in
Fig. S1. We choose the three tones with their frequency differences within 200
Hz because the used quantum phase modulation frequency is 500 Hz and the LIA
limited bandwidth is 200 Hz. The signal is modulated with a carrier of 9.5 kHz
so that the transmitted magnetic field is at frequencies around 10 kHz, which
can be detected by the sensor using the Hahn-echo sequence with
$T_{\phi}=50\mathop{}\\!\mathrm{\mu s}$. Then, the QPSD readouts of the three
tones will be at frequencies $f_{r1}=23$Hz, $f_{r2}=87$Hz, $f_{r3}=159$Hz. In
order to reconstruct the original signal, the readout is mixed with a 500 Hz
carrier. Three notch filters are applied to the mixed signal to remove the
sideband as the final step.

As for the speech signal, the signal is distributed over a frequency range
from DC to a few kHz. We cut off the spectrum from 4 kHz, shown in Fig. S2(a),
because signal amplitude after 4 kHz is much smaller than the components at
lower frequencies. In order to make the signal detectable by the sensor, the
signal is firstly compressed in the frequency domain by using interpolation
while keeping the same sampling rate for broadcasting. The signal size after
interpolation is 20 times the original signal size so that in the frequency
domain, the bandwidth is compressed from 4 kHz to 200 Hz, as shown in Fig.
S2(b). Then, the signal is mixed with a 10 kHz reference so that we can use
the same Hahn-echo sequence to measure. Figure S2(C) shows the spectrum in
which the audio signal spectrum distributes in symmetric to 10 kHz. Finally,
the QPSD readout detects the heterodyne signal, of which the frequency is
referenced to 10 kHz. To reconstruct the original sound, we compress the
readout in time domain, i.e., averaging the readout by every 20 data points.
In the frequency domain, the spectrum bandwidth is expanded from 200 Hz back
to 4 kHz.

S-FIG. 1:  (A) Signal waveform of the melody played for the detection. The
voltage signal is sent to an AWG for generating the magnetic field through a
loop antenna. The different colors represent different scales of the tones
that are labeled in the figure. Changing the amplitude indicates that the
loudness of the sound rises at the beginning and then falls in the end. (B)
The QPSD readout of the sensor. Different tones can be recognized from the
different heterodyne frequencies in the readout.

S-FIG. 2:  (A) Spectrum of the original speech signal. The signal is written
into an AWG in the unit of Volt. Therefore, the amplitudes of the spectra in
sub-figure A to C are all in Volt. (B) Frequency compression. The signal is
interpolated to a data size of 20 times larger than the original one.
Therefore, the signal is compressed in the frequency domain. (C) The signal
spectrum after mixing with a 10 kHz reference signal. (D) The QPSD readout
spectrum. The heterodyne signal is sampled every 20 data from the raw readout
to reconstruct the audio signal spectrum.