From 2e93517c19a3207daa703cdce7e752448e67c433 Mon Sep 17 00:00:00 2001 From: Ruslan Kabatsayev Date: Fri, 22 Sep 2023 22:51:22 +0400 Subject: [PATCH] Add Lens Distortion Estimator plugin --- CMakeLists.txt | 1 + .../LensDistortionEstimator/CMakeLists.txt | 5 + plugins/LensDistortionEstimator/COPYING | 340 ++ .../resources/LensDistortionEstimator.qrc | 6 + .../bt_LensDistortionEstimator_Off.png | Bin 0 -> 11622 bytes .../bt_LensDistortionEstimator_On.png | Bin 0 -> 22047 bytes .../resources/icons.svg | 4423 +++++++++++++++++ .../src/CMakeLists.txt | 73 + .../src/LensDistortionEstimator.cpp | 874 ++++ .../src/LensDistortionEstimator.hpp | 167 + .../src/gui/LensDistortionEstimatorDialog.cpp | 1822 +++++++ .../src/gui/LensDistortionEstimatorDialog.hpp | 173 + .../src/gui/lensDistortionEstimatorDialog.ui | 916 ++++ src/core/StelApp.cpp | 4 + 14 files changed, 8804 insertions(+) create mode 100644 plugins/LensDistortionEstimator/CMakeLists.txt create mode 100644 plugins/LensDistortionEstimator/COPYING create mode 100644 plugins/LensDistortionEstimator/resources/LensDistortionEstimator.qrc create mode 100644 plugins/LensDistortionEstimator/resources/bt_LensDistortionEstimator_Off.png create mode 100644 plugins/LensDistortionEstimator/resources/bt_LensDistortionEstimator_On.png create mode 100644 plugins/LensDistortionEstimator/resources/icons.svg create mode 100644 plugins/LensDistortionEstimator/src/CMakeLists.txt create mode 100644 plugins/LensDistortionEstimator/src/LensDistortionEstimator.cpp create mode 100644 plugins/LensDistortionEstimator/src/LensDistortionEstimator.hpp create mode 100644 plugins/LensDistortionEstimator/src/gui/LensDistortionEstimatorDialog.cpp create mode 100644 plugins/LensDistortionEstimator/src/gui/LensDistortionEstimatorDialog.hpp create mode 100644 plugins/LensDistortionEstimator/src/gui/lensDistortionEstimatorDialog.ui diff --git a/CMakeLists.txt b/CMakeLists.txt index 593e66d2b6dd4..20addd2c18b6e 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -549,6 +549,7 @@ ADD_PLUGIN(Satellites 1) ADD_PLUGIN(Scenery3d 1) ADD_PLUGIN(SolarSystemEditor 1) ADD_PLUGIN(Supernovae 1) +ADD_PLUGIN(LensDistortionEstimator 1) # Candidate to removing as an archaic plugin ADD_PLUGIN(TextUserInterface 1) ADD_PLUGIN(TelescopeControl 1) diff --git a/plugins/LensDistortionEstimator/CMakeLists.txt b/plugins/LensDistortionEstimator/CMakeLists.txt new file mode 100644 index 0000000000000..d8aaf5d484620 --- /dev/null +++ b/plugins/LensDistortionEstimator/CMakeLists.txt @@ -0,0 +1,5 @@ +SET(LENSDISTORTIONESTIMATOR_VERSION "0.1.0") +ADD_DEFINITIONS(-DLENSDISTORTIONESTIMATOR_PLUGIN_VERSION="${LENSDISTORTIONESTIMATOR_VERSION}") +ADD_DEFINITIONS(-DLENSDISTORTIONESTIMATOR_PLUGIN_LICENSE="GNU GPLv2 or later") + +ADD_SUBDIRECTORY( src ) diff --git a/plugins/LensDistortionEstimator/COPYING b/plugins/LensDistortionEstimator/COPYING new file mode 100644 index 0000000000000..b35f35c99338e --- /dev/null +++ b/plugins/LensDistortionEstimator/COPYING @@ -0,0 +1,340 @@ + GNU GENERAL PUBLIC LICENSE + Version 2, June 1991 + + Copyright (C) 1989, 1991 Free Software Foundation, Inc. + 51 Franklin Street, Suite 500, Boston, MA 02110-1335 USA + Everyone is permitted to copy and distribute verbatim copies + of this license document, but changing it is not allowed. + + Preamble + + The licenses for most software are designed to take away your +freedom to share and change it. By contrast, the GNU General Public +License is intended to guarantee your freedom to share and change free +software--to make sure the software is free for all its users. This +General Public License applies to most of the Free Software +Foundation's software and to any other program whose authors commit to +using it. (Some other Free Software Foundation software is covered by +the GNU Library General Public License instead.) You can apply it to +your programs, too. + + When we speak of free software, we are referring to freedom, not +price. Our General Public Licenses are designed to make sure that you +have the freedom to distribute copies of free software (and charge for +this service if you wish), that you receive source code or can get it +if you want it, that you can change the software or use pieces of it +in new free programs; and that you know you can do these things. + + To protect your rights, we need to make restrictions that forbid +anyone to deny you these rights or to ask you to surrender the rights. +These restrictions translate to certain responsibilities for you if you +distribute copies of the software, or if you modify it. + + For example, if you distribute copies of such a program, whether +gratis or for a fee, you must give the recipients all the rights that +you have. You must make sure that they, too, receive or can get the +source code. And you must show them these terms so they know their +rights. + + We protect your rights with two steps: (1) copyright the software, and +(2) offer you this license which gives you legal permission to copy, +distribute and/or modify the software. + + Also, for each author's protection and ours, we want to make certain +that everyone understands that there is no warranty for this free +software. If the software is modified by someone else and passed on, we +want its recipients to know that what they have is not the original, so +that any problems introduced by others will not reflect on the original +authors' reputations. + + Finally, any free program is threatened constantly by software +patents. We wish to avoid the danger that redistributors of a free +program will individually obtain patent licenses, in effect making the +program proprietary. To prevent this, we have made it clear that any +patent must be licensed for everyone's free use or not licensed at all. + + The precise terms and conditions for copying, distribution and +modification follow. + + GNU GENERAL PUBLIC LICENSE + TERMS AND CONDITIONS FOR COPYING, DISTRIBUTION AND MODIFICATION + + 0. This License applies to any program or other work which contains +a notice placed by the copyright holder saying it may be distributed +under the terms of this General Public License. The "Program", below, +refers to any such program or work, and a "work based on the Program" +means either the Program or any derivative work under copyright law: +that is to say, a work containing the Program or a portion of it, +either verbatim or with modifications and/or translated into another +language. (Hereinafter, translation is included without limitation in +the term "modification".) Each licensee is addressed as "you". + +Activities other than copying, distribution and modification are not +covered by this License; they are outside its scope. The act of +running the Program is not restricted, and the output from the Program +is covered only if its contents constitute a work based on the +Program (independent of having been made by running the Program). +Whether that is true depends on what the Program does. + + 1. You may copy and distribute verbatim copies of the Program's +source code as you receive it, in any medium, provided that you +conspicuously and appropriately publish on each copy an appropriate +copyright notice and disclaimer of warranty; keep intact all the +notices that refer to this License and to the absence of any warranty; +and give any other recipients of the Program a copy of this License +along with the Program. + +You may charge a fee for the physical act of transferring a copy, and +you may at your option offer warranty protection in exchange for a fee. + + 2. You may modify your copy or copies of the Program or any portion +of it, thus forming a work based on the Program, and copy and +distribute such modifications or work under the terms of Section 1 +above, provided that you also meet all of these conditions: + + a) You must cause the modified files to carry prominent notices + stating that you changed the files and the date of any change. + + b) You must cause any work that you distribute or publish, that in + whole or in part contains or is derived from the Program or any + part thereof, to be licensed as a whole at no charge to all third + parties under the terms of this License. + + c) If the modified program normally reads commands interactively + when run, you must cause it, when started running for such + interactive use in the most ordinary way, to print or display an + announcement including an appropriate copyright notice and a + notice that there is no warranty (or else, saying that you provide + a warranty) and that users may redistribute the program under + these conditions, and telling the user how to view a copy of this + License. (Exception: if the Program itself is interactive but + does not normally print such an announcement, your work based on + the Program is not required to print an announcement.) + +These requirements apply to the modified work as a whole. If +identifiable sections of that work are not derived from the Program, +and can be reasonably considered independent and separate works in +themselves, then this License, and its terms, do not apply to those +sections when you distribute them as separate works. But when you +distribute the same sections as part of a whole which is a work based +on the Program, the distribution of the whole must be on the terms of +this License, whose permissions for other licensees extend to the +entire whole, and thus to each and every part regardless of who wrote it. + +Thus, it is not the intent of this section to claim rights or contest +your rights to work written entirely by you; rather, the intent is to +exercise the right to control the distribution of derivative or +collective works based on the Program. + +In addition, mere aggregation of another work not based on the Program +with the Program (or with a work based on the Program) on a volume of +a storage or distribution medium does not bring the other work under +the scope of this License. + + 3. You may copy and distribute the Program (or a work based on it, +under Section 2) in object code or executable form under the terms of +Sections 1 and 2 above provided that you also do one of the following: + + a) Accompany it with the complete corresponding machine-readable + source code, which must be distributed under the terms of Sections + 1 and 2 above on a medium customarily used for software interchange; or, + + b) Accompany it with a written offer, valid for at least three + years, to give any third party, for a charge no more than your + cost of physically performing source distribution, a complete + machine-readable copy of the corresponding source code, to be + distributed under the terms of Sections 1 and 2 above on a medium + customarily used for software interchange; or, + + c) Accompany it with the information you received as to the offer + to distribute corresponding source code. (This alternative is + allowed only for noncommercial distribution and only if you + received the program in object code or executable form with such + an offer, in accord with Subsection b above.) + +The source code for a work means the preferred form of the work for +making modifications to it. For an executable work, complete source +code means all the source code for all modules it contains, plus any +associated interface definition files, plus the scripts used to +control compilation and installation of the executable. However, as a +special exception, the source code distributed need not include +anything that is normally distributed (in either source or binary +form) with the major components (compiler, kernel, and so on) of the +operating system on which the executable runs, unless that component +itself accompanies the executable. + +If distribution of executable or object code is made by offering +access to copy from a designated place, then offering equivalent +access to copy the source code from the same place counts as +distribution of the source code, even though third parties are not +compelled to copy the source along with the object code. + + 4. You may not copy, modify, sublicense, or distribute the Program +except as expressly provided under this License. Any attempt +otherwise to copy, modify, sublicense or distribute the Program is +void, and will automatically terminate your rights under this License. +However, parties who have received copies, or rights, from you under +this License will not have their licenses terminated so long as such +parties remain in full compliance. + + 5. You are not required to accept this License, since you have not +signed it. However, nothing else grants you permission to modify or +distribute the Program or its derivative works. These actions are +prohibited by law if you do not accept this License. Therefore, by +modifying or distributing the Program (or any work based on the +Program), you indicate your acceptance of this License to do so, and +all its terms and conditions for copying, distributing or modifying +the Program or works based on it. + + 6. Each time you redistribute the Program (or any work based on the +Program), the recipient automatically receives a license from the +original licensor to copy, distribute or modify the Program subject to +these terms and conditions. You may not impose any further +restrictions on the recipients' exercise of the rights granted herein. +You are not responsible for enforcing compliance by third parties to +this License. + + 7. If, as a consequence of a court judgment or allegation of patent +infringement or for any other reason (not limited to patent issues), +conditions are imposed on you (whether by court order, agreement or +otherwise) that contradict the conditions of this License, they do not +excuse you from the conditions of this License. If you cannot +distribute so as to satisfy simultaneously your obligations under this +License and any other pertinent obligations, then as a consequence you +may not distribute the Program at all. For example, if a patent +license would not permit royalty-free redistribution of the Program by +all those who receive copies directly or indirectly through you, then +the only way you could satisfy both it and this License would be to +refrain entirely from distribution of the Program. + +If any portion of this section is held invalid or unenforceable under +any particular circumstance, the balance of the section is intended to +apply and the section as a whole is intended to apply in other +circumstances. + +It is not the purpose of this section to induce you to infringe any +patents or other property right claims or to contest validity of any +such claims; this section has the sole purpose of protecting the +integrity of the free software distribution system, which is +implemented by public license practices. Many people have made +generous contributions to the wide range of software distributed +through that system in reliance on consistent application of that +system; it is up to the author/donor to decide if he or she is willing +to distribute software through any other system and a licensee cannot +impose that choice. + +This section is intended to make thoroughly clear what is believed to +be a consequence of the rest of this License. + + 8. If the distribution and/or use of the Program is restricted in +certain countries either by patents or by copyrighted interfaces, the +original copyright holder who places the Program under this License +may add an explicit geographical distribution limitation excluding +those countries, so that distribution is permitted only in or among +countries not thus excluded. In such case, this License incorporates +the limitation as if written in the body of this License. + + 9. The Free Software Foundation may publish revised and/or new versions +of the General Public License from time to time. Such new versions will +be similar in spirit to the present version, but may differ in detail to +address new problems or concerns. + +Each version is given a distinguishing version number. If the Program +specifies a version number of this License which applies to it and "any +later version", you have the option of following the terms and conditions +either of that version or of any later version published by the Free +Software Foundation. If the Program does not specify a version number of +this License, you may choose any version ever published by the Free Software +Foundation. + + 10. If you wish to incorporate parts of the Program into other free +programs whose distribution conditions are different, write to the author +to ask for permission. For software which is copyrighted by the Free +Software Foundation, write to the Free Software Foundation; we sometimes +make exceptions for this. Our decision will be guided by the two goals +of preserving the free status of all derivatives of our free software and +of promoting the sharing and reuse of software generally. + + NO WARRANTY + + 11. BECAUSE THE PROGRAM IS LICENSED FREE OF CHARGE, THERE IS NO WARRANTY +FOR THE PROGRAM, TO THE EXTENT PERMITTED BY APPLICABLE LAW. EXCEPT WHEN +OTHERWISE STATED IN WRITING THE COPYRIGHT HOLDERS AND/OR OTHER PARTIES +PROVIDE THE PROGRAM "AS IS" WITHOUT WARRANTY OF ANY KIND, EITHER EXPRESSED +OR IMPLIED, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF +MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. THE ENTIRE RISK AS +TO THE QUALITY AND PERFORMANCE OF THE PROGRAM IS WITH YOU. SHOULD THE +PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF ALL NECESSARY SERVICING, +REPAIR OR CORRECTION. + + 12. IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING +WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MAY MODIFY AND/OR +REDISTRIBUTE THE PROGRAM AS PERMITTED ABOVE, BE LIABLE TO YOU FOR DAMAGES, +INCLUDING ANY GENERAL, SPECIAL, INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING +OUT OF THE USE OR INABILITY TO USE THE PROGRAM (INCLUDING BUT NOT LIMITED +TO LOSS OF DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY +YOU OR THIRD PARTIES OR A FAILURE OF THE PROGRAM TO OPERATE WITH ANY OTHER +PROGRAMS), EVEN IF SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE +POSSIBILITY OF SUCH DAMAGES. + + END OF TERMS AND CONDITIONS + + How to Apply These Terms to Your New Programs + + If you develop a new program, and you want it to be of the greatest +possible use to the public, the best way to achieve this is to make it +free software which everyone can redistribute and change under these terms. + + To do so, attach the following notices to the program. It is safest +to attach them to the start of each source file to most effectively +convey the exclusion of warranty; and each file should have at least +the "copyright" line and a pointer to where the full notice is found. + + + Copyright (C) + + This program is free software; you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation; either version 2 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program; if not, write to the Free Software + Foundation, Inc., 51 Franklin Street, Suite 500, Boston, MA 02110-1335 USA + + +Also add information on how to contact you by electronic and paper mail. + +If the program is interactive, make it output a short notice like this +when it starts in an interactive mode: + + Gnomovision version 69, Copyright (C) year name of author + Gnomovision comes with ABSOLUTELY NO WARRANTY; for details type `show w'. + This is free software, and you are welcome to redistribute it + under certain conditions; type `show c' for details. + +The hypothetical commands `show w' and `show c' should show the appropriate +parts of the General Public License. Of course, the commands you use may +be called something other than `show w' and `show c'; they could even be +mouse-clicks or menu items--whatever suits your program. + +You should also get your employer (if you work as a programmer) or your +school, if any, to sign a "copyright disclaimer" for the program, if +necessary. Here is a sample; alter the names: + + Yoyodyne, Inc., hereby disclaims all copyright interest in the program + `Gnomovision' (which makes passes at compilers) written by James Hacker. + + , 1 April 1989 + Ty Coon, President of Vice + +This General Public License does not permit incorporating your program into +proprietary programs. If your program is a subroutine library, you may +consider it more useful to permit linking proprietary applications with the +library. If this is what you want to do, use the GNU Library General +Public License instead of this License. diff --git a/plugins/LensDistortionEstimator/resources/LensDistortionEstimator.qrc b/plugins/LensDistortionEstimator/resources/LensDistortionEstimator.qrc new file mode 100644 index 0000000000000..54ba44f5781a2 --- /dev/null +++ b/plugins/LensDistortionEstimator/resources/LensDistortionEstimator.qrc @@ -0,0 +1,6 @@ + + + bt_LensDistortionEstimator_Off.png + bt_LensDistortionEstimator_On.png + + diff --git a/plugins/LensDistortionEstimator/resources/bt_LensDistortionEstimator_Off.png b/plugins/LensDistortionEstimator/resources/bt_LensDistortionEstimator_Off.png new file mode 100644 index 0000000000000000000000000000000000000000..d2cb4c5cf067187ad6ca2a8b812c6d69c6f6f6cb GIT binary patch literal 11622 zcmb_?_d8te_w_L*j6PbBh#oyi2$Bd!^iBw(_vkHp8xt*}6QTn{7=Ywfl6drfr(GU7YL0059FDavYt&)ok$Faq%V!=UX0 z@B#IZQPP2dR{+c^68xXYRnfo$0F;LQ``~?*i#i2=WOys5|MsPe&0BACcWc1g+grfS z*}=oo+|^pZ#oabzTjCA?umDQ32pylytw6(Q#xdVz?{$;WleGrd9@Zycd=mJzsN_^C zq2}y(L&UFGo*LWhCN; z{H~yYs`i{XyyJn&*o_6;QGg^)2X%cVGa8F2HRV=HRF&aBlqlc0Xl2&Bo5vl6Vu`#A zZr$#KJ*xTenN_|@b&dQ%eet>W`Yw3 zK_}SS+t(NQQs?>QL?|q#dXf2hj``ZQYNEQL+ZH0N|N6Dxv+4l_Y<|Vas5sP$tJ?S- zEg%ZeHDQXr#T8{g_dj;@Svo&IuUzrXXj$;?7MFC$7wkMZ3lL~V6zMd^%syq@%vxmw zJS{C>;0wuOi;9X!k&lp~$3f4POGbxEM%QgB^fR~0IYK3l$1Nw~B+plGcQs{S?^VS7 zU8u!9?qX$GN|_U835C>$%@#bO+9rSiw9?z-6_zu^r>|dJ2`GVXCh$ z&Zm^vUhY_jvnhRbLE83Mgwm*^9!dk3sCO*VvDa3f;Rro^tnf5zprqKZaCps;rQUSm zt7WT-h%xWSvVABNsleoGwVfX)*+`rlHbp|(ZR*7@ zK)XBkoqmZlq zthwa5PBtaMZGlc7N13;n76m$K7E8Az@`;@Gp zpnypUt?jKgY(!{g?k#IORF0fAiEbaCN>+k@CpNFmsR|Du8FL=qNEuZKU&JSU5s zYu-JmORqDGJ5Qm}qT+OE_m_LaB?H7=3?A=d^zr`Dzq}%yb$uV}{6F<1*v$vGQ9;A? zF7@NMF{oMX>(bGQiHX89vdW1~Hpo!Lh3LTJ2H3VSzI3T}(c~f#3^+JA7{oi>F?#mw ze&g>Z0>TBmOhP6JXU%PoLXvkdX%)2gN(Z5_iODq-f8O(9e9%HQI}op{`<%(y*3xo8 z6h%1FdHTcBEZl%jy6JW4G?`1T-9NUxDp*sCz~As&xdyx+1tM_Z;tsNEqA!JKoSsz- zPHD*8iX*{823%jB2`6%p6%~Q|b-h_NSjHC5@hpJ~7t{0U2||E|_=%hrd4!lgcwG*q zwjS2JeY(XHHBTv{9SWy>w^l&-86V3Qs;O)@brRTC>x3D#x90EV4Hes_)a=B@)@E5X z&U6&QOg@Ps|NeJ2_o@aZRebLj7udxyJ@>1wZ*Fc5Yk3m?BPt>m`QgI{!7#lM`x?XS zmB&KyL7V2gV%nz*Qr8Plv}kPQw2a_Tb8;*z0#2@_OhGBLjMoSBP|LW(nEarH8adY# z4^60LwL-12dC+=qH1@$HQJxweA6 zHdTY+dd|+yLFCw={hwSO8WMXy5*$UJ3BuYR%Ji_YT2W%zSQoi#Sx3t|l9^aJ;%1k` z@0oOS^3uq)eGqIN58Ugd8w;2*WeuI&3k;q(xw0Pfbor=^9zB@}p7CVtF4o~4^Myd9 zEA#^>vB57~%NJHHHW~cyvvf{(`kFp~2#Kg&X#2dQ=2CBV-J9<$r;?GU+N>j`LRz9& z81Z+(Y#hjQ(Ze1$yZAGyY@_kyX_mU|$Bo8o5tLR#Xq)1xv{L;HZ@tD}5;_6JL5YJt z1}&N#0@o0@J2`eCkwY#zE4@Ij?>>?IJ2)F(RAnx=vQ<{uhArmiYUf4|V5XLl+zY&U z*;^@jb%D9L@TcdTZ%uGdgj&56Uo5SgR(`rg|fr@)_QJB`VXP zvn7Q%W~m@~-BB>}NGQO3Nf8*!yj_XGDkt1J7{FdlPPp2sTK-;9wJN7{n;^MYzf~}@8A<}(mnTaA~Zgj0*YUV zMxzUME|~F8tqX_GmOk?7IR~B2MN~+2)Dv-EujO2?Nr4CKAv2xytdr%jnL*0L8<9Sn zD!ae-kE^tcs?+-m!a}zCyRjxc{wt9vmapOFtCtl#=xQ*Xv}4oqdbjEfU1*n6I9xxF zV!;7M(R&Q9RsV4m9n0wF6!tr~h`P+%XiOYj84ii4E#0lp^nyOP)eNl zD_}stO(468WgGRbbCB$y^EE%!?rrv0kYk(0!m+%csChG zc%w9F7M* zo{D(P|9Br+5oi>dNAEv2C%uIP$w!QQzEs#Q`uaA5!4%U)h5=7Qe4yzB**v&7;K_x)4UrX7@hUYCLZ?F5h|IEcS z^R%ydUGuv1uc(+<$;tLqNm+ThUfRmu z+Hm1Y2UJ9@DhgXw5N(z{+>S&>zu6palcw^UY%ejrh1IFhkN)Jr*j>aDA3*E4LsI_| z5*|zfgDCu>P67e5IXtpw$eC{D!gLnR|004Qq%Z#+`Fl(aO4Ti6xYHdoSss{4A`#=7^;bc;8*v#r0zJg3C{ z4p%of$}vp7>i&pHs^1vjgXI-iz`*O4FmERd684-Av!#4CRTGgq_9A!cD8*|l7s6V6 zyWOUfY2@C+e^-7}BmJ_*fsI@`{&xqJ$L#^zmKz_TZ0*|J5kgFeWI%IYm5;eZ;C`(f z_(B7x>|G3R;o%)kTd>6EZ317)c z?S!vQO)+61%F%pqrB^ifWsav|$!R8>52>+}1&#eS!ZOkZDTNOi4QewM_RXde$&tq` z{wK~6-oHJsgixuC(b6oyyX*Vx?Ce&=2r3jC^2|7MEEZnJ5YEOj>jcea{5IGn|0O^w zrlx%Pz*jA0iPgmrWcts%x4)T|t?gc7W znoq<#ldXXl0_bWL*qW4soL1EE@NkiPLtg-_DYDp26rx92pVw~&AumdbM;gM6_Rrpw zyE4aE%+mAl3kV3Xyd5?kEA`FpK2LD-b%p^PN{MW(S|9+SeiobMj!_tap)CD1m`h2U z2_Lde+b8v#h%;id(wx3IoRv!>Tf#&jwcNZu?$afd`9=bDfB>n*kJi54^t8tDgJ@cRwpme+gC@7R+koQQ z+ov*E^}-H$v|9oSfG@AC6v7*~o$DMdd}=wI#qM9XdEwPhIoHunul6Lc4@WB8I&ePE z(#m)45F)n}dbgiz8D{&=qNA1mt^zV=BUwc4!&_%?$R5OzDpyoM`^bon6)@lEOroXq zoJi>Bkf(7SgGpMy&4);3wA&vpb1?o%{d9D#^PQ)ot9Hb>A&(+qF#5DRRs0UgdcW#d zPi*BGKKnzh7BgtVM-1z5G4S6)(QNbVvN^~CyI8haR?S6EqLn03n$G#V$xYVQxn5FB6*F@U+MAS0ck#PU}<``!E@-eY?KK~PmGOd$w zaYftE>-gZ2uOlLAK{jmUL_Gwo3AGY~mgK7L3V!f!23LLHCP&ZXdjh~4WO!}TnI(Ww8^IUC|bPp6r_2)_LB@06fJb0i59Zq;Wb8}v5lv>6d zZD1br`f!KW)kT85T9N}MGu}TLe6dw`Q&nIGvPD};j_HP*wqPhkim#%gLZ-ad;#_A~#oxWo=;cC0T$!KWR>oNtQXPu;%(DX3qSK0in>nW><2j)K3>r7swWY6)3o zIYTp{u8*446FK-FJg5nr!L0aB+gD8ltXj&1E_BD;bc@G55onxp9`iJgFq9ZPXbzAx zRYYsQaO10Uo@%UVZ8h;SP7O2jx?<51@nl>HQW7d!5I2_7VoR{E^l0Q#tPN5)l(J=MCIvndp=ar3#w1 zzW+1qhC2#lT1@S}CJ}@HzBV=m^qgi*Q_dJL{B|tISyA-ipBchgHyX>n7W0-A^E#w= zN!IA_qK|SN(hGQ^@;c*v0u;s_W}e-KY@@RY(IXVy%pzcHV6HbXLI2)}e3Q>N|GVM2 z9McxvX-*d)b%-NzQu)It$f~6MV_X%+gJ1tXQRMY z6dPhz>4!VG3nu*_SoR6uXWagIy|SR7$eHv9WTV;rHb*+VAN@8fy4BUy)W|;m1?f4p zoSVi}epctYOPoIJlx4pL=hQT6U!ywg_&7tY>d%``Aa$y#sj&u&jK1c}{a(|#Iv0>n z^@6L>hGqB|=e+Pj=xlInUUvMQ4jzrQL2SJt5$Zbu@;=szBKR3a<`%~Sw+DgXQFyzs z!bkO0%3TuDy!}s$Hph*-?H5~QVQ!)sYDE){%4lsRg4}_=zAx?2%f%!O0n^5}Ys-g+ z)mq#%Lk&g)Js?*GQ^q4wb|-py`v~YgL(;Bd6c1Y8z`$TT?6oM!T0J~Gs%&tgByNlE z5qObt68nl)k@7r`DWiLeVG_k>0B5!jlRUK$V=f0Py%R2f+g1Yt7T%J#^M>8?znF}d zM&FR$bhlQqNexdqi&05Eh{{8<5WCWb?4>`63ct^&fpz`Ds3A)mVi=gM?I$pq&EH$* zgxl!t?Y-XLe(@fwCiiEjyb=6Ia9$w$3LEvB}hw@!W)1 zIw-4!%pZ}`7?475O1Lo)83ZWarsf_LC4x-L$D1QN_rir~h1k{rMlG7sT|8N>Lm}8q zwNPVir#_T{7N-oGP_ZtY%9GOE51gTvO9&>A9TusWa^B&XOz!@gs*FYI+(wU(yRh_D z#^vHFVnb+1OT2Eq6JOB!XRq-|Hsz8Eh2?0~C-bk5mSuymEyGxGoo3PNy;lWhh89PS zu%@k%xJI$S6+K|e3H9qaDQYf*UV=GjM1XBo;q^Et{+5Xl5z!kqO^h0}*qg@pTwuB- zXc6uw=0FvC_emCw?Ow+{Go92B!WpHkkwAJ+0se5Sk+V~n$Sc{XZQ3LjmLwzGlO%FV zvk{InlfXVEq!zCa2EFmMl^1!fD89W>_!$!3mXjX?GLJ@|q*vd>-=qedwlX&mQH##+ zh|TRpXhAc38ltt^l83SiGBfDO6OjSYvh>s*Cb`M-q}WP+=APt5(JENdeizfxH;QJL zz*(x_2hBanI)hTKR1_$`Y&ur6dW&T!rRW=%LR3hgE7g?wz{YzEL1^^S;=S(oIwLet zbBiT1prws9m<^- zt94?(h;kK919u2x_w~U8xG*9e#3HK|{R6MQtk?0o8#FK`?A@!j5x2~q-`VW4oB3_i zhEk!~ZzI`h-0NT}(iq%U*OVYT0M(+Y*(DzHJkch3a+~px<%hTNNf0#p;BuP_k9Lz{ zFMQsj-T6)HPF%Y4q=P>cNZ{a+!8+IVwrd4DD73XY;SzW!pVKnh+Gb9f1t&w3Io{XQ znGVH$(Jtyy%+bD6o6#jlihze|@qSXm*26FkKh&__Pj;@P@qWK5p;a1UhoM4Fi4c_0 zs2wi66BXh_E&LoxXf!?qcoc=)^18s@pUvJhWtaFIwi9&HbfW|hjg!;5l!b?C+>Pl{ z6OeyQ8VuHBOadAAqb!sYIV!+esY@%baceAIaj*a^Nk3(E*v6hmwouTkD-i-&>V;t! zgx5)8GXWl8nIyG-aSf_M@~HEZsQPz_n%p$4#Iq%vQXN8#bb=<$B(ntQLrx?pz)cm6 zZ8$IooqqOR8MrriawRII%pMxw(9kfszHXnX9Us%k&g;W`fn2R&rw^xO4PEi&xAV_* zsz4}m)HZ@b8m%0qsqfOqFm{DlMlD@WEqG`WpbmokVcybwob3E}W?_Je}BkN#Qk1(j_dOnvFTI6N`4Ji*M_HXMjxHYuhooz>-piD z6yfd)sCkj1lc!(IuD;GoPT@h{S6<^pzfFQ`?@R1tTN^BW{P?jS6oFY(M?Oc0*u#fP zUNkemi*J?H)_pw0wV{OFJr4X2d*-^-u#$zD1 z-KLb=Q}qe%5L!0%nG355lR@y@{!qXjr7J&zxoETI5+TK}v+^>UM2u4Ez9k4PM*C@oTvQa5LVQf(~hDng_$GP1Jzm%&#G zq^h8Y*WPj_fwEXpb=gK&&?v_TlIKO$Q;k;=AIE=hc^MCj7(V}W&r~S?H5AdsO^bFa^E%(-3_$Zyn(E?{Ivi-Wj6IArakLcA&PYhO z$A!CaflnUm+Om316S;EwX975EZ%+&ZG)qfa=S z848r7g(et3>2})`!t-dzv^MyAC%lRh__G-O0!k`=7cp}daxnlh-``*@J$a#CSQRyk zb}*NY)Us31-whDY%^tO|*;)sdJ&oh&KReX?ql2AqIMC~`Ep>35(sgyPJXZepDfyo3 zNSM%f!k=H#*L<_r#MO0&MOxz!>N%25m+N|%#)#ZRUnxVHTD`Vrzh&xkMS9tjL2%<_ zzit~y=9GhB!v`Co^dMdxtGerUjWhN>W!g=QA6?D@OKP@{j<{CxSzeUtEq)fM%*4|tZdh5n~C04(BXoF1JBX2A^5cGMzYgp&>VG!k5Ar>|^ zHcV|V%;ZSU>K^fbM9yXM%uVv#fzZx3l^H5<|H5c;S5tod!YhC1DPQf0aah0{tnjn? z5Z$&iAV2Z9Ycauu<=*rqa_oi9S9!)nJs>X$xP=Tl-kf-0{xWqzd^DR~sCU@I?2lW; zGxkjgjnPeK<+qps$s2AAf-}hgD}5d+%;7jT{?3M!0quu}hZVcwwN5T(EXE7ml?e7t z03aRvuNMINMD$o0+mbE9{@EdYxn);bZTR;0K%(KTe@sVGNh#!jtxlGx+HfWPp(p2` z8S}B~pHdQW_qyN4AHC=OW&Ahxe$Jo*af<6L>XUSmh#jupaC2EIVld%T5eQBM z72JP6n_Jj+4&_hn<(vB$mDgukpSNH#$Id-TRg_awyiUq#S|10%L3x53Qt8*}hQ`40 zdO(+kco#IwbV=oZ6GlaZAqPM=lQ@+Z0@--HsH(JI(2XSnax(U__4 z_1TKSlhX^`@j0h30c5PT0i^}aht6Txd9hp3?DoHZ|KtUMuM`F;B#mVUe|(aV$XqMO zP-Kk!0@M7Np4uNl$$H)R8Hub1|Fez_l~3L30W8W(u9NY zW(BiB+m;G_$gaE=(Yg~{U0dt52Kq4%`RcIw8d%U)c*TD^g7Go9dPJnrzpQ>-KfA^5 zKX%XzwiS@+>YlZJa=kp=>-f+6by9tGO9fp1Fuh|Ky@Kz*pl2^%zNB2if$s?adKI@t z2rOpITo%0MzdcBV+JzwLV63K;(798l0XsugM;cXF8dqol8Foxo>&Md6fka&t>dOIrUyjed?4m^6r=h3$7i0Je5rKMA8V(PRp~+#62h2sEJSw1u&PBhg?N zn$%#tunVU&NRW}~(YuAfetr@0pF;@xi{+wsVzjF!qM_?=(!Cf%@Vv|IeR_65;gD~_ z(W|H7QAYoO<#-6|@otYQXTD0JhzOY?&sN#DZ$b;MjI{G|fb!1LJy3xDRo8k{x&s3H z6N40r=9C!(4i5o)Qq2A9&)S-us)oW;KaJdz!x)!a0L~L17_XJ0Ixq-O5Lso-+1`Vb1yk1X2y*(7CVN zQ3`&0O>wTD%kxXzy8&G71LHBvqo~;DD%UgDTQk=jS?awCaZFk1IXOAmX=!Po?tVu9 zE=fnGhe4)ClhBYxuBS5EV`TuO%vreHez*mY>)Llv8|om+FiF=VItaQP%mlJD>pH2@jJ4|IO>3jpl=Wi zqJaT*+7F!pk~xB*3p5s}$XWA#{0?%eL=ME5XWxK)7XTo2c#nRZfIhGHEyUMj;wF$l z9ZkEX&mFdcDk1^dIL7 zvShRvByyuMF)?6&gU-bGxLGW=?G)vNKM?|cZdE6SR^Wzsf&yQL`mB>Lf0^DpT}y#?F{-efhSWp|2c{DWhR}ee1POpur35=@+_mrm7yYzS;fMu{wj54Q>+ar7)bXnip+#FCo;dtI}ljKWx;Ks}bXr);52jors)8;Nf&906QAVT`Js&W2FxxWNZ%BXu6shK z@^LU?i>fm-o9j7C6`!%-<4MlM9(Sqpmz!UFaYp0DCfMB9mVQo7Tl))O|IqE=c$v?x zP{a>?+O~dMf$hM3w;b)&XUwe=<$yuSC@~NY0}BM7Ezo|Wf3A0@qBo{q{rLGNn_%nd z{0EJQV}+^-es*rPUJDNe-7s4#j|Pn(ELM>N!VN-XBnf!5&&8l~1mA8y?3JVEeD%qS z`Y0tEiszlz;S*<>gX@O;ZVb9Qw=HkpkLsX{*P^jr2Z3|0WGk+*e~Y27i#1YHpm1WFayQ}og1%$eWOdQ{^?_!qUT z{8^Jez2v#)`Uv${nehc3smyH#pHU|ADbffxW!QO!2iQX~a&qNduKfh;;1~s)3=Fo^ z{#C87neRW3L|&3~vF#UtS|vvdXpFanBa!z>6$hJ^BSu}tWnvKoEK)}!FIV91Ji=-o z{gvMykaBE>JpH`yqYvwGF8DCZX>jCsW%CJn7I19`<*4<5MJ zI;Da^HzO4Xml?=>k`D^GcDLFhA7=jEM?gubU&zGcq3(oJx^;k(tsdpqPJbAzr`v8) zkTfX*X|kzWp)bK~Vo>_Z#X|i01ncnQp#(n^O5oCfkNwX%ZLP0Ay@Lrcr`Ir>O@F*T z7kA!|t+l1gB0(c6xO;-*1J71Qq>r) ze%Wd9FmCqq_FogKJtEmzs3Eb-IEgh){nlCQD{n#-aV@xkzJ8>5M0&wb(&KM&V$Mi$ z1RO8*kH~|e0T3JoEbINQUc6JHj2~2#_NqTzu)%v2ItOh3^9d%_Z)m<@}hdAq0LB4)e`Pt>5fU>}jN#>u|% zi8;=zfeybzWo(Km*y60p+g~#!&-)HG4Xox*7NpQnLL>nb3j+@tg4TgUSrbvG~^Y9L=-Zo3?xe)k*>CKhCT%OwL}#4uTrn3cL_GyE4@jwo;y#Xe7~%7I za6k6TakOD$ownK&qKvnbSw>2vMD*9wweeq4*CcZJthN>KPP|(RwjNX003CghA1rfSo-hIKo9Oy z#yl>9543nB+KK_ZMKHJ$!2cQj4DIj$fS&kwr>fFR{tZ6K7pP|&XnDsy@UC-!8*umT zU3pJmZ#>S~&rSYLfJfms4SoO+0nn&xRv|@S+K$rX`u(r)x7b+UH5>abn024yO7^o` zog_D-L}!5vs-^oG!k@#5vh0T!|1gEo)4Y!6%b>axL%Z|hsf6Wox2*g4(floPO+Zp> zv@0q&IOwO{tHBaYDi-&O2SYW(10NKl*nzu&2LBJ+0iU5&^b^v=c%g8__t@jJ!;bR+ zp#A6bd2XOnw!ZZ-0LaVBa{--~qJA*{JUiL#zAqeA)&3*(PUp{`O>4XVc~1v+I-PD` zm#U`%B2IUlqPJ7eHm-Dp7E_Ld&BhA)3pO?m>Rv7;PtF$M*PBibcUQ-v&+s89Utg}8 zMIT#+94_1y13KMbBmuyAVJ#v6ES)SKiUOVIgd>2Zzh{T(PvP+ZKmg9sYq+$Z*ZQpk ze$C&LrTfBX-|ugSFP#Ye*)-6`4hP#nCOb!Bo**! zaj`gZH}qDUV~O#Hzv9QO`{(?k*vO>ypI?Jh0N`+UDcydDA%pLQ`#}3Hv(nkR)0rE1 z2BI1MwPodBO8~zf0+K*OQ$ZTieKz$(!f_;I#=h6MpXH@}ow#Vt{QQ4Gf6pR6>>f@) z%K#wvoy6_Sg#o7vzLoNqE@|%;EsQ@v!Q#c-T_q0|`+HqtXJ==hIo6kv6n$U*=b9Ke z+G~(5^qIu&on&j&ZAV|5UvJZ-vKo0XZqV#m&=t5Cjw4zSsv2j!1P$htQyuDOzpo8) zt`}wV^6~;1cXzy+WX{hXo>CUuRDy0?d6V5i@E^T#A{d#EiyGgzNpD*4O&4ELO$!MM z;=Xg|j7{NxeD2Br z;cKuQnu@(QkfYCEj&nX}A#mc|lF4KQ%dR|+3>-HbboROXeidmCJ9LZq1EK$N!S~18d*{+5R)XK8?zW%!x_S;6 z|Gp>ur;`dobGCd=`xzx~w`tZ<2}Y}uRu;OqmUN?j*8RoYhwDy})?kB`->qX)S8si6 zR*76HtniE{VHsg-s%c>B6CfNwr&{I5<^+F12|#DQ82}*2jPHDC?^FQjRwz0h=2aP; z?8{3#wAPrM?5kQ@-FYaG((hSNLG-Rmg-d{^lO0h(g0Qcurd98#YqQ+DQ_zng$4{Qf&MK$Q+V`hx9S@8l9zj+*W8f11DG5LMG~$laI{S@#V=G5IT$Cn% z320F!eE1T4l^IBdbOqk1e2CYFliLpUy)T08WLE*-5=tfk!1u)`6`q}Hca@X=oLNfB zrLKk?6f`v(`<*&YDzESIUN{Gc^o4NoEc;iin)RR5|1k#o$K2}l+S?J5VY}Bb`L=o$ zB!&2xzC4|dNFu`4J96ppk!OK3G%<$oZEr8&)7N0FzRkt3ox0CXe}6jtJM-VbW4TI6 z#zg)8?)>LT$MF*(K!-->N8sP%`oH#V2NOoy^Gj#H*!OQZIOHw(&eqL(x6Ob3NWzN6 z-%5a7?)w15Bq8kKWFF)#lrjbmi#Nx>;qj2ga0xOOT!MhW;8;R-$q?#1^j2*8SxMHR z5Lu(Oj)dHeJnD$5k6M7zAo`A9kCDu=Lt=xU+xoYl*@xG>S@#ZT%r3~r`zg;nI5LLO znq!UgrGYqdZ9bM$HbU+bxk4eh9dGU00NvsLXxGp!@n z+#^WC(ZQtIe6Pzm156!fOanZ`4v9W)ngI_I4fiCIT(&4hlUPh{BNi@EPSWt#aOMC2 z_Jo=BT=!`LB#u{Z{IAEfkx8rhPiRCW&iF|1(z6+UkG?q4gk}E|zPOT;(0*hOY8C)O zHt^H95zXr+8qs}uilz)yM{Q6sC{7_BKoU?gR~TU~7&sG5fy_n1Ix>7_ABR2OU%b+n z=ZL9pJ9LNA2l@LWzz%hjq7)CA-3fB2v-$q*&2FFY>83$1&kS(i+-Npvd$BECF}~F1 zNseHDgMw2bYqdOi>g(YN9ngp!vhnu85ihvHZ^c9G86IIV{6IMb$p~}BphX{->XG0Q z&O|Co5ejvUQw)vFO<1K{?Zsjw$X3o8sCX+ngk8||c2YkR%#h3#51|FP&;m$lX1!}z z2(t80h=hgr+48&(KaV7=5>8zOkKbsdAU3LT2r+3ghxT=bQ9%a%#Xibow^eX(|K zvj%(9uILLA{byAKnUDT80kl-YCB9UVP}^BI(=W#WIg-Q zeO^00{P0`L)15TI22tOrA1O{%ccM?`hwPsOzjqCv!y7 zqaVj%Wqwf5bzJVnjwNT>Sy*F)J_Af#D7)T>CE#&E{x)klDUOr{K@LEu%SrLrN>hc4 zm+w?5nySz%p)GpTm)rUuI7O*xg&hpYkl#b^!6n@4dP~j@D*IVN{QVy} zM*g&Z)!v?)F09?0(u;&Wwv&G#8Gp03iB=W|Sc(j@5 z$)B~K9Y=*7M;)n%tlNePeWeKcyKind$UpdS7zwUsr%~reTYiJUfLFkr+I8#j^L=3% z?>hTBn~Cr{<)r4_@NBFXWV{dh#u;jl`VtRe6)ir0d0{--r~s!0#8i0lZ!n@7U7PK2 z#I2J>LwFrRzs^S1v(8#L>i6dc)_)(|?fxPYMw>twf|eCjfPdx2NRnL}XTaI9f@BZe za^KxZJL|W=R8s*~7(LXN+;?seWM2EO0FF3oA70n=>GtJhjdk3m-(V-OS`2rqvwra4 z0l1P*AGTdX;p(R6zc$;+@mwamQgwROv5GDtmMQZ5aMIV$8N5W#7Y~`$o7Dn(U61EpwHxu^D zmXSRYv(-ZV*!a&;HKCF-XX|IfuduBXkNQYcAFJJW)nZ*R0;T0p4U1pwi^gFw)&7EY zGav6mY7y>M7(U?rjQ!kU+qnIG@N2F=$pKIM`Wt>=862&9*E^1LUPu4c`t;`Rc_afk zcwLBz8W&7hs>(_H?A)T;TdC!wm8G+Sz@ zHss*Xu!YE-|CYjDwns6*&K8zV7yjy_hPx$x@mOA{ry1}^xX5K{OQuJOJ+{Im-HF;S zyi+^7(=0^s&_0>#FNr$dGQ}hz8Hfl;rZ#cFkW1>Si}~w2Q~)*moEYiPUx7QV{zs!% zcU!U16mU6z3BGAe>$AtrG0p&%XCYOse2|8au0OX4wi}aB+}|~`v$eLcNiu07J^IgT z^pC}JucM2OK6~i$-*bAG9g*kp``g`)lW3C}4cOz$eO1u1<{1p7h#oLO4l}@}*4HKC zAva&nUs#V9>dk>z_*@nMk-(P|M->yAX}hO!MOUhBxtIm3B^a8MmMQHpMSUSmOUcis-|4%-S*3HI5HF%l+!>%BJk>RD|6$ALCn+QAkS)i zt7GY$OKn%@7fa*%8->&dM00d?7fEBo_vKX=9U$gsn9s7v58C&M1?@I>j>A~sZuu!B zj)TbEL}yj(AqE}_-Y)wVns^yTY~Lfuv762OfcF?gLy(Z4p2J>Ean;+J#}kEj+CvXS zOVrXt^E|uKT@<03LBN+tgUni$Yv|nD2zCXn*I;{QvzGN_FtLoo5w z>D}-XjlnW~&;kc02x{0~-$gstfHQGw(=P(Th4^a0*T7u{fqpAWFUyg`dtQdy`iWBW zyN}J9p=av>o*<<9Tq__lC`vZY>Bio!vZtetVSx3JvX*-olUc->R2~k`xilyoL@K7x znHXecd{u0AN;l3`L+4{T2-3wGp%ayZEnfp#z7P?~FG6KjSgEl?0|Mz<-8KoM$8;78i0m$PL_dxKsbP45128aQ#gMrfmNeUK|AdIB@cGY&*f3-7LqwQN z-og>b?PWzn<->p9jXu2Vq<%%piC)E|tm9;1x1A41TnuM{p(qX_d!sGOMf;;_+G%A1 z{Rp5lb5+;5q>H$CGIe{Wg_!7Q`J)D2 zL7W^+I?mq^$e8&^-JOSZ?woZDXMI{k*{3o`o?^a@T;;gBuk~{60;S!^+T6dKqBCG`l0-EfD-mO3X6L{lp0Mn28 zim$;}iL0B^i%=Y!eIFEu>VfXv6GJDua8A1?Vuv!=uvpAPY*rtXbuQl`rT4Mt{I8eI zkHCh3I!jt3NRW4f-#B+Hod+-ODLAVhM(!^A`DF}zF0g*<+HDjTr8jGZfkXQ(?849d z#}q|`n-JHH;Cc+#%1Pp(nv^^Pd7hMBa9-J4S#fEwugi5!iFlgPw{MPNY&%?3qYTCz zE}p2S;rvD);pNpS3VnHJo9SnpU#s~~jzHB)_~f+ROzom(IZIy>uuLUVQFJDpbS<+@ zxAZ&vT$9YY7@Q8UXLy^BVgjb7r>D0Im;P+b|9k}k%)fJEQNIM#{jq&e&$^s76uhp< zv!Kt`yO3A$S#dTIv9Y~9at6%B|{V(J@B481B@{94rTW zDg;?cQhwOrg|doV+N~G^7jwQuX9#jR_>CS4RTtFLnrv9?Neeui*>3g$<&@W#Ehm*t z=9_&^^c1VGn7Y}W>997($gkc`G~oK!g%-@-2oL^1@nNiVy)9>u=Y0M(k&@Q#kQw#c~q%9rWv)@8kYUfRmNe!BDVXINWn^xrdH8%B)x zsnz+1U)nt!D}P-n&bKhY9lc|-~+7~WBI$x`Z zp6qHGsYoDnD-eNFt>euysp#X%=-0)*&FR`_zmu9~K=N{x^*CDSvB)8FBqg4zu;T64 z;QzL(fwc+9^fgt*j3Ye;_Eq1P0^qp^J8R;IyV8tW0Eh=t;52EgstmA*-WVZ(8iCBf zCXzTRNTOw@RLPsy#p%>SZL&2h3o6tsSQk&OVbq?N(*T#^{b01`L(SmCYp-wA3`MshWtiG@XTUf zDxH4lB-lw^8dL++g+SPb$K#cNHBJnWwRjt)mwK6T7(6Gd&?uzXU0mnFBEk!yc((RA>QN|{0> z5m6a3^&K7RKU{UXlx>G-DQR*b0~LCrf^=hFnMGIIt+Ux@q|6}wm!gOCKTWu$KV5j~ zxyPk=2q^V_4Qu;m^JLnk$^K%&Yw3}|2K;FridzUv=6et2taQ#RC^n)Kz_}3Fi|`^0 zh;*gABog)*oEIp6+vaWuTeD_o#lCF`6HKBUVli=S>*==y`!BdwL-LM*PEIk^EFwD~ z`tbD{)eF`NBm>6#Axd8M(N>6piV0pZ_pto|mR@A=0XCpK*Z?;SMSgE&9C|f9&h;ch zZ#c0JN*_m~tt}i1ijZ9~KJ)x27C2d&%;Kt>2#J#_RLX+r6ez(_^6k3slLoSE7-H>r zxR8b?r%}=Mx$+rfZxSFykH?J8=?qrGXhZ(|Zs=%fX(21EfPCvx+X41)F?{!Q@nDf1 zp<>ca&0J2pNlY4BWz^0AJ95(E3D)tgJxKN8@#d}13$U22AE0E%Nv7!Z-sHYFmCv)j z`*g=N_(N33615tQ&ahOx9~L7O5Sb+%PeACeB{0Ga>w?O$>`8AYlrv)prA=CVPj*j= z`Izx)&%u{~Q-#bG2AE(kZOT_h0D`kJne?+x%m`47qs>YTEQ#IFQ>U0tkU;@ONjM`= zs$MB%w`@3xI(F<9dfmw74=pf2g{xA0f6~5W+z;Xq8#uS5DnFYR?}^Zp+d;dLzd`x% zF&C3B;T`FVs`viOGNJWdS_@ihW)-wsbS6QnX)#|99uNgNYPhW={tf&o-@O+a!9oQ&ZK z`Y9*H;V0OiInjslf||$CTFtd1^mFK3sBtW_X7VC}Pt&`WjNs_LDTbjV0e2a*%O;cN zBfNm>sK)0JQ{yBc2BAN=%@(j{pxnn$P-tjC20>5mMpXz+UtZ#Nx8Lald+ym$L0s@F z_vJ4IcpjB{8fGdWrrC#`-TB3wpjxxB2Tj781(>5uJYlp2`=2QTLY9hZuO>zJnm~D{ zAT{4t`yHUXU~^T~<8SD0f5F?BavFhGO|#4yMmhfW<_xix_fQPil4!mgTub{{V%X^n z4H;MUHn2cBU6X(2+IpBm~XxSO%PK&+{DJAzSgqM2)GzimdGv>LT`Ab(euS^g+y^{Z-jafQ{S^Jg4*< zXl-TiNWzQ)>Y8#1Rq+sh;FZ@Cjw2G%^8Ji7Iqmim$tc5L)v1!i5}LyU5J?jZFmMOH z4u+xud@$AM7CAlH_*+=E7G(tQjxsDhu@#&B{K`}KZBh7Uge`&C;;ZV45^CBwNk`bD ztkl?^N~poG91z0orLy;oV-NA)P%Hd6er$uZ(%C(g>n@2%Y7)Y7OsyDb8UAlX0AErq z+dgw(M?G%Fs&l=PeZ91J5#uM>tbA7$1VOi`GxJRByv5LssI<(rhy5v~B$;BEtNU|A zABV^{#nGE#OdMFPx63bl-Gd>>m_tWs^yis;*{6nG1rgo!@rySTyZt{N-##xzSZ7FtF`>{a0EolMHxi>`LB6Zwj9`AY3` zQ5)&g?m9Iy0l5{e2$pV3<_elt5t?>s>np<|C9knawfcwqG54MQZeM|k-ZvlWvE*Z* z0|tKSW6*t2?rM{Z)h@DYy|Wm&8&NGk6IymZ!bVmv4IW_gvOH|I`1#cc;_l(OHl3U|smG*ZC5mrg0zHN(J=CI(k>G@V zZ?nM0l3Ik=F3{{i_f5Fz#rt?=`Y-=>m?;+w2aTbJ0XGYFjS|qh(AR(I0G4t`8lwgh z$9!&V2_JhVr*6QoJ%V!&FYAJwK0HNpP~~gZR&ykDi<+MaYnn0kN7x;7aVS_2y?qS` zx++vKRHC`)HAXm>rH9UggFk8le$Ye)J1@h%;+ESwoke?Gv0u1>o2|(TeGI#rMxrLz zq<4DKH9nCj7%d4~aBoa}?w&PV!feq;WSY-StHcB@A+qX}wJ+lKw5*JT3#OM|H}gORmH5t1wQ1sRWUM1MsRI}T~Q~O_rmn6Jmo8K zDm_zVjGiOXLc^YIWI#sBFHNd(=7k1=^0PSN?l(R943=&Ys$~XlgU(fc;RsYUMqdwQ zO3{urLejApCEf%udSR^hjpsYFhZ?k*7^!s7Nh^81`i7NLXqJ^$|MhjHrro61+ zXLF6De%GK{^j@E#o9|#KoFk)1l-#so$!^lYdTJwH=UV!D*dQgd4Nr4qn+_vC;_ncb zH2+;LIGBkqb0H>ME0|vt&iZbE(p%?VvlTVM+o`26v85T@FkO$&6Q-rg-mHxI2n6)K z7u3~yKx^F6d`Y18OX$Pia~Vu)2jMYHDMP~sUiyac-YMug0HkM~Z@q~buyo8uqhv3q zM2*PZni!CRde`|h&e&i)rPl}lks)x&tg~4N-)%dj551Q)<(-79VR=WAqbk>^bEC=F zxipzpWfrCyB3(P9m}BRgdw21xwhc*TrLug3XS>#g|Z|8T+c5n}Ck*r1{G{=|lgGI66b*v*Ea8 zW>8t4+Aw@Ubpzer74t{kCNv!1L_gR)1vi*;n{Q5?dF{&hpSSg-vc5Oie$cE19YBm{ zpIoqQ{Fuf?N{5+IN+T z=&faproA2r#ZA|Mly$2Ki}`-i@%Q8or1*M?Y~qo2Y&;~ce42S~2hq1&l(KlT+*UGT z*f0}J(kmX|S+hDu#cn>1ntX)pX3UItw=Df|hGg(0#)9F6;xL|wG$@HSBfU^zHl0VD z&LcnWb|dzu4g-ugg-&eocYmXgcWsltcVyRZBw_AY;qJj7FbQyDcB>NoYd3;&^^QxT zKYbV9oA%vP1(vkv)0BR#eW>{ z)OU30B!IrEIo(eR-$q&#v$zs*^1c5iBM;5~k4oXhx2co2xkZ_>lOV{3zfkUT87YEA zWgdz*=bN9+f5psQOw%06k3lm0_-8z*;q&W@AndryxgazKfnpy4)FU2l&_wMGN|}ND zSQyV+N+66=?h{EI>x2759XFo_jCGf8uZ*d&k4QhYTHboZ z$)tfz_bPGPOw7{J70eYmF+@5t=n^Z{B+rEwr&fs+fOEPf{lVf%^zZx8zv=Gh#bnV) z=U=n;PC)2$SXpMTx`Ehqx~l;Lf2X=uf77$rUZmKhX)K=kX?0(G#5l?Q^)9F-YyJ7& zZPt2-TRJQ|Y&Ro)7;JbI2z=p!qmnAgS8i};tbdrR)p>ODg>DalU;%2Qe=QRdIgue~ z#+taS8&qjE(Qv>|^j^I25ak__sGHBpma%$g@L+A=+`mU>e= z^+O+J66F4|!V}#1am3K?%WvEMNg9N*jCT}wt5Mma2@lL*NyYp3rW)JU71fUx_7`t5 zz*xGM`~43hm(Ed&mO$MZn&=lAIv2GnuoPk*-M)MB`S|#_Vc-LSA)Yao#x)C1T}H|( zWDKAZ*K^D-AN!ukbuQoR+oiru>7NUfXxv~sj~a79|9M!yfD49Et-FkAm^oA6$j_SA z$QUYJKMs6;y=_j%$#)_|II`Vl8j9~LvePt4SbdgTp#(|=6My~!SqBf={6{__TfMu+cP`7D_I^u=(6J`$+~A4rV}D!ZimkLl;-VdfexLUs{!8*;;0BjnM9>V%3g_aRATA;mN{c5q ziQ|hWxp8mU16*k5vkXv+Yk4ZFJRjg_XZ=A~9o@v;-!0Hn-p}uw2>YkXBkO^!TsT3{ z33xT*JV>lK9Jhb+=Bbj;d;zHRxA#f^T>U9iDo**tD7lsoKf3 z{F-8c^!InlutU92zL2%X8C}d^KhywDzDFaU@YXz`oCJECo(@LMr3O8TK1OH75aifH zI!?|vq#x45f!8*^2$-&{myeU!zB*eI$ITVvx%t$+3=RE!m*i+>h-mM3Z#{H+Kz0q< z6`J{Ms{Jts5rmeOpmG!SZlQd4%RNffZuqTvp0-0q!cpUlYTBxtbdqGI(J6C53_lH4 zre-pgX4sJW`JrHjk_cIl?yYz}6&j-#=R%zkMlO@qRf~niL$<2SCht-P(D3_^I5R%h zD-|St7)lt&(wC>%_6^Bk5PWM=dATV2v-#1q=woZlTt;e7DT!(E1a&B2=^JYE1>4T3 zXl@V?3}vT?NZmw{NhQT}IMZM61Ml=fSJ6L}I%P7airjP*FQx&8!!bIsx3n&xJ%u=t z$WrBLXE;Bj+~wZ*DXRw_y7FUh=5F6?utiY>L17pJe{cW7!|qLhE&vb8?HnXzk@bD8->0OEi}3BNYkuWDT(@`z9ksgH#Q--j{{uX68vj-%H3E zGfb)(4=DL`R(L9G61v;06ll7OI?mM&Yd>K#udm>=*)R2A_p}Uys zn7OaA+#eBA^ch3vc3z z+T8{=Sj5?5MrfiOT?pA}o7|z~sm-?Y7bNXLL5C^WUA^uYI1>qE8w&L!M#(%~C( z^wtd9<%|AOPSONacTcODgo^ME@RUi4_2|6}$BovmfnljFzv34KM^T2lMGB9wL#e9+ z(|r1<190#Q3}xza7vYJsJX|w~_=FAbeWeG?lv(t8)!v#& z?2LWLF871KIz2j@i1ngL`yIoi<`>M@{48wdss5$HKqxOgffWyGqGk@7aXb;^ywX{;pNC- zA|XyOhepDFmkg;3+ts?@dVR~`St8mv7A8??*Q2R`h;;_Zd zyLx1*8u-G&WVs3@4;e6EcpGE{eW}jNGj7;&pKHRooUgyNk%YDUWHjkrt z8UEPp1Ghmhl`sQ^R#a_nTKcQW4EC66>DsbO^75qrjbHd0-GNx#aC8K7PAo8KO3{x+ zf!)X!GYgf~M;#zSwz{c#KjM7cb5|6@TV<27ct~Qmwz5ZouEp(p$t*oC25#y+x|j0( z69$f3O}FZ!!RbtPWvuOrNpHd42Rs?q;9OVr9#XkJ7TK=YsyDRk^>)caDqgGq&2K>( zuY4mX@nJNz!}I2lX80^MD(AnL*da}Y4>Xq%`pkfh5XP;l*OScN3{ z!F3uwa3MQM9F!%K5<6M)D}9v*48>RXT)$iV81Oc@XADiFnTZ6L8DX7jQq?)LOq*;e+OlT0F7(?B06_WlN!MS-k_pBo5GNH|{#g88 zdvJ(&@{Olsnf%Sn3t;zJl1qpFmDA<0ZW=I+%5VLy&{ldG!R6`s@^URDgi z`@=SpTEAIE1luQD?>hX`ar>9B_)+L&$Zn{hY(QlG)(f#VVfg+hn6s=%Ni?2rMjZdV z_HnAb1UWi$=Go#&b6hAttS5=K{C)D%dp!7#-bZgv)plQ5aZYAC`j4@`tG$K>O$3l1 zr70RnPg-XcWKS*_T^X7!a*1tG{?GGBk5gph_Brg^=420l3sEY&Y0^}6#?Ej51%C9b zJU1g4E_?6v%VATOVY~fGHSxn^e{&KZzgCgE^fGNC4O5XIEZW!s(hW{H;UL6N=%{Vm?Zy- zzk*9NT;yB}|5eq!t#nrv4FfR0_b)?~t=I+l5tGH#UAppjl9o=|mcsAbUwU;KrTr^X zc%uBRa#EBI>KbscS~y`p7*XEa(hgg0nxu9t-&rK>zlC4q0_wrL^O z3ylm=q@!@xK+76t(qOb@O}-%GWjC^0HcaK^oc3c!LbJ1_=NO;EANapOGTbACJr$L6 zA1kwHsjBEJx(gK%z7+6zu!r$&D-YLGXS4Ys!9K;j>`Si#7m}VOwbb0p<+2oK?NhyZ zb`h)&7zQIN6nkRB^?DhrBE-gRll(DaGL5uKo`;1$k|AV0=elFH6KFa63xdoaB$h0H z_am3?FZlA}h)`rFqlTq4*&Y?BKI;mfew^1Z@Mg8^)0Z}GbgFfZv^x>M9s8w>{`ne< z&PAaQ_onQ_DGu2x4=IkYb=9mH8|vwOT~s$rG&eEJi>1%u2KX%-;v&>Ky3ID`YUeRD ze zlQ!BIF>{VkyRn4iLNQqUt@)W|A8$(9k#xtmcsLlRyL+@Ym|gddf~`&F{`xIo7B>=;>$hE1YJe4VgW z5liAmMs#yR>tV6vkkL=LiSy^m__VuN&uO}Hf*9G|n!&G@+tBAhjjiCM>!qTnW#yz7 zEQZ7(mVRu9LJ8^l*VkY&E^t!W9_~hXYp>t0VuDLWNaftDoNE5cJ0XFFtITVES9-Kk zQ7T@;jE{)<6IWV3hK@56ES4Oi+H(#iB5pn{<7hiXtYN-_F^^qM`3LVf=b!ATgE__r zeDaKmrF zf)ACtrYU{$Y%ms?y>%X)5A73SOSn}&Fag; z(SsK>4Y&mV-a!)xYfz!-1TZe~Y4`9G9}^Wofht&6sYaz*_?t{>2=M+{etU1T?S8L| zKgqq9%1rVbNjbH{mpZEuRImf2=TMZC9XJn`1D*2 zyG(Y-RwyTR+J`sBg6hHZlb6vaFTG1rcg*Pu(=jj5HKoZRT8wr!RjOpG2%EyEc2*9o z7?%RwUb2UB4ZD@L#5#U1?p#Hw^1YXH0w$45{e^`Bnr@vV?vsDlq_mHkNE`(f7B^Xf)cE#fM)^{{<6T(BRZPT>Fucb6g$ z-_CZN{Vdbj1EUP;htmZQhLWacPMNQ59U0?U&jqo;gewF2LnY2}> zEqIlxNNZrtYnMLVpa?BAyiE%XH&(BhUcKxxYOL>dFub*Bf#sy2MeI2&`Rjsk#0GWV z`bPEt_U+*k!r;l_2og?oE;rW~7u%KTE!*cl^VuGy?rP*If6Vt)2qhmALpjG;;nuZ@aX)ITVXR7b& zn5N1YGh`DNNm+iQYU*$2EPj(4^#P2{D=<^CZq@X6Qu}nWgI(w0Wcy?Bo2ph9tx^%s z*^rZ?mllXc_-eM87Y4?#JgoTE`LpHoQIx&Cy#da%&9ZbMFnoSz?U5nJmHLe?L$cbH z==-MEP(z(nhK{OA@uKuavPcymrt*Urb556EPX>__E)h6;Sa?Tx8K8)3Xet>S@>e;Ko#{}8p7tDl zmbh%VJPT zgoIh~W4zf`5q`2_0IUK8d-8QkX#zXex1``!?2dZmr>^2cH6>1&EDSVYq@%2Xh`{k2DVQb(v@E zjg?;ARE5?W>`oh_&r2EB%ED0m{S76~I$)lqdppxBCe4Z1yE z8ENy4Lpcz+d9xnD&9)Mp4%N{o705lcQhwzOGu^l$=)6e4*m7s$Wht>9XcMAO?D3#T zGtot?V!ZJgeAZ1{%QzjO*!>_Y5YB7J^k}=lblBl?QTd`cf>*4jVV)ZBYp@+wq2kZW ztQ#&J?B#Y< zaI&L*2_t_6Gm|(dGY?1@hMS=rsj!_Jeu!%;7dS9IVDl3zB zTloUgaK`>yDXrhY=w~f@1xid0)IhOzZk2Ufavxsl;-8IL4t?toCZBm=7Wch#*mi-) zCYboB@&uWg9cLE~SZAn8ZMjP@F!sTr=Bscq99j<$cNo>-{R=4t~$l#T>oUw@TDkFqu+8wfxU$E| z_G{+DA8)q)@WNP+fZWGaQ5B4xF^~!*tFN(cs0u$>T3Oh?yXiQp1|wTvPi3T1Rjvtj z$DDz_3q2UW|9S61s3(!p`3+7!r?|a8WC@Rbc94#1q5@nxmQF}%x`>&Nxx^zwwnUDp zlvD}W!uh?#<(*>Ycu|+_lc=Y$_59tv4?D#r?7~N;Q&%DkHfP(?iDw=4HXZdj0Rz_aoWU-rz zA>W(=iI*z=m`)s}aD}Kb*iy8IH3`C;C81VL(R=T&m_Wq*n2sHS;VMu-Njr1PF2%adT4nqLI7a zZx0xX$|R}C$Pj5O>JAyT7~x(#4SYD!&Q!|9{s2Zh*Zh^yGUJBiX4H-tqTFL;8`AWV z*A>6bbiT{(y**~m^lD}W$?wfVz@xyY;`SmB4#I~N400?hQLQ3#%1;n5*9XFpXDy&4 zt}fB=HE|#G;2qumq#u9#;I^0T03b5bz3b%)Ir zpI1Ta!I)l(WN@f>5Y8YwOxJ2wUN`gFO+bF#>JN%T!DZA>-+l^PWob=Tf>?y_ssOjV z>yq!we@!NHL2i1R%CL4TaO^kmSY$op|6<$u7A89)q(?8XMZ%>g7GI}{aWv~@%3$u? zuULL1O!NLx7ci-`t~cU-ui%|?`51hNK3d+`&}1x?F|)<*yoc>~d2LY`ZR@{$6&0XL z`2_YepO^Ee*-~b|!2F+Cic}2rX>J?4s6kiOG|)SkI-n!f z2mCis3QHyQ*Qz2TRYoQ|1qF24>gI5{W%zuvk$Snui_$y+(N$IP zHO?Gqyq~(6XfIU~A9TP9_;FC+=6SpB;l8ohS}GPPeQ>2DQc3xFEwh@t!-bL*k`m1} z7`_vq0V|_)Cb+j6aV-uXr@f{p%)pfGN6(@?Pw^U?8(*X+no2@bDOn)CDpz}W)_OQ- zHxU$D020Tj4q=d(6<+OnG^o&1UYo2sugkpsh^o9;I?=jS$NyG|h?q)UV#Ha*N>*m$ zLOBa!%|#=qL5Tas5AL0E2b*khe-#AsTGzF2e?~ozS@K6+?x3I%_C38T_#k7J=d6b~^bxmR5Zxc8<%JBp@*ed(XfBeXGB(q6NI@l<@Q|#Ii zVco-OQBm_UtRwXF=j<+sazQH+^dBwUl?5MEzZXBn1T$Ub%u&akG+FHwc;l9ueI5kd z$gxG+Ukr$`6?@nTRG${W2HM@OP6E~*RZugU7~lC|wHv$sf|_gm@(ywYyiAJ(G>eut zCN45ld?ck`05@gt{ldQPHJGz*5jg<8@iQ{{5B`7pIP*s+yZ4XJlp^bh@Q4^YA!|(8 zcd}&48f7gE8iwp-m+U5Gm+WiQ5F+bH7=-cIJ;>Mxh3uc>`v-jIr~9Y5@B5tlT<2WZ z^?tqir;CSVM)o7?R`vc-mM~PT9=J#T`rB~DAT&lEa2L{a>X&c+Q7;-okhamAR&rx$ zZ$^?9LQ`7T#rRv^e`-&N<%Wrt2m(@CPB~Z`hs82|Wq{Ys%@Z`L#&W2d3nWHaO z;;eeYwv}%LXSiv<2;nFD?a;dge67{J5 zAT!LRN~Vpio!(ndT8?*m1v6L5^V#eHSbtIS(YgAVtJ@>bNhO;fBVLulhUJT}*^$sFfby05M5m|tc!CI8XV9dL2 ztucz@Wpx$+KEdzcW7Zkv)g|I1_iwX_Ko;IK?an>v4;%I1kmd-A56LZ;wq_rl)0fU9 zsdYge?*Dwqy&J?}i4wm|P?SV_hChgf6v}lO8>pzY*Vre0-9MMq{x@u-JUtXzCxu^D zKKi5_@NDNL=;b`mSfaVfyFL2`sW+De4!c#H=F8->Q*y}YzM`DfB0YuGaW=ogvn>u; z@dy#Jr=B@(HUI6GK#q%(2}&Vhu>o6is)5IUw}+zHuo|TFS1P)?q`EAli~N!T1raGp zu_j0kqrq@2yZE|9Ojdew5xxp`;oO_!S1ORGTcRAty!~s;gVElxq8 z%GThdTfXzkwe@G|OF!p?goo&Onx&7lALg?#gEDkT`1f&)s)v&uuf}COcA;L6XXn-B z+*><75Aw5`0n%~~b?R_bkF9~kmzF}&3phW-EsChHQ@80>Qo}?E1aC)ULeu-}763 z3Zp4iSCSHhad{hMzzyRsW@lQ2VXsEn0x_yb{bULRoTeHnjc+4l1yz@)D2nty!+1t%v;J%R?idS?NLof< zsz+cNW8Q~8!6csqy$FDEY;&D$b18>X=c;<>@ZserZ2AJD6$49Fx-tHMS)5-64k?PT z<&vjLr*}4T)3Y9J~LgZJkjmgj{nMIC(`oueg8U@5==^mCCY{*fJ zC{>di(3a6*cqqurTpTH|76aCTWezuw`%P~O-hzOVfWm80@N4gVqsjY64fH8H;M?UC z!HYX@1u9uhY6Z6R=SXnv8l%7D5n;wHr-y5-!`oC;CSOrc+btW7x+`+ymrS>Opej~Y zR#_4L=eU4LR&x5uGFioqux$AW6|lqOaQ~`qPjOfZ;KK~_#k4v z!>{mQ9#1ifvCz<=Gv4m=@ClCdIj9$c9>mob`3Cn3%ucL!e~%6*c))&qD#oUO%#w)` z$ldq$>T-iw%4F|1@Tz!-`}O;tI1Ll38gOZDI#r{b>B@&Klqj*8nVEL*c-{oz)^ui6 zQIYc(Dd_#{qWTVKf3osgI{N*UL`gmv1Oj5_t5GmJ<0&IISJc140+EnWzOmJE|LkT3 zAW$(Hjm4R!c3_{E{eZDxR|yZ}{d}pG^1}sea2D^H1Vm0Td?y+vH^TY$hm1~Uynyx1 zt)|k;6XbpMzT0Z)xvNc=mJojQ#%D^}!nwiQv){+>EXZHr!52B_#13USww9Nd1C+A& z%o0`cVe-w!?&s@FyA0))eS?~Ks@%ZD>bu2JX5}WV_6~!6LnMnIgWfdK#4=YAU6>@` z4w!3M)l|m38xJ`iLdWwD1^T*a83J$ZuIXEkW$UzVyjZj4If?0EO21n8+VNy+R9?HM z0OQe@nVA`+Xy&zz#D$GT804+&=XOPSS?8cfEV;6=K#%jKu`#2%pe&6G#zmw71ao)#aj!eFCN-f)VUQZRfW^RBXr($u>CgPNM0ai z%6n;hm_O}fix*a;bmh9F$V2n_3%$aND*3}a>D<6Gz(G&dp*#15C zlAI&n2#$35LKD{3_rpi)YFu}`(x_cV=2_qNI>t%*ux=dP0PF320)ne%cXZSw`K9$( z*V(t`^6~gdo%XI%+f?cNMV;HC@-dfTNe>KPxg!qvv$VHtxomIMJ+@;AzBC&|5Nn`7 z9HekTW1(%YKB#ICXve23HWfdY^9x2_}(Q50|I%^2lXnP`8&>;Igd?z#qGaG!g5%Qz90*j2qRDmNW zrX?od>A$Q;Xo4AWLh5_e5Kdn zkVny~0Opee0ZjM@$F)(1(BnCWxgEU2L?CtI88#DW&(JB+s5ZER(F)`6YDfcd%LJEIh z&lyJ~05K%2V5D@dk@&CBJriT*cnS_E1QJSfjfWj?{*SL+B#dZ@D#X^6i8N|=9t5_e zGb(2Ge-xswGZ&q5kWz(6IA=&rj>xm)+3JaUcs82pl*@ulR1=l|O|z)S{(Zz=8uxeh zC8p}n3<_`F!zLY?gttV2cks*Mh0M9&Jdt`?{ZH3y-{lK!t*^(0oF0(G6YQBZ8L^H% zZ+rP?XJ@N}&rW>A3msB>A`m2i&k_)S!Pa}a(QIZsyeY$?rehx$m!WoJQ;k7#_F;dt z1k;CNk*KjB`Cd*&ilI%loT$8EDni5lSgorZo*!z{FVNpcX)ztlCL~%$`DHWpCyr(i z-%_IWcF4`tXH%6YLWX8?z?mRJ-vwQfwTe!xi)yT^yf*F{G8cN_gLmV-=rG`&0oPK! z^7c~beYggvR2N`D1IXUYn(T}>@#ex2=8||*JeynPD9AD^mwq<@Ih%}uTbiRluhjm@ z{%rl&kbmC7g1rp>rDRl1(R?pAPw($68n{Ryvr}a$1tYi4v=m;{zR#DOns=ESUXmP@ zxo-ME;CPM^;x*Y?r?D*Lp3H%dQZei;Ng9OJ8geHrw7Jr-H%Tk~;?#N{``V4>)?hO4 z_=eWLurVtn@9p63Go2oiB=IN9?zoO!8nAxa{( zd(B-N>!^VvvKh*)a%h`oZ zxVWkVTIHV4e16_X!UTJI+7Si@+7R8v!@_{8EGeU2ZH(%1nEz4hl}*X?Xv{l{?^RWq z$Wwr>$uLI7C2^0zG*yiRm}iq*EswaaE5L=7G}U5oE()>XHF>tIEQv8oYJD{)N`Bd! zf2V6Yj9J%GO9JY2_l_7-VNDlb5?|vm`|+(t%sP#dN}-aMiI!B;#a6S1B{gjGv+*d+ z+ml=emWYagSPdR+ddA6!RG7QR=W>mB)BNA66^fnVl^B9NNuOC$OV7d*<($bvuUxn zWzII!0stTvynH9h^~B#M6}{P|d6liw=Le`I_Ag%KX##JhN)W&3>5lXtk_rHCT8-mhrKjX{-m~g_9v)j59C&zLY)sZUb%qHNcS5VMHaYjeRtAT60Z>S>`pe4 zUg*tJo8gjD;gp1Ys#LxjOH?mmRdJsR^h@nJ;yEtZz;=g!|=$6DK0uKkv1+TfUByv=U`8 zM8!u(M?Y1{1%h3odjl}P>?953UQ8wf6CyTQYMV6*3Ufg?1mf2rF(k`F^zH}8X+$z( zo^j@rN>kU65U$w#IY?!iHO)@M$x;mOsRdwW+Pv<#^}KWZ2k8xbLC&|>uU{jkr>E1r zC@Cq0B2;#KH56V7wCFRYi5~!eD=Mo0uDn5mnSerKViIpChJo)t<7}&g;dpVefPWtB z_0+&BM>N5n9-`7d0xplYw$KJqi;>%tX0?$`bz|mr-It93=Y!~K8fsLj*@XQM;@_Ua literal 0 HcmV?d00001 diff --git a/plugins/LensDistortionEstimator/resources/icons.svg b/plugins/LensDistortionEstimator/resources/icons.svg new file mode 100644 index 0000000000000..4a0d06ffe2e07 --- /dev/null +++ b/plugins/LensDistortionEstimator/resources/icons.svg @@ -0,0 +1,4423 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + image/svg+xml + + + + + + + + + + + + Made with Inkscape, Export with 480dpi + Notes: + Following guidelines explained on Stellarium/data/gui/icons.svg + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/plugins/LensDistortionEstimator/src/CMakeLists.txt b/plugins/LensDistortionEstimator/src/CMakeLists.txt new file mode 100644 index 0000000000000..a2f027d2b7e7c --- /dev/null +++ b/plugins/LensDistortionEstimator/src/CMakeLists.txt @@ -0,0 +1,73 @@ +FIND_PACKAGE(exiv2) + +CPMFindPackage(NAME NLopt + URL https://github.com/stevengj/nlopt/archive/refs/tags/v2.7.1.tar.gz + URL_HASH SHA256=db88232fa5cef0ff6e39943fc63ab6074208831dc0031cf1545f6ecd31ae2a1a + EXCLUDE_FROM_ALL yes + ) + +INCLUDE_DIRECTORIES( + . + gui + ${CMAKE_BINARY_DIR}/plugins/LensDistortionEstimator/src + ${CMAKE_BINARY_DIR}/plugins/LensDistortionEstimator/src/gui +) + +LINK_DIRECTORIES(${CMAKE_BINARY_DIR}/src) + +SET(LensDistortionEstimator_SRCS + LensDistortionEstimator.hpp + LensDistortionEstimator.cpp + gui/LensDistortionEstimatorDialog.hpp + gui/LensDistortionEstimatorDialog.cpp +) + +SET(LensDistortionEstimatorDialog_UIS + gui/lensDistortionEstimatorDialog.ui +) + +SET(LensDistortionEstimator_RES + ../resources/LensDistortionEstimator.qrc +) + +IF (${QT_VERSION_MAJOR} EQUAL "5") + QT5_WRAP_UI(LensDistortionEstimator_UIS_H ${LensDistortionEstimator_UIS}) + QT5_ADD_RESOURCES(LensDistortionEstimator_RES_CXX ${LensDistortionEstimator_RES}) +ELSE() + QT_WRAP_UI(LensDistortionEstimator_UIS_H ${LensDistortionEstimator_UIS}) + QT_ADD_RESOURCES(LensDistortionEstimator_RES_CXX ${LensDistortionEstimator_RES}) +ENDIF() + +ADD_LIBRARY(LensDistortionEstimator-static STATIC ${LensDistortionEstimator_SRCS} ${LensDistortionEstimator_RES_CXX} ${LensDistortionEstimator_UIS_H}) +TARGET_LINK_LIBRARIES(LensDistortionEstimator-static Qt${QT_VERSION_MAJOR}::Core Qt${QT_VERSION_MAJOR}::Widgets Qt${QT_VERSION_MAJOR}::OpenGL Qt${QT_VERSION_MAJOR}::Charts) +SET_TARGET_PROPERTIES(LensDistortionEstimator-static PROPERTIES OUTPUT_NAME "LensDistortionEstimator") +SET_TARGET_PROPERTIES(LensDistortionEstimator-static PROPERTIES COMPILE_FLAGS "-DQT_STATICPLUGIN") +if(exiv2_FOUND) + MESSAGE(STATUS "Found exiv2 ${exiv2_VERSION_MAJOR}.${exiv2_VERSION_MINOR}.${exiv2_VERSION_PATCH}") + TARGET_COMPILE_OPTIONS(LensDistortionEstimator-static PRIVATE "-DHAVE_EXIV2") + TARGET_LINK_LIBRARIES(LensDistortionEstimator-static exiv2lib) +endif() + +TARGET_COMPILE_OPTIONS(LensDistortionEstimator-static PRIVATE "-DHAVE_NLOPT") +IF(NLopt_ADDED) + ADD_DEPENDENCIES(LensDistortionEstimator-static generate-cpp) + TARGET_INCLUDE_DIRECTORIES(LensDistortionEstimator-static + PRIVATE ${CPM_PACKAGE_NLopt_BINARY_DIR} ${CPM_PACKAGE_NLopt_SOURCE_DIR}/src/api) + TARGET_LINK_LIBRARIES(LensDistortionEstimator-static nlopt) +ELSE() + # NLOPT_LIBRARIES variable is not set by CPM (unlike find_package; see + # https://github.com/cpm-cmake/CPM.cmake/issues/132), and the imported target + # may appear to be NLopt::nlopt or NLopt::nlopt_cxx, and using any one of these + # is incompatible with versions of nlopt where the other would be found. So + # let's check which one we should use. + IF(NLopt_DIR MATCHES "nlopt_cxx$") + set(NLOPT_LIBRARIES NLopt::nlopt_cxx) + ELSEIF(NLopt_DIR MATCHES "nlopt$") + set(NLOPT_LIBRARIES NLopt::nlopt) + ENDIF() + TARGET_LINK_LIBRARIES(LensDistortionEstimator-static ${NLOPT_LIBRARIES}) +ENDIF() + +ADD_DEPENDENCIES(AllStaticPlugins LensDistortionEstimator-static) + +SET_TARGET_PROPERTIES(LensDistortionEstimator-static PROPERTIES FOLDER "plugins/LensDistortionEstimator") diff --git a/plugins/LensDistortionEstimator/src/LensDistortionEstimator.cpp b/plugins/LensDistortionEstimator/src/LensDistortionEstimator.cpp new file mode 100644 index 0000000000000..b6ea06b23ac5c --- /dev/null +++ b/plugins/LensDistortionEstimator/src/LensDistortionEstimator.cpp @@ -0,0 +1,874 @@ +/* + * Copyright (C) 2023 Ruslan Kabatsayev + * + * This program is free software; you can redistribute it and/or + * modify it under the terms of the GNU General Public License + * as published by the Free Software Foundation; either version 2 + * of the License, or (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program; if not, write to the Free Software + * Foundation, Inc., 51 Franklin Street, Suite 500, Boston, MA 02110-1335, USA. + */ + +#include "LensDistortionEstimator.hpp" +#include "LensDistortionEstimatorDialog.hpp" +#include "StelApp.hpp" +#include "StelGui.hpp" +#include "StelCore.hpp" +#include "StelSRGB.hpp" +#include "StelUtils.hpp" +#include "StelPainter.hpp" +#include "StelFileMgr.hpp" +#include "StelGuiItems.hpp" +#include "StelModuleMgr.hpp" +#include "StelObjectMgr.hpp" +#include "StelTextureMgr.hpp" + +#include + +namespace +{ +constexpr float POINT_MARKER_RADIUS = 7; +constexpr int FS_QUAD_COORDS_PER_VERTEX = 2; +constexpr int SKY_VERTEX_ATTRIB_INDEX = 0; +constexpr char defaultImageAxesColor[] = "1,0.5,0"; +constexpr char defaultPointMarkerColor[] = "1,0,0"; +constexpr char defaultSelectedPointMarkerColor[] = "0,0.75,0"; +constexpr char defaultProjectionCenterMarkerColor[] = "0.25,0.5,0.5"; +} + +StelModule* LensDistortionEstimatorStelPluginInterface::getStelModule() const +{ + return new LensDistortionEstimator(); +} + +StelPluginInfo LensDistortionEstimatorStelPluginInterface::getPluginInfo() const +{ + // Allow to load the resources when used as a static plugin + Q_INIT_RESOURCE(LensDistortionEstimator); + + StelPluginInfo info; + info.id = "LensDistortionEstimator"; + info.displayedName = "Lens distortion estimator"; + info.authors = "Ruslan Kabatsayev"; + info.contact = STELLARIUM_DEV_URL; + info.description = "Estimates lens distortion by letting the user match stars to a photo and computing a fit."; + info.version = LENSDISTORTIONESTIMATOR_PLUGIN_VERSION; + info.license = LENSDISTORTIONESTIMATOR_PLUGIN_LICENSE; + return info; +} + +LensDistortionEstimator::LensDistortionEstimator() + : conf_(StelApp::getInstance().getSettings()) +{ + setObjectName("LensDistortionEstimator"); +} + +LensDistortionEstimator::~LensDistortionEstimator() +{ +} + +double LensDistortionEstimator::getCallOrder(StelModuleActionName actionName) const +{ + if(actionName==StelModule::ActionDraw) + return StelApp::getInstance().getModuleMgr().getModule("StarMgr")->getCallOrder(actionName)-1; + if(actionName==StelModule::ActionHandleMouseClicks) + return -11; + return 0; +} + +Vec2d LensDistortionEstimator::screenPointToImagePoint(const int x, const int y) const +{ + const auto prj = StelApp::getInstance().getCore()->getProjection(StelCore::FrameAltAz, StelCore::RefractionOff); + Vec3d point = Vec3d(0,0,0); + prj->unProject(x, y, point); + + const auto viewDir = point; + const auto imageCenterShift = dialog_->imageCenterShift(); + const auto imageSmallerSideFoV = M_PI/180 * dialog_->imageSmallerSideFoV(); + const auto projectionCenterDir = dialog_->projectionCenterDir(); + const auto imageUpDir = dialog_->imageUpDir(); + const double W = imageTex_->width(); + const double H = imageTex_->height(); + const double minSize = std::min(W,H); + const Vec2d normCoordToTexCoordCoefs = (minSize-1)*Vec2d(1./W,1./H); + + using namespace std; + const auto imageRightDir = normalize(projectionCenterDir ^ imageUpDir); + const auto R = acos(clamp(dot(viewDir, projectionCenterDir), -1., 1.)); + const auto angle = atan2(dot(viewDir, imageUpDir), + dot(viewDir, imageRightDir)); + const auto normalizedR = tan(R) / tan(imageSmallerSideFoV/2); + const auto distortedR = dialog_->applyDistortion(normalizedR); + const Vec2d centeredTexCoord = normCoordToTexCoordCoefs * (distortedR * Vec2d(cos(angle), sin(angle)) + imageCenterShift); + const Vec2d posInImage((0.5 + 0.5 * centeredTexCoord[0]) * W - 0.5, + (0.5 - 0.5 * centeredTexCoord[1]) * H - 0.5); + return posInImage; +} + +void LensDistortionEstimator::handleMouseClicks(QMouseEvent* event) +{ + if(!dialog_ || !dialog_->initialized() || !imageTex_) + { + event->setAccepted(false); + return; + } + +#if QT_VERSION >= QT_VERSION_CHECK(6, 0, 0) + const double x = event->position().x(), y = event->position().y(); +#else + const double x = event->x(), y = event->y(); +#endif + const auto posInImage = screenPointToImagePoint(x,y); + if(dialog_->isPickingAPoint()) + { + draggingImage_ = false; // Just in case + + if(event->type() == QEvent::MouseButtonPress) + { + event->setAccepted(true); + return; + } + if(event->type() != QEvent::MouseButtonRelease) + { + event->setAccepted(false); + return; + } + + const auto objs = StelApp::getInstance().getStelObjectMgr().getSelectedObject(); + if(objs.isEmpty()) + { + qWarning() << "Picked a point while nothing was selected... dropping the result"; + event->setAccepted(false); + return; + } + if(!objs.last()) + { + qWarning() << "NULL object is selected... dropping the result of picking"; + event->setAccepted(false); + return; + } + dialog_->resetImagePointPicking(); + dialog_->registerImagePoint(*objs[0], posInImage); + event->setAccepted(true); + return; + } + + if(draggingImage_ && event->type() == QEvent::MouseButtonRelease) + { + draggingImage_ = false; + imgPointToRotateAbout_ = screenPointToImagePoint(x,y); + event->setAccepted(true); + return; + } + + // If the modifiers are not what we use for dragging, or the button is wrong, we have nothing to do + const auto modifiers = qApp->keyboardModifiers() & (Qt::ShiftModifier|Qt::ControlModifier|Qt::AltModifier); + if(modifiers != (Qt::ShiftModifier|Qt::ControlModifier)) + { + event->setAccepted(false); + return; + } + + // To start dragging the image, the cursor must be inside + if(posInImage[0] < -0.5 || posInImage[0] > imageTex_->width() - 0.5 || + posInImage[1] < -0.5 || posInImage[1] > imageTex_->height() - 0.5) + { + event->setAccepted(false); + return; + } + + if(event->button() == Qt::LeftButton) + { + dragStartProjCenterElevation_ = dialog_->projectionCenterElevation(); + dragStartProjCenterAzimuth_ = dialog_->projectionCenterAzimuth(); + dragStartImgFieldRotation_ = dialog_->imageFieldRotation(); + dragStartPoint_ = QPointF(x,y); + draggingImage_ = true; + rotatingImage_ = false; + event->setAccepted(true); + } + else if(event->button() == Qt::RightButton) + { + dragStartImgSmallerSideFoV_ = dialog_->imageSmallerSideFoV(); + dragStartProjCenterElevation_ = dialog_->projectionCenterElevation(); + dragStartProjCenterAzimuth_ = dialog_->projectionCenterAzimuth(); + dragStartImgFieldRotation_ = dialog_->imageFieldRotation(); + dragStartPoint_ = QPointF(x,y); + draggingImage_ = false; + rotatingImage_ = true; + event->setAccepted(true); + } +} + +void LensDistortionEstimator::dragImage(const Vec2d& src, const Vec2d& dst) +{ + // These are used for computations, so initialize them to the starting values + dialog_->setProjectionCenterElevation(dragStartProjCenterElevation_); + dialog_->setProjectionCenterAzimuth(dragStartProjCenterAzimuth_); + dialog_->setImageFieldRotation(dragStartImgFieldRotation_); + + using namespace std; + // Trying to preserve the screen-up direction of the image + const auto prj = StelApp::getInstance().getCore()->getProjection(StelCore::FrameAltAz, StelCore::RefractionOff); + Vec3d srcDir1(0,0,0), srcDir2(0,0,0); + Vec3d dstDir1(0,0,0), dstDir2(0,0,0); + prj->unProject(src[0],src[1] , srcDir1); + prj->unProject(src[0],src[1]-1, srcDir2); + prj->unProject(dst[0],dst[1] , dstDir1); + prj->unProject(dst[0],dst[1]-1, dstDir2); + + auto crossDst = dstDir1 ^ dstDir2; + const auto crossSrc = srcDir1 ^ srcDir2; + + // Make sure angle between destination vectors is the same as between the source ones + const Mat3d rotToMatchAngle = Mat4d::rotation(normalize(crossDst), + asin(crossSrc.norm())-asin(crossDst.norm())).upper3x3(); + dstDir2 = rotToMatchAngle * dstDir2; + // Recompute cross product + crossDst = dstDir1 ^ dstDir2; + + // Find the matrix to rotate two points to the desired directions in space + const auto preRotSrc = Mat3d( srcDir1[0], srcDir1[1], srcDir1[2], + srcDir2[0], srcDir2[1], srcDir2[2], + crossSrc[0],crossSrc[1],crossSrc[2]); + const auto preRotDst = Mat3d( dstDir1[0], dstDir1[1], dstDir1[2], + dstDir2[0], dstDir2[1], dstDir2[2], + crossDst[0], crossDst[1], crossDst[2]); + const auto rotator = preRotDst * preRotSrc.inverse(); + + const auto origUpDir = dialog_->imageUpDir(); + const auto origCenterDir = dialog_->projectionCenterDir(); + const auto newCenterDir = normalize(rotator * origCenterDir); + const auto elevation = asin(clamp(newCenterDir[2], -1.,1.)); + const auto azimuth = atan2(newCenterDir[1], -newCenterDir[0]); + dialog_->setProjectionCenterAzimuth(180/M_PI * azimuth); + dialog_->setProjectionCenterElevation(180/M_PI * elevation); + + const auto origFromCenterToTop = normalize(origCenterDir + 1e-8 * origUpDir); + const auto newFromCenterToTop = normalize(rotator * origFromCenterToTop); + // Desired up direction + const auto upDirNew = normalize(newFromCenterToTop - newCenterDir); + // Renewed up direction takes into account the new center direction but not the desired orientation yet + dialog_->setImageFieldRotation(0); + const auto renewedUpDir = dialog_->imageUpDir(); + const auto upDirCross = renewedUpDir ^ upDirNew; + const auto upDirDot = dot(renewedUpDir, upDirNew); + const auto dirSign = dot(upDirCross, newCenterDir) > 0 ? -1. : 1.; + const auto upDirSinAngle = dirSign * (dirSign * upDirCross).norm(); + const auto upDirCosAngle = upDirDot; + const auto fieldRotation = atan2(upDirSinAngle, upDirCosAngle); + dialog_->setImageFieldRotation(180/M_PI * fieldRotation); +} + +void LensDistortionEstimator::rotateImage(const Vec2d& dragSrcPoint, const Vec2d& dragDstPoint) +{ + // These are used for computations, so initialize them to the starting values + dialog_->setImageSmallerSideFoV(dragStartImgSmallerSideFoV_); + dialog_->setProjectionCenterElevation(dragStartProjCenterElevation_); + dialog_->setProjectionCenterAzimuth(dragStartProjCenterAzimuth_); + dialog_->setImageFieldRotation(dragStartImgFieldRotation_); + + using namespace std; + + const auto srcPointInImg = screenPointToImagePoint(dragSrcPoint[0],dragSrcPoint[1]); + auto srcDir = dialog_->computeImagePointDir(srcPointInImg); + const auto desiredFixedDir = dialog_->computeImagePointDir(imgPointToRotateAbout_); + const auto prj = StelApp::getInstance().getCore()->getProjection(StelCore::FrameAltAz, StelCore::RefractionOff); + Vec3d dstDir(0,0,0); + prj->unProject(dragDstPoint[0],dragDstPoint[1], dstDir); + + // Find the correct FoV to make angle between the two image points correct + const auto desiredAngleBetweenDirs = acos(clamp(dot(desiredFixedDir, dstDir), -1., 1.)); + double fovMax = 0.999*M_PI; + double fovMin = 0; + double fov = (fovMin + fovMax) / 2.; + for(int n = 0; n < 53 && fovMin != fovMax; ++n) + { + const auto angleBetweenImagePoints = dialog_->computeAngleBetweenImagePoints(imgPointToRotateAbout_, + srcPointInImg, fov); + if(angleBetweenImagePoints > desiredAngleBetweenDirs) + fovMax = fov; + else + fovMin = fov; + fov = (fovMin + fovMax) / 2.; + } + dialog_->setImageSmallerSideFoV(180/M_PI * fov); + + // Recompute the directions after FoV update + const auto scaledImageFixedDir = dialog_->computeImagePointDir(imgPointToRotateAbout_); + srcDir = dialog_->computeImagePointDir(srcPointInImg); + + // Find the matrix to rotate two image points to corresponding directions in space + const auto crossDst = desiredFixedDir ^ dstDir; + const auto crossSrc = scaledImageFixedDir ^ srcDir; + const auto preRotSrc = Mat3d(scaledImageFixedDir[0], scaledImageFixedDir[1], scaledImageFixedDir[2], + srcDir[0], srcDir[1], srcDir[2], + crossSrc[0], crossSrc[1], crossSrc[2]); + const auto preRotDst = Mat3d(desiredFixedDir[0], desiredFixedDir[1], desiredFixedDir[2], + dstDir[0], dstDir[1], dstDir[2], + crossDst[0], crossDst[1], crossDst[2]); + const auto rotator = preRotDst * preRotSrc.inverse(); + + const auto origUpDir = dialog_->imageUpDir(); + const auto origCenterDir = dialog_->projectionCenterDir(); + const auto newCenterDir = normalize(rotator * origCenterDir); + const auto elevation = asin(clamp(newCenterDir[2], -1.,1.)); + const auto azimuth = atan2(newCenterDir[1], -newCenterDir[0]); + dialog_->setProjectionCenterAzimuth(180/M_PI * azimuth); + dialog_->setProjectionCenterElevation(180/M_PI * elevation); + + const auto origFromCenterToTop = normalize(origCenterDir + 1e-8 * origUpDir); + const auto newFromCenterToTop = normalize(rotator * origFromCenterToTop); + // Desired up direction + const auto upDirNew = normalize(newFromCenterToTop - newCenterDir); + // Renewed up direction takes into account the new center direction but not the desired orientation yet + dialog_->setImageFieldRotation(0); + const auto renewedUpDir = dialog_->imageUpDir(); + const auto upDirCross = renewedUpDir ^ upDirNew; + const auto upDirDot = dot(renewedUpDir, upDirNew); + const auto dirSign = dot(upDirCross, newCenterDir) > 0 ? -1. : 1.; + const auto upDirSinAngle = dirSign * (dirSign * upDirCross).norm(); + const auto upDirCosAngle = upDirDot; + const auto fieldRotation = atan2(upDirSinAngle, upDirCosAngle); + dialog_->setImageFieldRotation(180/M_PI * fieldRotation); +} + +bool LensDistortionEstimator::handleMouseMoves(const int x, const int y, const Qt::MouseButtons buttons) +{ + if(!dialog_ || !dialog_->initialized() || !imageTex_) return false; + if(draggingImage_ && (buttons & Qt::LeftButton)) + { + dragImage(Vec2d(dragStartPoint_.x(),dragStartPoint_.y()), Vec2d(x,y)); + return true; + } + if(rotatingImage_ && (buttons & Qt::RightButton)) + { + rotateImage(Vec2d(dragStartPoint_.x(),dragStartPoint_.y()), Vec2d(x,y)); + return true; + } + + draggingImage_ = false; + rotatingImage_ = false; + return false; +} + +void LensDistortionEstimator::setupCurrentVAO() +{ + auto& gl = *QOpenGLContext::currentContext()->functions(); + vbo_->bind(); + gl.glVertexAttribPointer(0, FS_QUAD_COORDS_PER_VERTEX, GL_FLOAT, false, 0, 0); + vbo_->release(); + gl.glEnableVertexAttribArray(SKY_VERTEX_ATTRIB_INDEX); +} + +void LensDistortionEstimator::bindVAO() +{ + if(vao_->isCreated()) + vao_->bind(); + else + setupCurrentVAO(); +} + +void LensDistortionEstimator::releaseVAO() +{ + if(vao_->isCreated()) + { + vao_->release(); + } + else + { + auto& gl = *QOpenGLContext::currentContext()->functions(); + gl.glDisableVertexAttribArray(SKY_VERTEX_ATTRIB_INDEX); + } +} + +void LensDistortionEstimator::init() +{ + dialog_.reset(new LensDistortionEstimatorDialog(this)); + + if(!conf_->childGroups().contains("LensDistortionEstimator")) + restoreDefaultSettings(); + loadSettings(); + + addAction("action_LensDistortionEstimator_togglePointPicking", N_("Lens Distortion Estimator"), + N_("Toggle pick-a-point mode"), this, + [this]{ dialog_->togglePointPickingMode(); }, "Ctrl+V"); + addAction("actionShow_LensDistortionEstimator_dialog", N_("Lens Distortion Estimator"), + N_("Show lens distortion estimator dialog"), "dialogVisible", "Ctrl+E"); + // Make sure that opening/closing the dialog by other means than dialogVisible property is reflected on this property + connect(dialog_.get(), &StelDialog::visibleChanged, this, [this](const bool visible) + { + if(visible != dialogVisible) + { + dialogVisible = visible; + emit dialogToggled(visible); + } + }); + + toolbarButton_ = new StelButton(nullptr, + QPixmap(":/LensDistortionEstimator/bt_LensDistortionEstimator_On.png"), + QPixmap(":/LensDistortionEstimator/bt_LensDistortionEstimator_Off.png"), + QPixmap(":/graphicGui/miscGlow32x32.png"), + "actionShow_LensDistortionEstimator_dialog", + false, + nullptr); + if(const auto gui = dynamic_cast(StelApp::getInstance().getGui())) + { + gui->getButtonBar()->addButton(toolbarButton_, "065-pluginsGroup"); + } + + auto& gl = *QOpenGLContext::currentContext()->functions(); + + vbo_.reset(new QOpenGLBuffer(QOpenGLBuffer::VertexBuffer)); + vbo_->create(); + vbo_->bind(); + const GLfloat vertices[]= + { + // full screen quad + -1, -1, + 1, -1, + -1, 1, + 1, 1, + }; + gl.glBufferData(GL_ARRAY_BUFFER, sizeof vertices, vertices, GL_STATIC_DRAW); + vbo_->release(); + + vao_.reset(new QOpenGLVertexArrayObject); + vao_->create(); + bindVAO(); + setupCurrentVAO(); + releaseVAO(); +} + +bool LensDistortionEstimator::configureGui(const bool show) +{ + if(show) + { + dialog_->setVisible(show); + dialog_->showSettingsPage(); + } + return true; +} + +QCursor LensDistortionEstimator::getPointPickCursor() const +{ + const int size = 32; + const QPointF center(size/2-0.5,size/2-0.5); + QPixmap pixmap(size,size); + pixmap.fill(Qt::transparent); + + QPainter painter(&pixmap); + + const float scale = StelApp::getInstance().getDevicePixelsPerPixel(); + const float radius = 0.5 + std::min((size-1)/2.f, POINT_MARKER_RADIUS * scale); + + painter.setRenderHint(QPainter::Antialiasing); + const auto lineWidth = scale * std::clamp(radius/7., 1., 2.5); + + for(float lineScale : {2.f, 1.f}) + { + painter.setPen(QPen(lineScale==1 ? Qt::black : Qt::white, lineWidth*lineScale)); + painter.drawEllipse(center, radius, radius); + + const float x = center.x(), y = center.y(); + const QPointF endpoints[] = {{x-radius+lineWidth/2, y}, + {x-radius*0.7f, y}, + {x+radius-lineWidth/2, y}, + {x+radius*0.7f, y}, + {x, y-radius+lineWidth/2}, + {x, y-radius*0.7f}, + {x, y+radius-lineWidth/2}, + {x, y+radius*0.7f}}; + painter.drawLines(endpoints, std::size(endpoints)/2); + } + + return QCursor(pixmap, center.x(), center.y()); +} + +void LensDistortionEstimator::renderPointMarker(StelPainter& sPainter, const float x, const float y, + float radius, const Vec3f& color) const +{ + const auto scale = sPainter.getProjector()->getDevicePixelsPerPixel(); + radius *= scale; + + sPainter.setBlending(true); + sPainter.setLineSmooth(true); + sPainter.setLineWidth(scale * std::clamp(radius/7, 1.f, 2.5f)); + sPainter.setColor(color); + + sPainter.drawCircle(x, y, radius); + + sPainter.enableClientStates(true); + const float vertexData[] = {x-radius , y, + x-radius*0.5f, y, + x+radius , y, + x+radius*0.5f, y, + x, y-radius , + x, y-radius*0.5f, + x, y+radius , + x, y+radius*0.5f}; + const auto vertCount = std::size(vertexData) / 2; + sPainter.setVertexPointer(2, GL_FLOAT, vertexData); + sPainter.drawFromArray(StelPainter::Lines, vertCount, 0, false); + sPainter.enableClientStates(false); +} + +void LensDistortionEstimator::renderProjectionCenterMarker(StelPainter& sPainter, const float x, const float y, + float radius, const Vec3f& color) const +{ + const auto scale = sPainter.getProjector()->getDevicePixelsPerPixel(); + radius *= scale; + + sPainter.setBlending(true); + sPainter.setLineSmooth(true); + sPainter.setLineWidth(scale * std::clamp(radius/7, 1.f, 2.5f)); + sPainter.setColor(color); + + sPainter.drawCircle(x, y, radius); + + sPainter.enableClientStates(true); + const auto size = radius / std::sqrt(2.f); + const float vertexData[] = {x-size , y-size , + x-size*0.1f, y-size*0.1f, + x+size , y+size , + x+size*0.1f, y+size*0.1f, + x+size , y-size , + x+size*0.1f, y-size*0.1f, + x-size , y+size , + x-size*0.1f, y+size*0.1f}; + const auto vertCount = std::size(vertexData) / 2; + sPainter.setVertexPointer(2, GL_FLOAT, vertexData); + sPainter.drawFromArray(StelPainter::Lines, vertCount, 0, false); + sPainter.enableClientStates(false); +} + +void LensDistortionEstimator::draw(StelCore* core) +{ + if(!imageEnabled) return; + + if(dialog_->imageChanged()) + { + const auto image = dialog_->image(); + if(image.isNull()) return; + + imageTex_.reset(new QOpenGLTexture(image)); + imageTex_->setMinMagFilters(QOpenGLTexture::LinearMipMapLinear, QOpenGLTexture::Linear); + imageTex_->setWrapMode(QOpenGLTexture::ClampToBorder); + imgPointToRotateAbout_ = Vec2d(imageTex_->width(),imageTex_->height()) / 2.; + dialog_->resetImageChangeFlag(); + } + if(!imageTex_) return; + + imageTex_->bind(0); + + const auto projector = core->getProjection(StelCore::FrameAltAz, StelCore::RefractionOff); + + if(!renderProgram_ || !prevProjector_ || !projector->isSameProjection(*prevProjector_)) + { + renderProgram_.reset(new QOpenGLShaderProgram); + prevProjector_ = projector; + + const auto vert = + StelOpenGL::globalShaderPrefix(StelOpenGL::VERTEX_SHADER) + + R"( +ATTRIBUTE highp vec3 vertex; +VARYING highp vec3 ndcPos; +void main() +{ + gl_Position = vec4(vertex, 1.); + ndcPos = vertex; +} +)"; + bool ok = renderProgram_->addShaderFromSourceCode(QOpenGLShader::Vertex, vert); + if(!renderProgram_->log().isEmpty()) + qWarning().noquote() << "LensDistortionEstimator: warnings while compiling vertex shader:\n" << renderProgram_->log(); + if(!ok) return; + + const auto frag = + StelOpenGL::globalShaderPrefix(StelOpenGL::FRAGMENT_SHADER) + + projector->getUnProjectShader() + + makeSRGBUtilsShader() + + R"( +VARYING highp vec3 ndcPos; +uniform sampler2D imageTex; +uniform mat4 projectionMatrixInverse; +uniform vec2 normCoordToTexCoordCoefs; +uniform float distNormalizationCoef; +uniform vec3 projectionCenterDir; +uniform vec3 imageUpDir; +uniform bool imageAxesEnabled; +uniform vec3 imageAxesColor; +uniform vec2 imageCenterShift; +uniform vec4 bgToSubtract; +uniform float brightness; +#define MODEL_POLY3 0 +#define MODEL_POLY5 1 +#define MODEL_PTLENS 2 +uniform int distortionModel; +uniform float distortionTerm1; +uniform float distortionTerm2; +uniform float distortionTerm3; +uniform float maxUndistortedR; +void main() +{ + const float PI = 3.14159265; + vec4 winPos = projectionMatrixInverse * vec4(ndcPos, 1); + bool unprojOK = false; + vec3 modelPos = unProject(winPos.x, winPos.y, unprojOK).xyz; + vec3 viewDir = normalize(modelPos); + + vec3 imageRightDir = normalize(cross(projectionCenterDir, imageUpDir)); + float viewProjCenterDot = dot(viewDir, projectionCenterDir); + float viewProjCenterCross = length(cross(viewDir, projectionCenterDir)); + float angle = atan(dot(viewDir, imageUpDir), dot(viewDir, imageRightDir)); + // Undistorted distance + float ru = distNormalizationCoef * tan(atan(viewProjCenterCross,viewProjCenterDot)); + + float ru2 = ru*ru; + float ru3 = ru2*ru; + float ru4 = ru2*ru2; + float k1 = distortionTerm1, k2 = distortionTerm2; + float a = distortionTerm1, b = distortionTerm2, c = distortionTerm3; + float distortedR = 0.; + if(distortionModel == MODEL_POLY3) + distortedR = ru*(1.-k1+k1*ru*ru); + else if(distortionModel == MODEL_POLY5) + distortedR = ru*(1+k1*ru2+k2*ru4); + else if(distortionModel == MODEL_PTLENS) + distortedR = ru*(a*ru3+b*ru2+c*ru+1-a-b-c); + + vec2 centeredTexCoord = normCoordToTexCoordCoefs * (distortedR * vec2(cos(angle), sin(angle)) + imageCenterShift); + + vec3 texColor = texture2D(imageTex, 0.5 + 0.5 * centeredTexCoord).rgb; + texColor = linearToSRGB(brightness * (srgbToLinear(texColor) - srgbToLinear(bgToSubtract.rgb)) + / + max(1. - srgbToLinear(bgToSubtract.rgb), vec3(1./255.))); + bool oppositeToImgDir = viewProjCenterDot < 0.; + FRAG_COLOR = oppositeToImgDir ? vec4(0) : vec4(texColor, 1); + + if(imageAxesEnabled && !oppositeToImgDir) + { + vec2 sGrad = vec2(dFdx(centeredTexCoord.s), dFdy(centeredTexCoord.s)); + vec2 tGrad = vec2(dFdx(centeredTexCoord.t), dFdy(centeredTexCoord.t)); + if(abs(centeredTexCoord.s) <= 1. && abs(centeredTexCoord.t) <= 1. && + (centeredTexCoord.t*centeredTexCoord.t < dot(tGrad,tGrad) || + centeredTexCoord.s*centeredTexCoord.s < dot(sGrad,sGrad))) + { + FRAG_COLOR = vec4(imageAxesColor,1); + } + } + + if(!unprojOK || ru > maxUndistortedR) + { + FRAG_COLOR = vec4(0); + } +} +)"; + ok = renderProgram_->addShaderFromSourceCode(QOpenGLShader::Fragment, frag); + if(!renderProgram_->log().isEmpty()) + qWarning().noquote() << "LensDistortionEstimator: warnings while compiling fragment shader:\n" << renderProgram_->log(); + + if(!ok) return; + + renderProgram_->bindAttributeLocation("vertex", SKY_VERTEX_ATTRIB_INDEX); + + if(!StelPainter::linkProg(renderProgram_.get(), "LensDistortionEstimator image render program")) + return; + + renderProgram_->bind(); + shaderVars.imageTex = renderProgram_->uniformLocation("imageTex"); + shaderVars.imageUpDir = renderProgram_->uniformLocation("imageUpDir"); + shaderVars.projectionCenterDir = renderProgram_->uniformLocation("projectionCenterDir"); + shaderVars.distNormalizationCoef = renderProgram_->uniformLocation("distNormalizationCoef"); + shaderVars.normCoordToTexCoordCoefs = renderProgram_->uniformLocation("normCoordToTexCoordCoefs"); + shaderVars.imageAxesEnabled = renderProgram_->uniformLocation("imageAxesEnabled"); + shaderVars.imageAxesColor = renderProgram_->uniformLocation("imageAxesColor"); + shaderVars.imageCenterShift = renderProgram_->uniformLocation("imageCenterShift"); + shaderVars.bgToSubtract = renderProgram_->uniformLocation("bgToSubtract"); + shaderVars.imageBrightness = renderProgram_->uniformLocation("brightness"); + shaderVars.distortionModel = renderProgram_->uniformLocation("distortionModel"); + shaderVars.distortionTerm1 = renderProgram_->uniformLocation("distortionTerm1"); + shaderVars.distortionTerm2 = renderProgram_->uniformLocation("distortionTerm2"); + shaderVars.distortionTerm3 = renderProgram_->uniformLocation("distortionTerm3"); + shaderVars.maxUndistortedR = renderProgram_->uniformLocation("maxUndistortedR"); + shaderVars.projectionMatrixInverse = renderProgram_->uniformLocation("projectionMatrixInverse"); + renderProgram_->release(); + } + if(!renderProgram_ || !renderProgram_->isLinked()) + return; + + renderProgram_->bind(); + + const int imageTexSampler = 0; + imageTex_->bind(imageTexSampler); + renderProgram_->setUniformValue(shaderVars.imageTex, imageTexSampler); + renderProgram_->setUniformValue(shaderVars.imageUpDir, dialog_->imageUpDir().toQVector()); + renderProgram_->setUniformValue(shaderVars.projectionCenterDir, dialog_->projectionCenterDir().toQVector()); + const float W = imageTex_->width(); + const float H = imageTex_->height(); + const float minSize = std::min(W,H); + const double imageSmallerSideFoV = M_PI/180 * dialog_->imageSmallerSideFoV(); + renderProgram_->setUniformValue(shaderVars.distNormalizationCoef, + float(1. / (std::tan(imageSmallerSideFoV / 2.)))); + renderProgram_->setUniformValue(shaderVars.imageAxesEnabled, imageAxesEnabled); + renderProgram_->setUniformValue(shaderVars.imageAxesColor, imageAxesColor.toQVector()); + renderProgram_->setUniformValue(shaderVars.imageCenterShift, dialog_->imageCenterShift().toQVector()); + renderProgram_->setUniformValue(shaderVars.bgToSubtract, dialog_->bgToSubtract()); + renderProgram_->setUniformValue(shaderVars.imageBrightness, float(dialog_->imageBrightness())); + renderProgram_->setUniformValue(shaderVars.distortionModel, int(dialog_->distortionModel())); + renderProgram_->setUniformValue(shaderVars.normCoordToTexCoordCoefs, (minSize-1)*QVector2D(1./W,1./H)); + renderProgram_->setUniformValue(shaderVars.distortionTerm1, float(dialog_->distortionTerm1())); + renderProgram_->setUniformValue(shaderVars.distortionTerm2, float(dialog_->distortionTerm2())); + renderProgram_->setUniformValue(shaderVars.distortionTerm3, float(dialog_->distortionTerm3())); + renderProgram_->setUniformValue(shaderVars.maxUndistortedR, float(dialog_->maxUndistortedR())); + renderProgram_->setUniformValue(shaderVars.projectionMatrixInverse, projector->getProjectionMatrix().toQMatrix().inverted()); + projector->setUnProjectUniforms(*renderProgram_); + + auto& gl = *QOpenGLContext::currentContext()->functions(); + + gl.glEnable(GL_BLEND); + gl.glBlendFunc(GL_ONE,GL_ONE); // allow colored sky background + bindVAO(); + gl.glDrawArrays(GL_TRIANGLE_STRIP, 0, 4); + releaseVAO(); + gl.glDisable(GL_BLEND); + + if(!pointMarkersEnabled && !projectionCenterMarkerEnabled) + return; + + StelPainter sPainter(core->getProjection(StelCore::FrameAltAz, StelCore::RefractionOff)); + if(pointMarkersEnabled) + { + for(const auto& status : dialog_->imagePointDirections()) + { + const auto color = status.selected ? selectedPointMarkerColor : pointMarkerColor; + Vec3d win; + if(sPainter.getProjector()->project(status.direction, win)) + renderPointMarker(sPainter, win[0], win[1], POINT_MARKER_RADIUS, color); + } + } + + if(projectionCenterMarkerEnabled) + { + const auto dir = dialog_->projectionCenterDir(); + Vec3d win; + if(sPainter.getProjector()->project(dir, win)) + renderProjectionCenterMarker(sPainter, win[0], win[1], 12, projectionCenterMarkerColor); + } +} + +void LensDistortionEstimator::setImageAxesColor(const Vec3f& color) +{ + if(imageAxesColor == color) return; + + imageAxesColor = color; + conf_->setValue("LensDistortionEstimator/image_axes_color", color.toStr()); + emit imageAxesColorChanged(color); +} + +void LensDistortionEstimator::setPointMarkerColor(const Vec3f& color) +{ + if(pointMarkerColor == color) return; + + pointMarkerColor = color; + conf_->setValue("LensDistortionEstimator/point_marker_color", color.toStr()); + emit pointMarkerColorChanged(color); +} + +void LensDistortionEstimator::setSelectedPointMarkerColor(const Vec3f& color) +{ + if(selectedPointMarkerColor == color) return; + + selectedPointMarkerColor = color; + conf_->setValue("LensDistortionEstimator/selected_point_marker_color", color.toStr()); + emit selectedPointMarkerColorChanged(color); +} + +void LensDistortionEstimator::setProjectionCenterMarkerColor(const Vec3f& color) +{ + if(projectionCenterMarkerColor == color) return; + + projectionCenterMarkerColor = color; + conf_->setValue("LensDistortionEstimator/center_of_projection_marker_color", color.toStr()); + emit projectionCenterMarkerColorChanged(color); +} + +void LensDistortionEstimator::setPointMarkersEnabled(const bool enable) +{ + if(pointMarkersEnabled == enable) return; + + pointMarkersEnabled = enable; + conf_->setValue("LensDistortionEstimator/mark_picked_points", enable); + emit pointMarkersToggled(enable); +} + +void LensDistortionEstimator::setProjectionCenterMarkerEnabled(const bool enable) +{ + if(projectionCenterMarkerEnabled == enable) return; + + projectionCenterMarkerEnabled = enable; + conf_->setValue("LensDistortionEstimator/mark_center_of_projection", enable); + emit pointMarkersToggled(enable); +} + +void LensDistortionEstimator::setImageAxesEnabled(const bool enable) +{ + if(imageAxesEnabled == enable) return; + + imageAxesEnabled = enable; + conf_->setValue("LensDistortionEstimator/show_image_axes", enable); + emit imageAxesToggled(enable); +} + +void LensDistortionEstimator::setImageEnabled(const bool enable) +{ + if(imageEnabled == enable) return; + + imageEnabled = enable; + conf_->setValue("LensDistortionEstimator/show_image", enable); + emit imageToggled(enable); +} + +void LensDistortionEstimator::setDialogVisible(const bool enable) +{ + if(dialogVisible == enable) return; + + dialogVisible = enable; + dialog_->setVisible(enable); + if(enable) + dialog_->showComputationPage(); +} + +void LensDistortionEstimator::loadSettings() +{ + setImageEnabled(conf_->value("LensDistortionEstimator/show_image", true).toBool()); + setImageAxesEnabled(conf_->value("LensDistortionEstimator/show_image_axes", false).toBool()); + setPointMarkersEnabled(conf_->value("LensDistortionEstimator/mark_picked_points", true).toBool()); + setProjectionCenterMarkerEnabled(conf_->value("LensDistortionEstimator/mark_center_of_projection", false).toBool()); + imageAxesColor = Vec3f(conf_->value("LensDistortionEstimator/image_axes_color", defaultImageAxesColor).toString()); + pointMarkerColor = Vec3f(conf_->value("LensDistortionEstimator/point_marker_color", defaultPointMarkerColor).toString()); + selectedPointMarkerColor = Vec3f(conf_->value("LensDistortionEstimator/selected_point_marker_color", defaultSelectedPointMarkerColor).toString()); + projectionCenterMarkerColor = Vec3f(conf_->value("LensDistortionEstimator/center_of_projection_marker_color", defaultProjectionCenterMarkerColor).toString()); +} + +void LensDistortionEstimator::restoreDefaultSettings() +{ + // Remove the old values... + conf_->remove("LensDistortionEstimator"); + // ...load the default values... + loadSettings(); + // But this doesn't save the colors, so: + conf_->beginGroup("LensDistortionEstimator"); + conf_->setValue("image_axes_color", defaultImageAxesColor); + conf_->setValue("point_marker_color", defaultPointMarkerColor); + conf_->setValue("selected_point_marker_color", defaultSelectedPointMarkerColor); + conf_->setValue("center_of_projection_marker_color", defaultProjectionCenterMarkerColor); + conf_->endGroup(); +} diff --git a/plugins/LensDistortionEstimator/src/LensDistortionEstimator.hpp b/plugins/LensDistortionEstimator/src/LensDistortionEstimator.hpp new file mode 100644 index 0000000000000..4530068474b41 --- /dev/null +++ b/plugins/LensDistortionEstimator/src/LensDistortionEstimator.hpp @@ -0,0 +1,167 @@ +/* + * Copyright (C) 2023 Ruslan Kabatsayev + * + * This program is free software; you can redistribute it and/or + * modify it under the terms of the GNU General Public License + * as published by the Free Software Foundation; either version 2 + * of the License, or (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program; if not, write to the Free Software + * Foundation, Inc., 51 Franklin Street, Suite 500, Boston, MA 02110-1335, USA. + */ + +#ifndef LENSDISTORTIONESTIMATOR_HPP +#define LENSDISTORTIONESTIMATOR_HPP + +#include +#include "StelCore.hpp" +#include "StelModule.hpp" +#include +#include +#include +#include + +class QSettings; +class StelButton; +class QMouseEvent; +class LensDistortionEstimatorDialog; +class LensDistortionEstimator : public StelModule +{ + Q_OBJECT + Q_PROPERTY(bool dialogVisible READ isDialogVisible WRITE setDialogVisible NOTIFY dialogToggled) + Q_PROPERTY(bool imageEnabled READ isImageEnabled WRITE setImageEnabled NOTIFY imageToggled) + Q_PROPERTY(bool imageAxesEnabled READ areImageAxesEnabled WRITE setImageAxesEnabled NOTIFY imageAxesToggled) + Q_PROPERTY(bool pointMarkersEnabled READ arePointMarkersEnabled WRITE setPointMarkersEnabled NOTIFY pointMarkersToggled) + Q_PROPERTY(bool projectionCenterMarkerEnabled READ isProjectionCenterMarkerEnabled WRITE setProjectionCenterMarkerEnabled NOTIFY projectionCenterMarkerToggled) + Q_PROPERTY(Vec3f selectedPointMarkerColor READ getSelectedPointMarkerColor WRITE setSelectedPointMarkerColor NOTIFY selectedPointMarkerColorChanged ) + Q_PROPERTY(Vec3f pointMarkerColor READ getPointMarkerColor WRITE setPointMarkerColor NOTIFY pointMarkerColorChanged ) + Q_PROPERTY(Vec3f projectionCenterMarkerColor READ getProjectionCenterMarkerColor WRITE setProjectionCenterMarkerColor NOTIFY projectionCenterMarkerColorChanged ) + Q_PROPERTY(Vec3f imageAxesColor READ getImageAxesColor WRITE setImageAxesColor NOTIFY imageAxesColorChanged ) +public: + enum class DistortionModel + { + Poly3, + Poly5, + PTLens, + }; + + LensDistortionEstimator(); + ~LensDistortionEstimator() override; + + void init() override; + void draw(StelCore* core) override; + double getCallOrder(StelModuleActionName actionName) const override; + void handleMouseClicks(QMouseEvent* event) override; + bool handleMouseMoves(int x, int y, Qt::MouseButtons) override; + bool configureGui(bool show = true) override; + bool arePointMarkersEnabled() const { return pointMarkersEnabled; } + bool isProjectionCenterMarkerEnabled() const { return projectionCenterMarkerEnabled; } + bool areImageAxesEnabled() const { return imageAxesEnabled; } + bool isImageEnabled() const { return imageEnabled; } + bool isDialogVisible() const { return dialogVisible; } + QCursor getPointPickCursor() const; + Vec3f getProjectionCenterMarkerColor() const { return projectionCenterMarkerColor; } + Vec3f getSelectedPointMarkerColor() const { return selectedPointMarkerColor; } + Vec3f getPointMarkerColor() const { return pointMarkerColor; } + Vec3f getImageAxesColor() const { return imageAxesColor; } + void restoreDefaultSettings(); + +public slots: + void setImageEnabled(bool enable); + void setDialogVisible(bool enable); + void setImageAxesEnabled(bool enabled); + void setPointMarkersEnabled(bool enabled); + void setProjectionCenterMarkerEnabled(bool enabled); + void setProjectionCenterMarkerColor(const Vec3f& color); + void setSelectedPointMarkerColor(const Vec3f& color); + void setPointMarkerColor(const Vec3f& color); + void setImageAxesColor(const Vec3f& color); +signals: + void imageToggled(bool enabled); + void dialogToggled(bool enabled); + void imageAxesToggled(bool enabled); + void pointMarkersToggled(bool enabled); + void projectionCenterMarkerToggled(bool enabled); + void projectionCenterMarkerColorChanged(const Vec3f& color); + void selectedPointMarkerColorChanged(const Vec3f& color); + void pointMarkerColorChanged(const Vec3f& color); + void imageAxesColorChanged(const Vec3f& color); +private: + void loadSettings(); + void setupCurrentVAO(); + void bindVAO(); + void releaseVAO(); + void renderPointMarker(StelPainter& sPainter, float x, float y, float radius, const Vec3f& color) const; + void renderProjectionCenterMarker(StelPainter& sPainter, float x, float y, float radius, const Vec3f& color) const; + Vec2d screenPointToImagePoint(int x, int y) const; + void dragImage(const Vec2d& src, const Vec2d& dst); + void rotateImage(const Vec2d& dragSrcPoint, const Vec2d& dragDstPoint); +private: + struct + { + int imageTex; + int imageUpDir; + int projectionCenterDir; + int distNormalizationCoef; + int normCoordToTexCoordCoefs; + int imageAxesEnabled; + int imageAxesColor; + int imageCenterShift; + int bgToSubtract; + int imageBrightness; + int distortionModel; + int distortionTerm1; + int distortionTerm2; + int distortionTerm3; + int maxUndistortedR; + int projectionMatrixInverse; + } shaderVars; + + QSettings* conf_; + StelButton* toolbarButton_ = nullptr; + StelProjectorP prevProjector_; + std::unique_ptr vbo_; + std::unique_ptr vao_; + std::unique_ptr renderProgram_; + std::unique_ptr dialog_; + std::unique_ptr imageTex_; + QPointF dragStartPoint_{0,0}; + Vec2d imgPointToRotateAbout_{0,0}; + double dragStartProjCenterElevation_ = 0; + double dragStartProjCenterAzimuth_ = 0; + double dragStartImgFieldRotation_ = 0; + double dragStartImgSmallerSideFoV_ = 0; + bool draggingImage_ = false; + bool rotatingImage_ = false; + bool pointMarkersEnabled = true; + bool projectionCenterMarkerEnabled = false; + bool imageAxesEnabled = false; + bool imageEnabled = true; + bool dialogVisible = false; + Vec3f imageAxesColor = Vec3f(0,0,0); + Vec3f pointMarkerColor = Vec3f(0,0,0); + Vec3f selectedPointMarkerColor = Vec3f(0,0,0); + Vec3f projectionCenterMarkerColor = Vec3f(0,0,0); +}; + + +#include +#include "StelPluginInterface.hpp" + +class LensDistortionEstimatorStelPluginInterface : public QObject, public StelPluginInterface +{ + Q_OBJECT + Q_PLUGIN_METADATA(IID StelPluginInterface_iid) + Q_INTERFACES(StelPluginInterface) +public: + StelModule* getStelModule() const override; + StelPluginInfo getPluginInfo() const override; +}; + +#endif diff --git a/plugins/LensDistortionEstimator/src/gui/LensDistortionEstimatorDialog.cpp b/plugins/LensDistortionEstimator/src/gui/LensDistortionEstimatorDialog.cpp new file mode 100644 index 0000000000000..3684a77fe64ef --- /dev/null +++ b/plugins/LensDistortionEstimator/src/gui/LensDistortionEstimatorDialog.cpp @@ -0,0 +1,1822 @@ +/* + * Stellarium: Lens distortion estimator plugin + * Copyright (C) 2023 Ruslan Kabatsayev + * + * This program is free software; you can redistribute it and/or + * modify it under the terms of the GNU General Public License + * as published by the Free Software Foundation; either version 2 + * of the License, or (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program; if not, write to the Free Software + * Foundation, Inc., 51 Franklin Street, Suite 500, Boston, MA 02110-1335, USA. +*/ + +#include "LensDistortionEstimatorDialog.hpp" +#include "LensDistortionEstimator.hpp" +#include "StarMgr.hpp" +#include "StelApp.hpp" +#include "StelGui.hpp" +#include "StelCore.hpp" +#include "StelUtils.hpp" +#include "StelMainView.hpp" +#include "StelModuleMgr.hpp" +#include "StelObjectMgr.hpp" +#include "StelMovementMgr.hpp" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#ifdef HAVE_EXIV2 +# include +#endif +#ifdef HAVE_NLOPT +# include +#endif +#if QT_VERSION < QT_VERSION_CHECK(6,0,0) +using namespace QtCharts; +#endif +#include "ui_lensDistortionEstimatorDialog.h" + +namespace +{ +struct Column +{ +enum +{ + ObjectName, + x, + y, +}; +}; +struct Role +{ +enum +{ + ObjectEnglishName = Qt::UserRole, +}; +}; +struct Page +{ +enum +{ + ImageLoading, + Fitting, + Settings, + About, +}; +}; +template auto sqr(T x) { return x*x; } +template auto cube(T x) { return x*x*x; } +using DistortionModel = LensDistortionEstimator::DistortionModel; + +//! Finds real roots of the polynomial ax^2+bx+c, sorted in ascending order +std::vector solveQuadratic(const double a, const double b, const double c) +{ + if(a==0 && b==0) return {}; + const auto D = b*b - 4*a*c; + if(D < 0) return {}; + const auto sqrtD = std::sqrt(D); + const auto signB = b < 0 ? -1. : 1.; + const auto u = -b - signB * sqrtD; + const auto x1 = u / (2*a); + const auto x2 = 2*c / u; + if(!std::isfinite(x1)) return {x2}; // a~0, the equation is linear + return x1 < x2 ? std::vector{x1, x2} : std::vector{x2, x1}; +} + +// Locates the root t in [tMin, tMax] of the equation t^3+Bt^2+Ct-1=0. The root must be unique. +double locateRootInCubicPoly(const double B, const double C, double tMin, double tMax) +{ + for(int n = 0; n < 100 && tMin != tMax; ++n) + { + const auto t = (tMin + tMax) / 2; + const auto poly = -1+t*(C+t*(B+t)); + if(poly > 0) + tMax = t; + else + tMin = t; + } + return (tMin + tMax) / 2; +} + +//! Finds real roots of the polynomial Ax^3+ax^2+bx+c, sorted in ascending order +std::vector solveCubic(double A, double a, double b, double c) +{ + // The algorithm here is intended to handle the hard cases like + // A~1e-17, a~1, b~1, c~1 well, without introducing arbitrary epsilons. + + using namespace std; + + if(abs(A*A*A) < std::numeric_limits::min()) + return solveQuadratic(a,b,c); + + if(c == 0) + { + // Easy special case: x*(Ax^2+ax+b)=0. It must be handled to + // avoid division by zero in the subsequent general solution. + auto roots = solveQuadratic(A,a,b); + roots.push_back(0); + sort(roots.begin(), roots.end()); + return roots; + } + + // Bring the equation to the form x^3+ax^2+bx+c=0 + a /= A; + b /= A; + c /= A; + A = 1; + + // Let's bracket the first root following https://www.johndcook.com/blog/2022/04/05/cubic/ + const auto t2x = -cbrt(c); // Changing variable x = t2x*t + const auto B = -a/cbrt(c); + const auto C = b/sqr(cbrt(c)); + // New equation: f(t)=0, i.e. t^3+Bt^2+Ct-1=0. + // f(0) = -1 + const auto f1 = B+C; // f(1) + double t; + if(f1==0) + { + t = 1; + } + else if(f1 > 0) + { + // We have a root t in (0,1) + t = locateRootInCubicPoly(B,C, 0,1); + } + else + { + // We have a root t in (1,inf). Change variable t=1/p: + // p^3-Cp^2-Bp-1=0. + const auto p = locateRootInCubicPoly(-C,-B, 0,1); + t = 1/p; + } + const auto x1 = t2x*t; + + // Long polynomial division of the LHS of the equation by (x-x1) allows us to reach this decomposition: + // (x - x1) (A x^2 + (A x1 + a) x + (A x1^2 + a x1 + b)) + (A x1^3 + a x1^2 + b x1 + c) = A x^3 + a x^2 + b x + c + // Now, assuming that x1 is a root, we can drop the remainder (A x1^3 + a x1^2 + b x1 + c), + // since it must be zero, so we get a quadratic equation for the remaining two roots: + // A x^2 + (A x1 + a) x + (A x1^2 + a x1 + b) = 0 + auto roots = solveQuadratic(A, A*x1 + a, A*x1*x1 + a*x1 + b); + roots.push_back(x1); + sort(roots.begin(), roots.end()); + return roots; +} + +} + +LensDistortionEstimatorDialog::LensDistortionEstimatorDialog(LensDistortionEstimator* lde) + : StelDialog("LensDistortionEstimator") + , lde_(lde) + , ui_(new Ui_LDEDialog{}) + , errorsChart_(new QChart) +{ + warningCheckTimer_.setInterval(1000); + connect(&warningCheckTimer_, &QTimer::timeout, this, &LensDistortionEstimatorDialog::periodicWarningsCheck); + + // Same as in AstroCalcChart + errorsChart_->setBackgroundBrush(QBrush(QColor(86, 87, 90))); + errorsChart_->setTitleBrush(QBrush(Qt::white)); + errorsChart_->setMargins(QMargins(2, 1, 2, 1)); // set to 0/0/0/0 for max space usage. This is between the title/axis labels and the enclosing QChartView. + errorsChart_->layout()->setContentsMargins(0, 0, 0, 0); + errorsChart_->setBackgroundRoundness(0); // remove rounded corners +} + +LensDistortionEstimatorDialog::~LensDistortionEstimatorDialog() +{ +} + +void LensDistortionEstimatorDialog::retranslate() +{ + if(dialog) + { + ui_->retranslateUi(dialog); + setAboutText(); + } +} + +void LensDistortionEstimatorDialog::createDialogContent() +{ + ui_->setupUi(dialog); + ui_->tabs->setCurrentIndex(0); + + // Kinetic scrolling + kineticScrollingList << ui_->about; + StelGui* gui= dynamic_cast(StelApp::getInstance().getGui()); + if(gui) + { + enableKineticScrolling(gui->getFlagUseKineticScrolling()); + connect(gui, &StelGui::flagUseKineticScrollingChanged, this, &LensDistortionEstimatorDialog::enableKineticScrolling); + } + + connect(&StelApp::getInstance(), &StelApp::languageChanged, this, &LensDistortionEstimatorDialog::retranslate); + connect(ui_->titleBar, &TitleBar::closeClicked, this, &StelDialog::close); + connect(ui_->titleBar, &TitleBar::movedTo, this, &LensDistortionEstimatorDialog::handleMovedTo); + + core_ = StelApp::getInstance().getCore(); + starMgr_ = GETSTELMODULE(StarMgr); + objectMgr_ = GETSTELMODULE(StelObjectMgr); + Q_ASSERT(objectMgr_); + connect(objectMgr_, &StelObjectMgr::selectedObjectChanged, this, &LensDistortionEstimatorDialog::handleSelectedObjectChange); + handleSelectedObjectChange(StelModule::ReplaceSelection); + + // Main tabs + clearWarnings(); + connect(ui_->imageFilePath, &QLineEdit::editingFinished, this, &LensDistortionEstimatorDialog::updateImage); + connect(ui_->imageFilePath, &QLineEdit::textChanged, this, &LensDistortionEstimatorDialog::updateImagePathStatus); + connect(ui_->imageBrowseBtn, &QAbstractButton::clicked, this, &LensDistortionEstimatorDialog::browseForImage); + connect(ui_->pickImagePointBtn, &QAbstractButton::clicked, this, &LensDistortionEstimatorDialog::startImagePointPicking); + connect(ui_->removePointBtn, &QAbstractButton::clicked, this, &LensDistortionEstimatorDialog::removeImagePoint); + connect(ui_->imagePointsPicked, &QTreeWidget::itemSelectionChanged, this, &LensDistortionEstimatorDialog::handlePointSelectionChange); + connect(ui_->imagePointsPicked, &QTreeWidget::itemDoubleClicked, this, &LensDistortionEstimatorDialog::handlePointDoubleClick); + connect(ui_->repositionBtn, &QAbstractButton::clicked, this, &LensDistortionEstimatorDialog::repositionImageByPoints); + connect(ui_->fixWarningBtn, &QPushButton::clicked, this, &LensDistortionEstimatorDialog::fixWarning); + connect(ui_->exportPointsBtn, &QPushButton::clicked, this, &LensDistortionEstimatorDialog::exportPoints); + connect(ui_->importPointsBtn, &QPushButton::clicked, this, &LensDistortionEstimatorDialog::importPoints); + connect(ui_->distortionModelCB, qOverload(&QComboBox::currentIndexChanged), this, + &LensDistortionEstimatorDialog::updateDistortionCoefControls); + connect(ui_->resetPlacementBtn, &QPushButton::clicked, this, &LensDistortionEstimatorDialog::resetImagePlacement); + connect(ui_->resetDistortionBtn, &QPushButton::clicked, this, &LensDistortionEstimatorDialog::resetDistortion); +#ifdef HAVE_NLOPT + connect(ui_->optimizeBtn, &QAbstractButton::clicked, this, &LensDistortionEstimatorDialog::optimizeParameters); +#else + ui_->optimizeBtn->hide(); +#endif + connect(ui_->generateLensfunModelBtn, &QPushButton::clicked, this, &LensDistortionEstimatorDialog::generateLensfunModel); + connect(ui_->lensMake, &QLineEdit::textChanged, this, &LensDistortionEstimatorDialog::setupGenerateLensfunModelButton); + connect(ui_->lensModel, &QLineEdit::textChanged, this, &LensDistortionEstimatorDialog::setupGenerateLensfunModelButton); + connect(ui_->lensMount, &QLineEdit::textChanged, this, &LensDistortionEstimatorDialog::setupGenerateLensfunModelButton); + connect(ui_->imageFilePath, &QLineEdit::textChanged, this, &LensDistortionEstimatorDialog::setupGenerateLensfunModelButton); + setupGenerateLensfunModelButton(); + + ui_->errorsChartView->setChart(errorsChart_); + ui_->errorsChartView->setRenderHint(QPainter::Antialiasing); + ui_->imagePointsPicked->header()->setSectionResizeMode(QHeaderView::ResizeToContents); + ui_->listChartSplitter->setStretchFactor(0, 1); + ui_->listChartSplitter->setStretchFactor(1, 4); + updateDistortionCoefControls(); + + // Settings tab + connectBoolProperty(ui_->imageEnabled, "LensDistortionEstimator.imageEnabled"); + connectBoolProperty(ui_->imageAxesEnabled, "LensDistortionEstimator.imageAxesEnabled"); + connectBoolProperty(ui_->markPickedPoints, "LensDistortionEstimator.pointMarkersEnabled"); + connectBoolProperty(ui_->markProjectionCenter, "LensDistortionEstimator.projectionCenterMarkerEnabled"); + ui_->imageAxesColor->setup("LensDistortionEstimator.imageAxesColor", "LensDistortionEstimator/image_axes_color"); + ui_->pointMarkerColor->setup("LensDistortionEstimator.pointMarkerColor", "LensDistortionEstimator/point_marker_color"); + ui_->selectedPointMarkerColor->setup("LensDistortionEstimator.selectedPointMarkerColor", + "LensDistortionEstimator/selected_point_marker_color"); + ui_->projectionCenterMarkerColor->setup("LensDistortionEstimator.projectionCenterMarkerColor", + "LensDistortionEstimator/center_of_projection_marker_color"); + connect(ui_->restoreDefaultsBtn, &QPushButton::clicked, this, &LensDistortionEstimatorDialog::restoreDefaults); + + // About tab + setAboutText(); + if(gui) + { + ui_->about->document()->setDefaultStyleSheet(QString(gui->getStelStyle().htmlStyleSheet)); + } + init(); + + initialized_ = true; +} + +void LensDistortionEstimatorDialog::restoreDefaults() +{ + if(askConfirmation()) + { + qDebug() << "[LensDistortionEstimator] restore defaults..."; + GETSTELMODULE(LensDistortionEstimator)->restoreDefaultSettings(); + } + else + { + qDebug() << "[LensDistortionEstimator] restore defaults is canceled..."; + } +} + +void LensDistortionEstimatorDialog::emitWarning(const QString& text, const bool autoFixable) +{ + ui_->warnings->setStyleSheet("margin-left: 1px; border-radius: 5px; background-color: #ff5757;"); + const auto oldText = ui_->warnings->text(); + if(oldText.isEmpty()) + ui_->warnings->setText(text); + else + ui_->warnings->setText(text+"\n"+oldText); + ui_->warnings->show(); + if(autoFixable) + { + ui_->fixWarningBtn->show(); + ui_->fixWarningBtn->setToolTip(q_("Sets simulation date/time to image EXIF date/time, and freezes simulation time")); + } +} + +void LensDistortionEstimatorDialog::clearWarnings() +{ + ui_->warnings->hide(); + ui_->warnings->setText(""); + ui_->fixWarningBtn->hide(); +} + +double LensDistortionEstimatorDialog::computeExifTimeDifference() const +{ + const auto jd = core_->getJD(); + const auto stelTZ = core_->getUTCOffset(jd); // in hours + const auto jdc = jd + stelTZ/24.; // UTC -> local TZ + int year, month, day, hour, minute, second; + StelUtils::getDateTimeFromJulianDay(jdc, &year, &month, &day, &hour, &minute, &second); + const bool exifHasLocalTime = exifDateTime_.timeSpec() == Qt::LocalTime; + const auto currentDateTime = exifHasLocalTime ? + QDateTime(QDate(year, month, day), QTime(hour, minute, second, 0), Qt::LocalTime) : + QDateTime(QDate(year, month, day), QTime(hour, minute, second, 0), QTimeZone(stelTZ*3600)); + return exifDateTime_.msecsTo(currentDateTime); +} + +void LensDistortionEstimatorDialog::periodicWarningsCheck() +{ + if(!exifDateTime_.isValid()) return; + clearWarnings(); + const bool exifHasLocalTime = exifDateTime_.timeSpec() == Qt::LocalTime; + const auto timeDiff = computeExifTimeDifference(); + if(std::abs(timeDiff) > 5000) + emitWarning(q_("Stellarium time differs from image EXIF time %1 by %2 seconds") + .arg(exifHasLocalTime ? exifDateTime_.toString("yyyy-MM-dd hh:mm:ss") + : exifDateTime_.toString("yyyy-MM-dd hh:mm:ss UTC+ttt")) + .arg(timeDiff/1000.) + .replace("UTC+UTCUTCUTC", "UTC+0") /* XXX: a Qt 6.4.1 bug I suppose */, + true); +} + +void LensDistortionEstimatorDialog::fixWarning() +{ + if(!exifDateTime_.isValid()) return; + + const auto timeDiff = computeExifTimeDifference(); + const auto jd = core_->getJD(); + core_->setJD(jd - timeDiff/86400e3); + core_->setZeroTimeSpeed(); + + clearWarnings(); +} + +void LensDistortionEstimatorDialog::exportPoints() const +{ + const auto path = QFileDialog::getSaveFileName(&StelMainView::getInstance(), q_("Save points to file"), {}, + q_("Comma-separated values (*.csv)")); + if(path.isNull()) return; + QFile file(path); + if(!file.open(QFile::WriteOnly)) + { + QMessageBox::critical(&StelMainView::getInstance(), q_("Error saving data"), + q_("Failed to open output file for writing:\n%1").arg(file.errorString())); + return; + } + QTextStream str(&file); + str << "Object name,Localized object name,Image position X,Image position Y,Exact object azimuth,Exact object elevation\n"; + const auto nMax = ui_->imagePointsPicked->topLevelItemCount(); + for(int n = 0; n < nMax; ++n) + { + const auto item = ui_->imagePointsPicked->topLevelItem(n); + const auto objectName = item->data(Column::ObjectName, Role::ObjectEnglishName).toString(); + const auto object = findObjectByName(objectName); + double trueAzimuth = NAN, trueElevation = NAN; + if(!object) + { + qWarning() << "Object" << objectName << "not found, will export its sky position as NaN"; + } + else + { + using namespace std; + const auto objDir = normalize(object->getAltAzPosApparent(core_)); + trueElevation = 180/M_PI * asin(clamp(objDir[2], -1.,1.)); + trueAzimuth = 180/M_PI * atan2(objDir[1], -objDir[0]); + } + str << objectName << "," + << item->data(Column::ObjectName, Qt::DisplayRole).toString() << "," + << item->data(Column::x, Qt::DisplayRole).toDouble() << "," + << item->data(Column::y, Qt::DisplayRole).toDouble() << "," + << (trueAzimuth < 0 ? 360 : 0) + trueAzimuth << "," + << trueElevation << "\n"; + } + str.flush(); + if(!file.flush()) + { + QMessageBox::critical(&StelMainView::getInstance(), q_("Error saving data"), + q_("Failed to write the output file:\n%1").arg(file.errorString())); + } +} + +void LensDistortionEstimatorDialog::importPoints() +{ + if(ui_->imagePointsPicked->topLevelItemCount()) + { + const auto ret = QMessageBox::warning(&StelMainView::getInstance(), q_("Confirm deletion of old points"), + q_("There are some points defined. They will be deleted " + "and replaced with the ones from the file. Proceed?"), + QMessageBox::Yes|QMessageBox::Cancel, QMessageBox::Cancel); + if(ret != QMessageBox::Yes) return; + } + + const auto path = QFileDialog::getOpenFileName(&StelMainView::getInstance(), q_("Save points to file"), {}, + q_("Comma-separated values (*.csv)")); + if(path.isNull()) return; + QFile file(path); + if(!file.open(QFile::ReadOnly)) + { + QMessageBox::critical(&StelMainView::getInstance(), q_("Error reading data"), + q_("Failed to open input file for reading:\n%1").arg(file.errorString())); + return; + } + + // Keep new items out of the widget till the end, so that, if parsing fails, we didn't destroy the original widget contents + std::vector> newItems; + + QTextStream str(&file); + for(int lineNum = 1;; ++lineNum) + { + if(str.atEnd()) break; + + const auto entries = str.readLine().split(","); + constexpr int minEntryCount = 4; + if(entries.size() < minEntryCount) + { + QMessageBox::critical(&StelMainView::getInstance(), q_("File format error"), + q_("Can't import file: expected to find at least %1 entries at line %2, but got %3.") + .arg(minEntryCount).arg(lineNum).arg(entries.size())); + return; + } + if(lineNum==1) continue; + const auto objectName = entries[0]; + const auto object = StelApp::getInstance().getStelObjectMgr().searchByName(objectName); + if(!object) + { + QMessageBox::critical(&StelMainView::getInstance(), q_("File format error"), + q_("Can't import file: object \"%1\" not found at line %2.").arg(objectName).arg(lineNum)); + return; + } + const auto locObjName = entries[1]; + bool ok = false; + const auto imgPosX = entries[2].toDouble(&ok); + if(!ok) + { + QMessageBox::critical(&StelMainView::getInstance(), q_("File format error"), + q_("Can't import file: failed to parse image position X at line %1.").arg(lineNum)); + return; + } + const auto imgPosY = entries[3].toDouble(&ok); + if(!ok) + { + QMessageBox::critical(&StelMainView::getInstance(), q_("File format error"), + q_("Can't import file: failed to parse image position Y at line %1.").arg(lineNum)); + return; + } + const auto item = new QTreeWidgetItem({locObjName, + QString("%1").arg(imgPosX), + QString("%1").arg(imgPosY)}); + item->setFlags(item->flags() & ~Qt::ItemIsDropEnabled); + item->setData(Column::ObjectName, Role::ObjectEnglishName, objectName); + newItems.emplace_back(item); + } + + for(int n = ui_->imagePointsPicked->topLevelItemCount() - 1; n >= 0; --n) + delete ui_->imagePointsPicked->topLevelItem(n); + + for(auto& item : newItems) + ui_->imagePointsPicked->addTopLevelItem(item.release()); + + updateRepositionButtons(); + updateErrorsChart(); + ui_->exportPointsBtn->setEnabled(true); +} + +void LensDistortionEstimatorDialog::removeImagePoint() +{ + const auto item = ui_->imagePointsPicked->currentItem(); + if(!item) return; + delete item; + + updateRepositionButtons(); + updateErrorsChart(); + if(ui_->imagePointsPicked->topLevelItemCount() == 0) + ui_->exportPointsBtn->setDisabled(true); +} + +void LensDistortionEstimatorDialog::updateRepositionButtons() +{ + ui_->repositionBtn->setEnabled(ui_->imagePointsPicked->topLevelItemCount() >= 2 && !image_.isNull()); + ui_->optimizeBtn->setEnabled(ui_->imagePointsPicked->topLevelItemCount() >= 3 && !image_.isNull()); +} + +void LensDistortionEstimatorDialog::handlePointSelectionChange() +{ + ui_->removePointBtn->setEnabled(!!ui_->imagePointsPicked->currentItem()); +} + +void LensDistortionEstimatorDialog::handlePointDoubleClick(QTreeWidgetItem*const item, int) +{ + const auto objectName = item->data(Column::ObjectName, Role::ObjectEnglishName).toString(); + const auto object = findObjectByName(objectName); + if(object) objectMgr_->setSelectedObject(object); +} + +void LensDistortionEstimatorDialog::handleSelectedObjectChange(StelModule::StelModuleSelectAction) +{ + if(optimizing_) return; + + if(objectMgr_->getWasSelected() && !image_.isNull()) + { + ui_->pickImagePointBtn->setEnabled(true); + } + else + { + ui_->pickImagePointBtn->setChecked(false); + ui_->pickImagePointBtn->setDisabled(true); + QApplication::restoreOverrideCursor(); + isPickingAPoint_ = false; + } +} + +void LensDistortionEstimatorDialog::startImagePointPicking(const bool buttonChecked) +{ + if(!buttonChecked) + { + QApplication::restoreOverrideCursor(); + isPickingAPoint_ = false; + return; + } + QApplication::setOverrideCursor(lde_->getPointPickCursor()); + isPickingAPoint_ = true; +} + +void LensDistortionEstimatorDialog::togglePointPickingMode() +{ + if(!initialized_ || !ui_->pickImagePointBtn->isEnabled()) return; + const bool mustBeChecked = !ui_->pickImagePointBtn->isChecked(); + ui_->pickImagePointBtn->setChecked(mustBeChecked); + startImagePointPicking(mustBeChecked); +} + +void LensDistortionEstimatorDialog::resetImagePointPicking() +{ + QApplication::restoreOverrideCursor(); + isPickingAPoint_ = false; + ui_->pickImagePointBtn->setChecked(false); +} + +void LensDistortionEstimatorDialog::registerImagePoint(const StelObject& object, const Vec2d& posInImage) +{ + const auto englishName = getObjectName(object); + auto name = object.getNameI18n(); + if(name.isEmpty()) + name = englishName; + const auto item = new QTreeWidgetItem({name, + QString("%1").arg(posInImage[0]), + QString("%1").arg(posInImage[1])}); + item->setFlags(item->flags() & ~Qt::ItemIsDropEnabled); + item->setData(Column::ObjectName, Role::ObjectEnglishName, englishName); + ui_->imagePointsPicked->addTopLevelItem(item); + ui_->imagePointsPicked->setCurrentItem(item); + + updateRepositionButtons(); + updateErrorsChart(); + ui_->exportPointsBtn->setEnabled(true); +} + +void LensDistortionEstimatorDialog::setupChartAxisStyle(QValueAxis& axis) +{ + // Same as in AstroCalcChart + static const QPen pen(Qt::white, 1, Qt::SolidLine); + static const QPen gridPen(Qt::white, 0.5, Qt::SolidLine); + static const QPen minorGridPen(Qt::white, 0.35, Qt::DotLine); + + axis.setLinePen(pen); + axis.setGridLinePen(gridPen); + axis.setMinorGridLinePen(minorGridPen); + axis.setTitleBrush(Qt::white); + axis.setLabelsBrush(Qt::white); +} + +double LensDistortionEstimatorDialog::roundToNiceAxisValue(const double x) +{ + if(!std::isfinite(x)) return x; + + const bool negative = x<0; + // Going through serialization to avoid funny issues with imprecision. + // This function is not supposed to be called in performance-critical tasks. + const auto str = QString::number(std::abs(x), 'e', 0); + auto absRounded = str.toDouble(); + + if(absRounded >= std::abs(x)) + return negative ? -absRounded : absRounded; + + // We want to round away from zero, so need to tweak the result + const int firstDigit = str[0].toLatin1() - '0'; + assert(firstDigit >= 0 && firstDigit <= 9); + absRounded *= double(firstDigit+1) / firstDigit; + return negative ? -absRounded : absRounded; +} + +void LensDistortionEstimatorDialog::updateErrorsChart() +{ + errorsChart_->removeAllSeries(); + for(const auto axis : errorsChart_->axes()) + { + errorsChart_->removeAxis(axis); + delete axis; + } + + if(image_.isNull()) return; + + using namespace std; + constexpr auto radian = 180/M_PI; + + const auto series = new QScatterSeries; + const auto nMax = ui_->imagePointsPicked->topLevelItemCount(); + double xMin = +INFINITY, xMax = -INFINITY, yMin = +INFINITY, yMax = -INFINITY; + for(int n = 0; n < nMax; ++n) + { + const auto item = ui_->imagePointsPicked->topLevelItem(n); + const auto objectName = item->data(Column::ObjectName, Role::ObjectEnglishName).toString(); + const auto object = findObjectByName(objectName); + if(!object) + { + qWarning() << "Object" << objectName << "not found, chart creation fails"; + delete series; + return; + } + const auto objectDirTrue = normalize(object->getAltAzPosApparent(core_)); + const Vec2d objectXY(item->data(Column::x, Qt::DisplayRole).toDouble(), + item->data(Column::y, Qt::DisplayRole).toDouble()); + const auto objectDirByImg = computeImagePointDir(objectXY); + const auto centerDir = projectionCenterDir(); + + const auto distFromCenterTrue = acos(clamp(dot(objectDirTrue, centerDir), -1.,1.)); + + const auto x = distFromCenterTrue * radian; + const auto y = acos(clamp(dot(objectDirByImg,objectDirTrue), -1.,1.)) * radian; + if(x < xMin) xMin = x; + if(x > xMax) xMax = x; + if(y < yMin) yMin = y; + if(y > yMax) yMax = y; + series->append(x, y); + } + errorsChart_->addSeries(series); + + const auto xAxis = new QValueAxis; + xAxis->setTitleText(q_(u8"Distance from center, \u00b0")); + setupChartAxisStyle(*xAxis); + const auto xRange = xMax - xMin; + xAxis->setRange(0, max({xMax + 0.05 * xRange + 1, min(90., imageSmallerSideFoV()/2 * 1.5)})); + + const auto yAxis = new QValueAxis; + yAxis->setTitleText(q_(u8"Error estimate, \u00b0")); + setupChartAxisStyle(*yAxis); + const auto yRange = yMax - yMin; + const auto yRangeToUse = roundToNiceAxisValue(max(abs(yMin), abs(yMax)) + 0.1 * yRange); + if(yMin >= 0) + yAxis->setRange(0, yRangeToUse); + else if(yMax <= 0) + yAxis->setRange(-yRangeToUse, 0); + else + yAxis->setRange(-yRangeToUse, yRangeToUse); + + errorsChart_->addAxis(xAxis, Qt::AlignBottom); + errorsChart_->addAxis(yAxis, Qt::AlignLeft); + + series->setBorderColor(Qt::transparent); + series->setMarkerSize(4 * StelApp::getInstance().getDevicePixelsPerPixel()); + series->attachAxis(xAxis); + series->attachAxis(yAxis); + + errorsChart_->legend()->hide(); +} + +double LensDistortionEstimatorDialog::applyDistortion(const double ru) const +{ + const auto model = distortionModel(); + const auto ru2 = ru*ru; + const auto ru3 = ru2*ru; + const auto ru4 = ru2*ru2; + switch(model) + { + case DistortionModel::Poly3: + { + const auto k1 = distortionTerm1(); + return ru*(1-k1+k1*ru2); + } + case DistortionModel::Poly5: + { + const auto k1 = distortionTerm1(); + const auto k2 = distortionTerm2(); + return ru*(1+k1*ru2+k2*ru4); + } + case DistortionModel::PTLens: + { + const auto a = distortionTerm1(); + const auto b = distortionTerm2(); + const auto c = distortionTerm3(); + return ru*(a*ru3+b*ru2+c*ru+1-a-b-c); + } + } + qWarning() << "Unknown distortion model chosen: " << static_cast(model); + return ru; +} + +double LensDistortionEstimatorDialog::maxUndistortedR() const +{ + // The limit value is defined as the point where the distortion model ceases to be monotonic. + + const auto model = distortionModel(); + // The value for the case when the model is monotonic for all ru in [0,inf) + constexpr double unlimited = std::numeric_limits::max(); + + using namespace std; + switch(model) + { + case DistortionModel::Poly3: + { + const auto k1 = distortionTerm1(); + const auto square = (k1-1)/(3*k1); + if(square < 0) return unlimited; + return sqrt(square); + } + case DistortionModel::Poly5: + { + const auto k1 = distortionTerm1(); + const auto k2 = distortionTerm2(); + if(k2 == 0 && k1 < 0) return 1/sqrt(-3*k1); + const auto D = 9*k1*k1-20*k2; + if(D < 0) return unlimited; + const auto square = -(3*k1+sqrt(D))/(10*k2); + if(square < 0) return unlimited; + return sqrt(square); + } + case DistortionModel::PTLens: + { + const auto a = distortionTerm1(); + const auto b = distortionTerm2(); + const auto c = distortionTerm3(); + const auto roots = solveCubic(4*a, 3*b, 2*c, 1-a-b-c); + for(auto root : roots) + if(root > 0) return root; + return unlimited; + } + } + qWarning() << "Unknown distortion model chosen: " << static_cast(model); + return unlimited; +} + +Vec2d LensDistortionEstimatorDialog::imagePointToNormalizedCoords(const Vec2d& imagePoint) const +{ + // This is how it works in lensfun: i increases to the right, j increases towards the bottom of the image, and + // normPos = [(2*i - (width -1)) / min(width-1, height-1) - lensfun.center[0], + // (2*j - (height-1)) / min(width-1, height-1) - lensfun.center[1]] + // Our coordinates are inverted in y with respect to lensfun, so that they are easier + // to handle, avoiding too many (1,-1) multiplers scattered around the code. + const double W = image_.width(); + const double H = image_.height(); + const double w = W - 1; + const double h = H - 1; + const auto distortedPos = Vec2d(1,-1) * (2. * imagePoint - Vec2d(w, h)) / std::min(w, h) - imageCenterShift(); + + using namespace std; + const auto rd = distortedPos.norm(); + double ruMin = 0; + double ruMax = 1.3 * hypot(double(image_.width()),image_.height())/min(image_.width(),image_.height()); // FIXME: make a more sound choice + for(int n = 0; n < 53 && ruMin != ruMax; ++n) + { + const auto ru = (ruMin + ruMax) / 2; + const auto rdCurr = applyDistortion(ru); + if(rdCurr < rd) + ruMin = ru; + else + ruMax = ru; + } + const auto ru = (ruMin + ruMax) / 2; + + return distortedPos / rd * ru; +} + +Vec3d LensDistortionEstimatorDialog::computeImagePointDir(const Vec2d& imagePoint) const +{ + using namespace std; + const auto centeredPoint = imagePointToNormalizedCoords(imagePoint); + const auto upDir = imageUpDir(); + const auto centerDir = projectionCenterDir(); + const auto rightDir = centerDir ^ upDir; + const auto smallerSideFoV = M_PI/180 * this->imageSmallerSideFoV(); + + const auto x = centeredPoint[0] * tan(smallerSideFoV/2); + const auto y = centeredPoint[1] * tan(smallerSideFoV/2); + const auto z = 1.; + return normalize(x*rightDir + y*upDir + z*centerDir); +} + +double LensDistortionEstimatorDialog::computeAngleBetweenImagePoints(const Vec2d& pointA, const Vec2d& pointB, const double smallerSideFoV) const +{ + using namespace std; + + const auto centeredPointA = imagePointToNormalizedCoords(pointA); + const auto centeredPointB = imagePointToNormalizedCoords(pointB); + const auto tanFoV2 = tan(smallerSideFoV/2); + + auto pointA3D = Vec3d(centeredPointA[0] * tanFoV2, + centeredPointA[1] * tanFoV2, + 1); + pointA3D.normalize(); + + auto pointB3D = Vec3d(centeredPointB[0] * tanFoV2, + centeredPointB[1] * tanFoV2, + 1); + pointB3D.normalize(); + + return acos(clamp(dot(pointA3D, pointB3D), -1., 1.)); +} + +void LensDistortionEstimatorDialog::repositionImageByPoints() +{ + Vec2d object0xy, object1xy; + Vec3d object0dirTrue, object1dirTrue; + using namespace std; + + if(optimizing_) + { + if(objectPositions_.size() < 2) + { + qWarning() << "Reposition called while number of points less than 2. Ignoring."; + return; + } + object0xy = objectPositions_[0].imgPos; + object1xy = objectPositions_[1].imgPos; + object0dirTrue = objectPositions_[0].skyDir; + object1dirTrue = objectPositions_[1].skyDir; + } + else + { + if(ui_->imagePointsPicked->topLevelItemCount() < 2) + { + qWarning() << "Reposition called while number of points less than 2. Ignoring."; + return; + } + const auto item0 = ui_->imagePointsPicked->topLevelItem(0); + const auto object0Name = item0->data(Column::ObjectName, Role::ObjectEnglishName).toString(); + const auto object0 = findObjectByName(object0Name); + if(!object0) + { + qWarning() << "Object" << object0Name << "not found, repositioning fails"; + return; + } + + const auto item1 = ui_->imagePointsPicked->topLevelItem(1); + const auto object1Name = item1->data(Column::ObjectName, Role::ObjectEnglishName).toString(); + const auto object1 = findObjectByName(object1Name); + if(!object1) + { + qWarning() << "Object" << object1Name << "not found, repositioning fails"; + return; + } + + object0dirTrue = normalize(object0->getAltAzPosApparent(core_)); + object0xy = Vec2d(item0->data(Column::x, Qt::DisplayRole).toDouble(), + item0->data(Column::y, Qt::DisplayRole).toDouble()); + + object1dirTrue = normalize(object1->getAltAzPosApparent(core_)); + object1xy = Vec2d(item1->data(Column::x, Qt::DisplayRole).toDouble(), + item1->data(Column::y, Qt::DisplayRole).toDouble()); + } + auto object0dirByImg = computeImagePointDir(object0xy); + auto object1dirByImg = computeImagePointDir(object1xy); + + // Find the correct FoV to make angle between the two image points correct + const auto angleBetweenObjectsTrue = acos(clamp(dot(object0dirTrue, object1dirTrue), -1., 1.)); + double fovMax = 0.999*M_PI; + double fovMin = 0; + double fov = (fovMin + fovMax) / 2.; + for(int n = 0; n < 53 && fovMin != fovMax; ++n) + { + const auto angleBetweenObjectsByImg = computeAngleBetweenImagePoints(object0xy, object1xy, fov); + if(angleBetweenObjectsByImg > angleBetweenObjectsTrue) + fovMax = fov; + else + fovMin = fov; + fov = (fovMin + fovMax) / 2.; + } + setImageSmallerSideFoV(180/M_PI * fov); + + // Recompute the directions after FoV update + object0dirByImg = computeImagePointDir(object0xy); + object1dirByImg = computeImagePointDir(object1xy); + + // Find the matrix to rotate two image points to corresponding directions in space + const auto crossTrue = object0dirTrue ^ object1dirTrue; + const auto crossByImg = object0dirByImg ^ object1dirByImg; + const auto preRotByImg = Mat3d(object0dirByImg[0], object0dirByImg[1], object0dirByImg[2], + object1dirByImg[0], object1dirByImg[1], object1dirByImg[2], + crossByImg[0], crossByImg[1], crossByImg[2]); + const auto preRotTrue = Mat3d(object0dirTrue[0], object0dirTrue[1], object0dirTrue[2], + object1dirTrue[0], object1dirTrue[1], object1dirTrue[2], + crossTrue[0], crossTrue[1], crossTrue[2]); + const auto rotator = preRotTrue * preRotByImg.inverse(); + + const auto origUpDir = imageUpDir(); + const auto origCenterDir = projectionCenterDir(); + const auto newCenterDir = normalize(rotator * origCenterDir); + const auto elevation = asin(clamp(newCenterDir[2], -1.,1.)); + const auto azimuth = atan2(newCenterDir[1], -newCenterDir[0]); + setProjectionCenterAzimuth(180/M_PI * azimuth); + setProjectionCenterElevation(180/M_PI * elevation); + + const auto origFromCenterToTop = normalize(origCenterDir + 1e-8 * origUpDir); + const auto newFromCenterToTop = normalize(rotator * origFromCenterToTop); + // Desired up direction + const auto upDirNew = normalize(newFromCenterToTop - newCenterDir); + // Renewed up direction takes into account the new center direction but not the desired orientation yet + setImageFieldRotation(0); + const auto renewedUpDir = imageUpDir(); + const auto upDirCross = renewedUpDir ^ upDirNew; + const auto upDirDot = dot(renewedUpDir, upDirNew); + const auto dirSign = dot(upDirCross, newCenterDir) > 0 ? -1. : 1.; + const auto upDirSinAngle = dirSign * (dirSign * upDirCross).norm(); + const auto upDirCosAngle = upDirDot; + const auto fieldRotation = atan2(upDirSinAngle, upDirCosAngle); + setImageFieldRotation(180/M_PI * fieldRotation); + + if(!optimizing_) + updateErrorsChart(); +} + +void LensDistortionEstimatorDialog::resetImagePlacement() +{ + if(image_.isNull()) + { + setImageFieldRotation(0); + setProjectionCenterAzimuth(0); + setProjectionCenterElevation(0); + setImageSmallerSideFoV(60); + return; + } + + const auto ret = QMessageBox::warning(&StelMainView::getInstance(), q_("Confirm resetting image placement"), + q_("This will move the image to the center of the screen. Proceed?"), + QMessageBox::Yes|QMessageBox::Cancel, QMessageBox::Cancel); + if(ret != QMessageBox::Yes) return; + + placeImageInCenterOfScreen(); +} + +void LensDistortionEstimatorDialog::resetDistortion() +{ + const auto ret = QMessageBox::warning(&StelMainView::getInstance(), q_("Confirm resetting distortion"), + q_("This will zero out all distortion and shift coefficients. Proceed?"), + QMessageBox::Yes|QMessageBox::Cancel, QMessageBox::Cancel); + if(ret != QMessageBox::Yes) return; + + setImageCenterShiftX(0); + setImageCenterShiftY(0); + setDistortionTerm1(0); + setDistortionTerm2(0); + setDistortionTerm3(0); +} + +void LensDistortionEstimatorDialog::optimizeParameters() +{ +#ifdef HAVE_NLOPT + const auto funcToMinimize = [](const unsigned n, const double *v, double*, void*const params) + { + const auto dialog = static_cast(params); + + const double xShift = v[0]; + const double yShift = v[1]; + const double imageFieldRot = v[2]; + const double centerAzim = v[3]; + const double centerElev = v[4]; + const double a = v[5]; + const double b = dialog->distortionTerm2inUse_ ? v[6] : 0.; + const double c = dialog->distortionTerm3inUse_ ? v[7] : 0.; + + dialog->setImageCenterShiftX(xShift); + dialog->setImageCenterShiftY(yShift); + dialog->setImageFieldRotation(imageFieldRot); + dialog->setProjectionCenterAzimuth(centerAzim); + dialog->setProjectionCenterElevation(centerElev); + dialog->setDistortionTerm1(a); + dialog->setDistortionTerm2(b); + dialog->setDistortionTerm3(c); + + dialog->repositionImageByPoints(); + + double errorMeasure = 0; + for(const auto& obj : dialog->objectPositions_) + { + const auto objectDirByImg = dialog->computeImagePointDir(obj.imgPos); + + using namespace std; + const auto error = acos(clamp(dot(objectDirByImg,obj.skyDir), -1.,1.)); + errorMeasure += sqr(error); + } + return errorMeasure; + }; + + objectPositions_.clear(); + for(int n = 0; n < ui_->imagePointsPicked->topLevelItemCount(); ++n) + { + const auto item = ui_->imagePointsPicked->topLevelItem(n); + const auto objectName = item->data(Column::ObjectName, Role::ObjectEnglishName).toString(); + const auto object = findObjectByName(objectName); + if(!object) + { + qWarning() << "Object" << objectName << "not found, skipping it during optimization"; + continue; + } + const auto objectDirTrue = normalize(object->getAltAzPosApparent(core_)); + const Vec2d objectXY(item->data(Column::x, Qt::DisplayRole).toDouble(), + item->data(Column::y, Qt::DisplayRole).toDouble()); + objectPositions_.push_back({objectXY, objectDirTrue}); + } + + distortionTerm1_ = distortionTerm1(); distortionTerm1min_ = ui_->distortionTerm1->minimum(); distortionTerm1max_ = ui_->distortionTerm1->maximum(); + distortionTerm2_ = distortionTerm2(); distortionTerm2min_ = ui_->distortionTerm2->minimum(); distortionTerm2max_ = ui_->distortionTerm2->maximum(); + distortionTerm3_ = distortionTerm3(); distortionTerm3min_ = ui_->distortionTerm3->minimum(); distortionTerm3max_ = ui_->distortionTerm3->maximum(); + + imageCenterShiftX_ = imageCenterShiftX(); + imageCenterShiftXmin_ = ui_->imageCenterShiftX->minimum(); + imageCenterShiftXmax_ = ui_->imageCenterShiftX->maximum(); + + imageCenterShiftY_ = imageCenterShiftY(); + imageCenterShiftYmin_ = ui_->imageCenterShiftY->minimum(); + imageCenterShiftYmax_ = ui_->imageCenterShiftY->maximum(); + + imageSmallerSideFoV_ = ui_->imageSmallerSideFoV->value(); + + imageFieldRotation_ = imageFieldRotation(); + imageFieldRotationMin_ = ui_->imageFieldRotation->minimum(); + imageFieldRotationMax_ = ui_->imageFieldRotation->maximum(); + + projectionCenterAzimuth_ = projectionCenterAzimuth(); + projectionCenterAzimuthMin_ = ui_->projectionCenterAzimuth->minimum(); + projectionCenterAzimuthMax_ = ui_->projectionCenterAzimuth->maximum(); + + projectionCenterElevation_ = projectionCenterElevation(); + projectionCenterElevationMin_ = ui_->projectionCenterElevation->minimum(); + projectionCenterElevationMax_ = ui_->projectionCenterElevation->maximum(); + + distortionModel_ = distortionModel(); + optimizing_ = true; + + repositionImageByPoints(); + + const unsigned numVars = 5 + + distortionTerm1inUse_ + + distortionTerm2inUse_ + + distortionTerm3inUse_; + std::vector values{ + imageCenterShiftX_, + imageCenterShiftY_, + imageFieldRotation_, + projectionCenterAzimuth_, + projectionCenterElevation_, + distortionTerm1_, + distortionTerm2_, + distortionTerm3_, + }; + std::vector lowerBounds{ + imageCenterShiftXmin_, + imageCenterShiftYmin_, + imageFieldRotationMin_, + projectionCenterAzimuthMin_, + projectionCenterElevationMin_, + distortionTerm1min_, + distortionTerm2min_, + distortionTerm3min_, + }; + std::vector upperBounds{ + imageCenterShiftXmax_, + imageCenterShiftYmax_, + imageFieldRotationMax_, + projectionCenterAzimuthMax_, + projectionCenterElevationMax_, + distortionTerm1max_, + distortionTerm2max_, + distortionTerm3max_, + }; + // Remove trailing values if needed + values.resize(numVars); + lowerBounds.resize(numVars); + upperBounds.resize(numVars); + + nlopt::opt minimizer(nlopt::algorithm::LN_NELDERMEAD, numVars); + minimizer.set_lower_bounds(lowerBounds); + minimizer.set_upper_bounds(upperBounds); + minimizer.set_min_objective(funcToMinimize, this); + minimizer.set_ftol_rel(1e-3); + + try + { + double minf = INFINITY; + const auto result = minimizer.optimize(values, minf); + if(result > 0) + qDebug() << "Found minimum, error measure:" << minf; + else + qCritical() << "NLOpt failed with result " << minimizer.get_errmsg(); + } + catch(const std::exception& ex) + { + qCritical() << "NLOpt failed:" << ex.what(); + } + + optimizing_ = false; + setDistortionTerm1(distortionTerm1_); + setDistortionTerm2(distortionTerm2_); + setDistortionTerm3(distortionTerm3_); + setImageCenterShiftX(imageCenterShiftX_); + setImageCenterShiftY(imageCenterShiftY_); + setImageFieldRotation(imageFieldRotation_); + setImageSmallerSideFoV(imageSmallerSideFoV_); + setProjectionCenterAzimuth(projectionCenterAzimuth_); + setProjectionCenterElevation(projectionCenterElevation_); + + updateErrorsChart(); +#endif +} + +void LensDistortionEstimatorDialog::updateImagePathStatus() +{ + if(QFileInfo(ui_->imageFilePath->text()).exists()) + { + ui_->imageFilePath->setStyleSheet(""); + } + else + { + ui_->imageFilePath->setStyleSheet("color:red;"); + } +} + +void LensDistortionEstimatorDialog::computeColorToSubtract() +{ + // Compute the (rounded down) median of all the pixels + // First get the histogram + std::array histogramR{}; + std::array histogramG{}; + std::array histogramB{}; + const uchar*const data = image_.bits(); + const auto stride = image_.bytesPerLine(); + for(int j = 0; j < image_.height(); ++j) + { + for(int i = 0; i < image_.width(); ++i) + { + const auto r = data[j*stride+4*i+0]; + const auto g = data[j*stride+4*i+1]; + const auto b = data[j*stride+4*i+2]; + ++histogramR[r]; + ++histogramG[g]; + ++histogramB[b]; + } + } + // Now use the histogram to find the (rounded down) median + const auto totalR = std::accumulate(histogramR.begin(), histogramR.end(), size_t(0)); + const auto totalG = std::accumulate(histogramG.begin(), histogramG.end(), size_t(0)); + const auto totalB = std::accumulate(histogramB.begin(), histogramB.end(), size_t(0)); + size_t sumR = 0, sumG = 0, sumB = 0; + int midR = -1, midG = -1, midB = -1; + for(unsigned n = 0; n < histogramR.size(); ++n) + { + sumR += histogramR[n]; + sumG += histogramG[n]; + sumB += histogramB[n]; + if(midR < 0 && sumR > totalR/2) midR = n; + if(midG < 0 && sumG > totalG/2) midG = n; + if(midB < 0 && sumB > totalB/2) midB = n; + } + ui_->bgRed->setValue(midR); + ui_->bgGreen->setValue(midG); + ui_->bgBlue->setValue(midB); +} + +void LensDistortionEstimatorDialog::placeImageInCenterOfScreen() +{ + using namespace std; + + const auto mvtMan = core_->getMovementMgr(); + const auto vFoV = mvtMan->getCurrentFov(); + setImageSmallerSideFoV(clamp(0.98*vFoV, 0.,170.)); + + const auto j2000 = mvtMan->getViewDirectionJ2000(); + const auto altAz = core_->j2000ToAltAz(j2000, StelCore::RefractionOff); + const auto elevation = asin(clamp(altAz[2], -1.,1.)); + const auto azimuth = atan2(altAz[1], -altAz[0]); + setProjectionCenterElevation(180/M_PI * elevation); + setProjectionCenterAzimuth(180/M_PI * azimuth); + setImageFieldRotation(0); + + if(!mvtMan->getEquatorialMount()) + return; + + const auto proj = core_->getProjection(StelCore::FrameAltAz, StelCore::RefractionMode::RefractionOff); + const auto imageCenterDir = projectionCenterDir(); + const auto imageUpperDir = normalize(imageCenterDir + 1e-3 * imageUpDir()); + Vec3d imgCenterWin, imgUpWin; + proj->project(imageCenterDir, imgCenterWin); + proj->project(imageUpperDir, imgUpWin); + const auto angle = -atan2(imgUpWin[1]-imgCenterWin[1], imgUpWin[0]-imgCenterWin[0]) + M_PI/2; + setImageFieldRotation(180/M_PI * angle); +} + +void LensDistortionEstimatorDialog::updateImage() +{ + const auto path = ui_->imageFilePath->text(); + image_ = QImage(path).convertToFormat(QImage::Format_RGBA8888).mirrored(false, true); + updateRepositionButtons(); + + if(image_.isNull()) + { + ui_->imageFilePath->setStyleSheet("color:red;"); + return; + } + + ui_->imageFilePath->setStyleSheet(""); + imageChanged_ = true; + setupGenerateLensfunModelButton(); + computeColorToSubtract(); + placeImageInCenterOfScreen(); + +#ifdef HAVE_EXIV2 + try + { + const auto image = Exiv2::ImageFactory::open(path.toStdString()); + if(!image.get()) return; + image->readMetadata(); + const auto& exif = image->exifData(); + + QString date; + for(const auto& key : {"Exif.Photo.DateTimeOriginal", + "Exif.Image.DateTimeOriginal", + "Exif.Image.DateTime"}) + { + const auto it=exif.findKey(Exiv2::ExifKey(key)); + if(it!=exif.end()) + { + date = QString::fromStdString(it->toString()); + break; + } + } + + int timeZone = INT_MIN; + for(const auto& key : {"Exif.CanonTi.TimeZone"}) + { + const auto it=exif.findKey(Exiv2::ExifKey(key)); + if(it != exif.end() && it->typeId() == Exiv2::signedLong) + { +#if EXIV2_MINOR_VERSION < 28 + const auto num = it->toLong(); +#else + const auto num = it->toInt64(); +#endif + timeZone = num * 60; // save as seconds + break; + } + } + + QString gpsTime; + for(const auto& key : {"Exif.GPSInfo.GPSTimeStamp"}) + { + const auto it=exif.findKey(Exiv2::ExifKey(key)); + if(it!=exif.end()) + { + if(it->count() != 3) + continue; + const auto hour = it->toFloat(0); + const auto min = it->toFloat(1); + const auto sec = it->toFloat(2); + if(hour < 0 || hour > 59 || hour - std::floor(hour) != 0 || + min < 0 || min > 59 || min - std::floor(min ) != 0) + continue; + gpsTime = QString("%1:%2:%3").arg(int(hour), 2, 10, QChar('0')) + .arg(int(min), 2, 10, QChar('0')) + .arg(double(sec), 2, 'f', 0, QChar('0')); + break; + } + } + + QString gpsDate; + for(const auto& key : {"Exif.GPSInfo.GPSDateStamp"}) + { + const auto it=exif.findKey(Exiv2::ExifKey(key)); + if(it!=exif.end()) + { + gpsDate = QString::fromStdString(it->toString()); + break; + } + } + + // If GPS date and time are present, take them as more reliable. + // At the very least, time zone is present there unconditionally. + if(!gpsDate.isEmpty() && !gpsTime.isEmpty()) + { + date = gpsDate + " " + gpsTime; + timeZone = 0; + } + + exifDateTime_ = {}; + clearWarnings(); + if(date.isEmpty()) + { + emitWarning(q_("Can't get EXIF date from the image, please make " + "sure current date and time setting is correct.")); + } + else + { + exifDateTime_ = QDateTime::fromString(date, "yyyy:MM:dd HH:mm:ss"); + if(timeZone != INT_MIN) + exifDateTime_.setTimeZone(QTimeZone(timeZone)); + if(!exifDateTime_.isValid()) + { + emitWarning(q_("Failed to parse EXIF date/time, please make sure " + "current date and time setting is correct.")); + } + else + periodicWarningsCheck(); + warningCheckTimer_.start(); + } + + for(const auto& key : {"Exif.Photo.FocalLength", "Exif.Image.FocalLength"}) + { + const auto it=exif.findKey(Exiv2::ExifKey(key)); + if(it==exif.end() || it->typeId() != Exiv2::unsignedRational) + continue; + const auto [num,denom] = it->toRational(); + if(denom==0) continue; + ui_->lensFocalLength->setValue(double(num)/denom); + break; + } + + for(const auto& key : {"Exif.Photo.LensModel"}) + { + const auto it=exif.findKey(Exiv2::ExifKey(key)); + if(it != exif.end() && it->typeId() == Exiv2::signedLong) + { + const auto model = it->toString(); + if(model.empty()) continue; + ui_->lensMake->setText(model.c_str()); + break; + } + } + } + catch(Exiv2::Error& e) + { + qDebug() << "exiv2 error:" << e.what(); + } +#endif +} + +void LensDistortionEstimatorDialog::browseForImage() +{ + const auto path = QFileDialog::getOpenFileName(&StelMainView::getInstance(), q_("Open image file"), {}, + q_("Image files (*.png *.bmp *.jpg *.jpeg *.tif *.tiff *.webm *.pbm *.pgm *.ppm *.xbm *.xpm)")); + if(path.isNull()) return; + ui_->imageFilePath->setText(path); + updateImage(); +} + +void LensDistortionEstimatorDialog::init() +{ +} + +Vec3d LensDistortionEstimatorDialog::imageUpDir() const +{ + const auto az = M_PI/180 * projectionCenterAzimuth(); + const auto el = M_PI/180 * projectionCenterElevation(); + // This unrotated direction is just a derivative of projectionCenterDir with respect to elevation + const auto sinEl = std::sin(el); + const auto cosEl = std::cos(el); + const auto cosAz = std::cos(az); + const auto sinAz = std::sin(az); + const Vec3d unrotated(sinEl * cosAz, -sinEl * sinAz, cosEl); + // The rotated direction that we return is the one that takes imageFieldRotation into account. + const Vec3d centerDir(-cosEl * cosAz, cosEl * sinAz, sinEl); // same as projectionCenterDir(), but reusing sines & cosines + const Mat3d rotator = Mat4d::rotation(centerDir, -M_PI/180 * imageFieldRotation()).upper3x3(); + const Vec3d rotated = rotator * unrotated; + return Vec3d(rotated[0], rotated[1], rotated[2]); +} + +Vec3d LensDistortionEstimatorDialog::projectionCenterDir() const +{ + const auto az = M_PI/180 * projectionCenterAzimuth(); + const auto el = M_PI/180 * projectionCenterElevation(); + // Following FrameAltAz + return Vec3d(-std::cos(el) * std::cos(az), + std::cos(el) * std::sin(az), + std::sin(el)); +} + +Vec2d LensDistortionEstimatorDialog::imageCenterShift() const +{ + return 2. * Vec2d(-imageCenterShiftX(), imageCenterShiftY()) + / + std::min(image_.width()-1, image_.height()-1); +} + +QColor LensDistortionEstimatorDialog::bgToSubtract() const +{ + if(ui_->subtractBG->isChecked()) + return QColor(ui_->bgRed->value(), ui_->bgGreen->value(), ui_->bgBlue->value()); + return QColor(0,0,0); +} + +double LensDistortionEstimatorDialog::imageBrightness() const +{ + return ui_->imageBrightness->value() / 100.; +} + +double LensDistortionEstimatorDialog::distortionTerm1() const +{ + if(optimizing_) return distortionTerm1_; + return ui_->distortionTerm1->value(); +} + +double LensDistortionEstimatorDialog::distortionTerm2() const +{ + if(optimizing_) return distortionTerm2_; + return ui_->distortionTerm2->value(); +} + +double LensDistortionEstimatorDialog::distortionTerm3() const +{ + if(optimizing_) return distortionTerm3_; + return ui_->distortionTerm3->value(); +} + +double LensDistortionEstimatorDialog::imageCenterShiftX() const +{ + if(optimizing_) return imageCenterShiftX_; + return ui_->imageCenterShiftX->value(); +} + +double LensDistortionEstimatorDialog::imageCenterShiftY() const +{ + if(optimizing_) return imageCenterShiftY_; + return ui_->imageCenterShiftY->value(); +} + +double LensDistortionEstimatorDialog::imageFieldRotation() const +{ + if(optimizing_) return imageFieldRotation_; + return ui_->imageFieldRotation->value(); +} + +double LensDistortionEstimatorDialog::projectionCenterAzimuth() const +{ + if(optimizing_) return projectionCenterAzimuth_; + return ui_->projectionCenterAzimuth->value(); +} + +double LensDistortionEstimatorDialog::projectionCenterElevation() const +{ + if(optimizing_) return projectionCenterElevation_; + return ui_->projectionCenterElevation->value(); +} + +double LensDistortionEstimatorDialog::imageSmallerSideFoV() const +{ + if(optimizing_) return imageSmallerSideFoV_; + return ui_->imageSmallerSideFoV->value(); +} + +DistortionModel LensDistortionEstimatorDialog::distortionModel() const +{ + if(optimizing_) return distortionModel_; + return static_cast(ui_->distortionModelCB->currentIndex()); +} + +void LensDistortionEstimatorDialog::setDistortionTerm1(const double a) +{ + if(optimizing_) distortionTerm1_ = a; + else ui_->distortionTerm1->setValue(a); +} + +void LensDistortionEstimatorDialog::setDistortionTerm2(const double b) +{ + if(optimizing_) distortionTerm2_ = b; + else ui_->distortionTerm2->setValue(b); +} + +void LensDistortionEstimatorDialog::setDistortionTerm3(const double c) +{ + if(optimizing_) distortionTerm3_ = c; + else ui_->distortionTerm3->setValue(c); +} + +void LensDistortionEstimatorDialog::setImageCenterShiftX(const double xShift) +{ + if(optimizing_) imageCenterShiftX_ = xShift; + else ui_->imageCenterShiftX->setValue(xShift); +} + +void LensDistortionEstimatorDialog::setImageCenterShiftY(const double yShift) +{ + if(optimizing_) imageCenterShiftY_ = yShift; + else ui_->imageCenterShiftY->setValue(yShift); +} + +void LensDistortionEstimatorDialog::setImageFieldRotation(const double imageFieldRot) +{ + if(optimizing_) imageFieldRotation_ = imageFieldRot; + else ui_->imageFieldRotation->setValue(imageFieldRot); +} + +void LensDistortionEstimatorDialog::setProjectionCenterAzimuth(const double centerAzim) +{ + if(optimizing_) projectionCenterAzimuth_ = centerAzim; + else ui_->projectionCenterAzimuth->setValue(centerAzim); +} + +void LensDistortionEstimatorDialog::setProjectionCenterElevation(const double centerElev) +{ + if(optimizing_) projectionCenterElevation_ = centerElev; + else ui_->projectionCenterElevation->setValue(centerElev); +} + +void LensDistortionEstimatorDialog::setImageSmallerSideFoV(const double fov) +{ + if(optimizing_) imageSmallerSideFoV_ = fov; + else ui_->imageSmallerSideFoV->setValue(fov); +} + +auto LensDistortionEstimatorDialog::imagePointDirections() const -> std::vector +{ + std::vector statuses; + const auto nMax = ui_->imagePointsPicked->topLevelItemCount(); + for(int n = 0; n < nMax; ++n) + { + const auto item = ui_->imagePointsPicked->topLevelItem(n); + const Vec2d objectXY(item->data(Column::x, Qt::DisplayRole).toDouble(), + item->data(Column::y, Qt::DisplayRole).toDouble()); + statuses.push_back({computeImagePointDir(objectXY), item->isSelected()}); + } + return statuses; +} + +QString LensDistortionEstimatorDialog::getObjectName(const StelObject& object) const +{ + auto name = object.getEnglishName(); + if(!name.isEmpty()) return name; + + double raJ2000, decJ2000; + StelUtils::rectToSphe(&raJ2000, &decJ2000, object.getJ2000EquatorialPos(core_)); + return QString("unnamed{RA=%1;DEC=%2}").arg(raJ2000 * 12/M_PI,0,'g',17).arg(decJ2000 * 180/M_PI,0,'g',17); +} + +StelObjectP LensDistortionEstimatorDialog::findObjectByName(const QString& name) const +{ + const QLatin1String unnamedObjectPrefix("unnamed{RA="); + if(!name.startsWith(unnamedObjectPrefix)) + { + const auto obj = objectMgr_->searchByName(name); + if(!obj) return nullptr; + // Sometimes we get a completely wrong object when the correct one isn't found + if(obj->getEnglishName() != name) return nullptr; + return obj; + } + + // Parse the unnamed object + if(!name.endsWith("}")) return nullptr; + const auto raDecStr = name.mid(unnamedObjectPrefix.size(), name.size() - unnamedObjectPrefix.size() - 1); + const auto raDecStrings = raDecStr.split(";DEC="); + if(raDecStrings.size() != 2 || raDecStrings[0].isEmpty() || raDecStrings[1].isEmpty()) + { + qWarning() << "Failed to split unnamed object name" << raDecStrings; + return nullptr; + } + bool okRA = false, okDE = false; + const auto ra = raDecStrings[0].toDouble(&okRA); + const auto dec = raDecStrings[1].toDouble(&okDE); + if(!okRA || !okDE) + { + qWarning() << "Failed to parse RA" << raDecStrings[0] << "or DE" << raDecStrings[1]; + return nullptr; + } + + // There's no specific string handle for unnamed objects, so they have to be searched by coordinates + Vec3d coords(0,0,0); + StelUtils::spheToRect(ra*M_PI/12, dec*M_PI/180, coords); + + QList stars = starMgr_->searchAround(coords, 1e-4, core_); + StelObjectP object = nullptr; + double minDist = INFINITY; + for (const auto& currObj : stars) + { + const auto dist = coords.angle(currObj->getJ2000EquatorialPos(core_)); + if (dist < minDist) + { + minDist = dist; + object = currObj; + } + } + return object; +} + +void LensDistortionEstimatorDialog::updateDistortionCoefControls() +{ + const auto model = distortionModel(); + distortionTerm1inUse_ = true; + // NOTE: simply hide/show is not sufficient to query this status when the tab is not at + // the page with these controls, so we use enable/disable in addition to visibility. + switch(model) + { + case DistortionModel::Poly3: + ui_->distortionTerm2->hide(); distortionTerm2inUse_ = false; + ui_->distortionTerm3->hide(); distortionTerm3inUse_ = false; + break; + case DistortionModel::Poly5: + ui_->distortionTerm2->show(); distortionTerm2inUse_ = true; + ui_->distortionTerm3->hide(); distortionTerm3inUse_ = false; + break; + case DistortionModel::PTLens: + ui_->distortionTerm2->show(); distortionTerm2inUse_ = true; + ui_->distortionTerm3->show(); distortionTerm3inUse_ = true; + break; + } +} + +void LensDistortionEstimatorDialog::setupGenerateLensfunModelButton() +{ + const bool anyEmpty = ui_->lensMake->text().isEmpty() || + ui_->lensModel->text().isEmpty() || + ui_->lensMount->text().isEmpty() || + image_.isNull(); + ui_->generateLensfunModelBtn->setEnabled(!anyEmpty); + if(anyEmpty) + ui_->generateLensfunModelBtn->setToolTip(q_("Make, model and mount fields must be filled")); + else + ui_->generateLensfunModelBtn->setToolTip(""); +} + +void LensDistortionEstimatorDialog::generateLensfunModel() +{ + const auto lensfunCenter = -imageCenterShift(); + const auto model = distortionModel(); + QString distModel; + switch(model) + { + case DistortionModel::Poly3: + distModel = QString(R"()") + .arg(ui_->lensFocalLength->value(), 0, 'g', 7) + .arg(distortionTerm1()); + break; + case DistortionModel::Poly5: + distModel = QString(R"()") + .arg(ui_->lensFocalLength->value(), 0, 'g', 7) + .arg(distortionTerm1()) + .arg(distortionTerm2()); + break; + case DistortionModel::PTLens: + distModel = QString(R"()") + .arg(ui_->lensFocalLength->value(), 0, 'g', 7) + .arg(distortionTerm1()) + .arg(distortionTerm2()) + .arg(distortionTerm3()); + break; + default: + ui_->lensfunModelXML->setText(q_("Internal error: bad distortion model chosen")); + return; + } + ui_->lensfunModelXML->setText(QString(1+R"( + + %1 + %2 + %3 + %4 + %5 +
+ + %8 + + +)").arg(ui_->lensMake->text()) + .arg(ui_->lensModel->text()) + .arg(ui_->lensMount->text()) + .arg(computeLensCropFactor(), 0, 'g', 7) + .arg(double(image_.width())/image_.height(), 0, 'g', 7) + .arg(lensfunCenter[0], 0, 'f', 7) + .arg(lensfunCenter[1], 0, 'f', 7) + .arg(distModel)); +} + +double LensDistortionEstimatorDialog::computeLensCropFactor() const +{ + const double f = ui_->lensFocalLength->value(); + const double alpha = M_PI/180 * imageSmallerSideFoV(); + const double w = 2*f * std::tan(alpha/2); + const double aspect = double(image_.width())/image_.height(); + const double h = w/aspect; + return 36 / std::hypot(w,h); +} + +void LensDistortionEstimatorDialog::showSettingsPage() +{ + ui_->tabs->setCurrentIndex(Page::Settings); +} + +void LensDistortionEstimatorDialog::showComputationPage() +{ + const auto current = ui_->tabs->currentIndex(); + if(current != Page::ImageLoading && current != Page::Fitting) + ui_->tabs->setCurrentIndex(Page::ImageLoading); +} + +void LensDistortionEstimatorDialog::setAboutText() +{ + QString html = "" + "

" + q_("Lens Distortion Estimator Plug-in") + "

" + "" + "" + "" + "" + "
" + q_("Version") + ":" + LENSDISTORTIONESTIMATOR_PLUGIN_VERSION + "
" + q_("License") + ":" + LENSDISTORTIONESTIMATOR_PLUGIN_LICENSE + "
" + q_("Author") + ":Ruslan Kabatsayev <b7.10110111+stellarium@gmail.com>
"; + // Overview + html += "

" + q_("Overview") + "

"; + + html += "

" + q_("This plugin lets you estimate distortion of the objective lens of your camera and generate " + "a distortion model for use with the lensfun library.") + "

"; + html += "

" + q_("The process is as follows.") + "

"; + html += "
    " + "
  1. " + q_("Take a photo of the night sky so that you can easily find enough stars scattered from the " + "center of the frame to the border. Some care must be taken to get a useful shot, see more " + "details in Preparing a photo section.\n") + "
    2. " + "
  2. " + q_("On the Image Loading tab choose the image to load and use the available controls (or " + "hotkeys) to position it so that it approximately matches the direction where the camera " + "looked and corresponding orientation. You need to be able to tell which point in the image " + "corresponds to which simulated star/planet.") + "
    3. " + "
  3. " + q_("Then you switch to the Fitting tab, and pick some points in the image for a set " + "of celestial objects. To do this,\n" + "
    1. Select an object as you'd normally do.
      2. \n" + "
    2. Click Pick an image point for selected object... button (or use a hotkey).
      3. \n" + "
    3. Click the point in the image corresponding to this object (hide the simulated stars " + "if they are too confusing at this point).
    ") + "
    4. " + "
  4. " + q_("Finally, when enough points have been picked to fill the different distances from image " + "center, click Optimize button. After a few moments the image should appear aligned " + "with the stars. If it's not, try picking more points from the ones which don't match, and " + "retry optimization.") + "
    5. " + "
  5. " + q_("On the Model for lensfun tab you can generate the XML entry for lensfun database. " + "Enter make and model for the lens, and the mount (the name must match one of the mounts " + "supported by your camera as listed in the database).") + "
    6. " + "
"; + + html += "

" + q_("Preparing a photo") + "

"; + html += q_("To make a useful shot, you should follow some guidelines:\n"); + html += "
    " + "
  • " + q_("Try to keep away from the horizon, because star elevation there is very dependent on " + "weather conditions due to atmospheric refraction.") + "
    • " + "
  • " + q_("Be sure to disable any geometry distortion correction, because failure to disable " + "distortion correction may make good fitting impossible. In particular:") + + "
      " + "
    • " + q_("If you're shooting RAW, lens distortion correction must be turned off in RAW development " + "software like Lightroom or Darktable that you use for processing of your shots.") + "
      • " + "
    • " + q_("If you're shooting JPEG, lens distortion correction must be turned off in the camera " + "settings before taking the photo.") + "
      • " + "
    " + "
    • " + "
  • " + q_("Disable automatic image rotation according to its orientation saved in EXIF tags. While " + "such rotation may make the image look better, the shifts of the center of projection " + "that will be computed during the fitting process will have wrong directions if such " + "rotation was applied.") + "
    • " + "
"; + html += "

" + q_("The image that you input into this estimator should contain as much EXIF information from the " + "original RAW as possible, this will let the UI tell you if e.g. current time and date set in " + "the simulation are not the same as specified in the image.") + "

"; + html += "

" + q_("One method that will take an input RAW image and output a result following the above " + "guidelines (except those about shooting) is to use dcraw (the original one from " + "Dave Coffin) or dcraw_emu (from libRaw package), plus " + "exiftool as in the following commands (adjust the file names to your needs):") + "

" + "
dcraw_emu -t 0 -T -W -Z 20230921_235104-dcraw.tiff 20230921_235104.dng\n"
+	             "exiftool -tagsfromfile 20230921_235104.jpg -Time:all -DateTimeOriginal -gps:all "
+	             "-wm cg 20230921_235104-dcraw.tiff
" + "

" + q_("Here the -t 0 option for dcraw_emu prevents automatic rotation, " + "-T generates a TIFF file instead of the default PPM, -W prevents " + "automatic brightening which may make stars not as small as they are in the photo, and " + "-Z gives an explicit name to the file.") + "

" + "

" + q_("The exiftool command copies all time-related and GPS tags from the JPEG sidecar " + "(which on my Samsung Galaxy S8 appeared to contain more complete EXIF info than the DNG) into " + "the TIFF. Your case may differ, you may have to copy from your RAW file instead.") + "

"; + + html += "

" + q_("Hot Keys") + "

"; + html += "
    " + "
  • " + q_("Ctrl+Shift & left mouse drag: move the image around the sky") + "
    • " + "
  • " + q_("Ctrl+Shift & right mouse drag: scale and rotate the image around the last drag point") + "
    • " + "
  • " + q_("Keyboard hotkeys can be configured in the Keyboard shortcuts editor (F7).") + "
    • " + "
"; + + html += StelApp::getInstance().getModuleMgr().getStandardSupportLinksInfo("Lens Distortion Estimator plugin"); + html += ""; + + ui_->about->setHtml(html); +} diff --git a/plugins/LensDistortionEstimator/src/gui/LensDistortionEstimatorDialog.hpp b/plugins/LensDistortionEstimator/src/gui/LensDistortionEstimatorDialog.hpp new file mode 100644 index 0000000000000..4e80e40c12fa6 --- /dev/null +++ b/plugins/LensDistortionEstimator/src/gui/LensDistortionEstimatorDialog.hpp @@ -0,0 +1,173 @@ +/* + * Stellarium: Lens distortion estimator plugin + * Copyright (C) 2023 Ruslan Kabatsayev + * + * This program is free software; you can redistribute it and/or + * modify it under the terms of the GNU General Public License + * as published by the Free Software Foundation; either version 2 + * of the License, or (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program; if not, write to the Free Software + * Foundation, Inc., 51 Franklin Street, Suite 500, Boston, MA 02110-1335, USA. +*/ + +#ifndef LDEDIALOG_HPP +#define LDEDIALOG_HPP + +#include +#include "VecMath.hpp" +#include "StelDialog.hpp" +#include "StelObjectType.hpp" +#include +#include +#include +#include +#include +#include "LensDistortionEstimator.hpp" + +class StarMgr; +class Ui_LDEDialog; + +class LensDistortionEstimatorDialog : public StelDialog +{ + Q_OBJECT + +public: + LensDistortionEstimatorDialog(LensDistortionEstimator* lde); + ~LensDistortionEstimatorDialog() override; + void init(); + double distortionTerm1() const; + double distortionTerm2() const; + double distortionTerm3() const; + Vec3d imageUpDir() const; + Vec3d projectionCenterDir() const; + Vec2d imageCenterShift() const; + QColor bgToSubtract() const; + double imageBrightness() const; + double imageSmallerSideFoV() const; + LensDistortionEstimator::DistortionModel distortionModel() const; + Vec3d computeImagePointDir(const Vec2d& imagePoint) const; + double applyDistortion(const double ru) const; + double maxUndistortedR() const; + QImage image() const { return image_; } + bool imageChanged() const { return imageChanged_; } + void resetImageChangeFlag() { imageChanged_ = false; } + bool isPickingAPoint() const { return isPickingAPoint_; } + void resetImagePointPicking(); + void registerImagePoint(const StelObject& object, const Vec2d& posInImage); + void repositionImageByPoints(); + void setupGenerateLensfunModelButton(); + double imageFieldRotation() const; + double projectionCenterAzimuth() const; + double projectionCenterElevation() const; + void setImageFieldRotation(double imageFieldRot); + void setProjectionCenterAzimuth(double centerAzim); + void setProjectionCenterElevation(double centerElev); + void setImageSmallerSideFoV(double fov); + double computeAngleBetweenImagePoints(const Vec2d& pointA, const Vec2d& pointB, double smallerSideFoV) const; + void togglePointPickingMode(); + void generateLensfunModel(); + void showSettingsPage(); + void showComputationPage(); + bool initialized() const { return initialized_; } + + struct ImagePointStatus + { + Vec3d direction; + bool selected; + }; + std::vector imagePointDirections() const; + +protected: + void createDialogContent() override; + void emitWarning(const QString& text, bool autoFixable = false); + void clearWarnings(); + +public slots: + void retranslate() override; + +private: + double imageCenterShiftX() const; + double imageCenterShiftY() const; + + void setDistortionTerm1(double a); + void setDistortionTerm2(double b); + void setDistortionTerm3(double c); + void setImageCenterShiftX(double xShift); + void setImageCenterShiftY(double yShift); + + void fixWarning(); + void updateImage(); + void exportPoints() const; + void importPoints(); + void setAboutText(); + void browseForImage(); + void restoreDefaults(); + void resetDistortion(); + void removeImagePoint(); + void updateErrorsChart(); + void optimizeParameters(); + void resetImagePlacement(); + void periodicWarningsCheck(); + void updateImagePathStatus(); + void computeColorToSubtract(); + void updateRepositionButtons(); + void handlePointSelectionChange(); + void placeImageInCenterOfScreen(); + void updateDistortionCoefControls(); + double computeLensCropFactor() const; + double computeExifTimeDifference() const; + void startImagePointPicking(bool buttonChecked); + Vec3d computeObjectDir(const Vec2d& imagePoint) const; + QString getObjectName(const StelObject& object) const; + StelObjectP findObjectByName(const QString& name) const; + void handlePointDoubleClick(QTreeWidgetItem* item, int column); + Vec2d imagePointToNormalizedCoords(const Vec2d& imagePoint) const; + void handleSelectedObjectChange(StelModule::StelModuleSelectAction); + static double roundToNiceAxisValue(double x); + static void setupChartAxisStyle(QValueAxis& axis); + +private: + LensDistortionEstimator* lde_ = nullptr; + StelCore* core_ = nullptr; + StarMgr* starMgr_ = nullptr; + StelObjectMgr* objectMgr_ = nullptr; + std::unique_ptr ui_; + QImage image_; + bool initialized_ = false; + bool imageChanged_ = false; + bool isPickingAPoint_ = false; + QDateTime exifDateTime_; + QTimer warningCheckTimer_; + QChart* errorsChart_ = nullptr; + + // These variables hold the state used during optimization + struct ObjectPos + { + Vec2d imgPos; + Vec3d skyDir; + }; + std::vector objectPositions_; + double distortionTerm1_ = 0, distortionTerm1min_ = 0, distortionTerm1max_ = 0; + double distortionTerm2_ = 0, distortionTerm2min_ = 0, distortionTerm2max_ = 0; + double distortionTerm3_ = 0, distortionTerm3min_ = 0, distortionTerm3max_ = 0; + double imageCenterShiftX_ = 0, imageCenterShiftXmin_ = 0, imageCenterShiftXmax_ = 0; + double imageCenterShiftY_ = 0, imageCenterShiftYmin_ = 0, imageCenterShiftYmax_ = 0; + double imageFieldRotation_ = 0, imageFieldRotationMin_ = 0, imageFieldRotationMax_ = 0; + double projectionCenterAzimuth_ = 0, projectionCenterAzimuthMin_ = 0, projectionCenterAzimuthMax_ = 0; + double projectionCenterElevation_ = 0, projectionCenterElevationMin_ = 0, projectionCenterElevationMax_ = 0; + LensDistortionEstimator::DistortionModel distortionModel_ = LensDistortionEstimator::DistortionModel::PTLens; + double imageSmallerSideFoV_ = 0; + bool optimizing_ = false; + bool distortionTerm1inUse_ = true; + bool distortionTerm2inUse_ = false; + bool distortionTerm3inUse_ = false; +}; + +#endif // LDEDIALOG_HPP diff --git a/plugins/LensDistortionEstimator/src/gui/lensDistortionEstimatorDialog.ui b/plugins/LensDistortionEstimator/src/gui/lensDistortionEstimatorDialog.ui new file mode 100644 index 0000000000000..2877e0ab536a0 --- /dev/null +++ b/plugins/LensDistortionEstimator/src/gui/lensDistortionEstimatorDialog.ui @@ -0,0 +1,916 @@ + + + LDEDialog + + + + 0 + 0 + 300 + 300 + + + + Lens Distortion Estimator + + + + 0 + + + 0 + + + 0 + + + 0 + + + 0 + + + + + Lens Distortion Estimator + + + + + + + 0 + + + false + + + + Image Loading + + + + + + Qt::Vertical + + + + 20 + 40 + + + + + + + + Image placement + + + + + + Distortion model: + + + + + + + + + px + + + -500.000000000000000 + + + 500.000000000000000 + + + 0.000000000000000 + + + + + + + + 0 + 0 + + + + , + + + + + + + px + + + -500.000000000000000 + + + 500.000000000000000 + + + 0.000000000000000 + + + + + + + + + ° + + + 7 + + + -360.000000000000000 + + + 360.000000000000000 + + + 0.100000000000000 + + + + + + + Field of view in smaller direction: + + + imageSmallerSideFoV + + + + + + + ° + + + 7 + + + 0.001000000000000 + + + 179.000000000000000 + + + 60.000000000000000 + + + + + + + Image rotation angle: + + + imageFieldRotation + + + + + + + ° + + + 7 + + + -360.000000000000000 + + + 360.000000000000000 + + + 0.100000000000000 + + + + + + + Center of projection elevation: + + + projectionCenterElevation + + + + + + + 2 + + + + 3rd order polynomial + + + + + 5th order polynomial + + + + + PTLens + + + + + + + + Image center shift: + + + + + + + Center of projection azimuth: + + + projectionCenterAzimuth + + + + + + + 7 + + + -0.100000000000000 + + + 0.100000000000000 + + + 0.001000000000000 + + + + + + + Distortion coefficients: + + + + + + + ° + + + 7 + + + -90.000000000000000 + + + 90.000000000000000 + + + 0.100000000000000 + + + + + + + 7 + + + -0.100000000000000 + + + 0.100000000000000 + + + 0.001000000000000 + + + + + + + 7 + + + -0.100000000000000 + + + 0.100000000000000 + + + 0.001000000000000 + + + + + + + + 0 + 0 + + + + Reset +placement + + + + + + + + 0 + 0 + + + + Reset +distortion + + + + + + + + + + + + Image to load: + + + + + + + + + + Browse for an image... + + + + + + + :/graphicGui/uibtFolder.png + :/graphicGui/uibtFolder-disabled.png:/graphicGui/uibtFolder.png + + + + + + + + + + + + 0 + 0 + + + + + + + true + + + Qt::LinksAccessibleByMouse|Qt::TextSelectableByKeyboard|Qt::TextSelectableByMouse + + + + + + + Fix + + + + + + + + + + + Subtract background RGB: + + + + + + + 0 + + + 255 + + + + + + + 0 + + + 255 + + + + + + + 0 + + + 255 + + + + + + + Qt::Horizontal + + + + 40 + 20 + + + + + + + + Image brightness: + + + + + + + % + + + 1 + + + 100.000000000000000 + + + 100.000000000000000 + + + + + + + + + + Fitting + + + + + + + + false + + + Reposition image by first two points + + + + + + + Pick an image point for selected object... + + + true + + + + + + + false + + + + 0 + 0 + + + + Optimize + + + + + + + + + Qt::Horizontal + + + + + + + + 0 + 0 + + + + true + + + QAbstractItemView::InternalMove + + + Qt::MoveAction + + + + Object + + + + + x, px + + + + + y, px + + + + + + + + Import list... + + + + + + + false + + + Remove point + + + + + + + false + + + Export list... + + + + + + + + + 0 + 0 + + + + + + + + + + Model for lensfun + + + + + + true + + + false + + + + + + + mm + + + 1.000000000000000 + + + 9999.000000000000000 + + + 0.100000000000000 + + + 45.000000000000000 + + + + + + + Make: + + + + + + + + + + Mount: + + + + + + + Focal length: + + + + + + + false + + + Generate + + + + + + + Model: + + + + + + + + + + + + + + + + + Settings + + + + + + Show image + + + true + + + true + + + + + + + + + + + + + Qt::Horizontal + + + + 40 + 20 + + + + + + + + normal + + + + + + + Qt::Vertical + + + + 20 + 40 + + + + + + + + Mark image axes + + + + + + + selected + + + + + + + Mark picked points + + + + + + + + + + + + + + + + + + + + + Mark center of projection + + + + + + + + + + + + + + + + + Restore defaults + + + + + + + + About + + + + 9 + + + + + true + + + true + + + + + + + + + + + + TitleBar + QFrame +
Dialog.hpp
+ 1 + + + ColorButton + QToolButton +
gui/ColorButton.hpp
+ + + QChartView + QGraphicsView +
QChartView
+ + + + tabs + imageFilePath + imageBrowseBtn + fixWarningBtn + projectionCenterAzimuth + projectionCenterElevation + imageFieldRotation + imageSmallerSideFoV + imageCenterShiftX + imageCenterShiftY + distortionModelCB + distortionTerm1 + distortionTerm2 + distortionTerm3 + subtractBG + bgRed + bgGreen + bgBlue + imageBrightness + repositionBtn + optimizeBtn + pickImagePointBtn + imagePointsPicked + removePointBtn + importPointsBtn + exportPointsBtn + errorsChartView + imageEnabled + markPickedPoints + pointMarkerColor + selectedPointMarkerColor + markProjectionCenter + projectionCenterMarkerColor + imageAxesEnabled + imageAxesColor + restoreDefaultsBtn + about + lensMake + lensModel + lensMount + lensFocalLength + generateLensfunModelBtn + lensfunModelXML + + + + diff --git a/src/core/StelApp.cpp b/src/core/StelApp.cpp index 8575eb31fed1a..5da376498d9de 100644 --- a/src/core/StelApp.cpp +++ b/src/core/StelApp.cpp @@ -150,6 +150,10 @@ Q_IMPORT_PLUGIN(TelescopeControlStelPluginInterface) Q_IMPORT_PLUGIN(SolarSystemEditorStelPluginInterface) #endif +#ifdef USE_STATIC_PLUGIN_LENSDISTORTIONESTIMATOR +Q_IMPORT_PLUGIN(LensDistortionEstimatorStelPluginInterface) +#endif + #ifdef USE_STATIC_PLUGIN_METEORSHOWERS Q_IMPORT_PLUGIN(MeteorShowersStelPluginInterface) #endif