diff --git a/LICENSE b/LICENSE index a8aa44255..29ebfa545 100644 --- a/LICENSE +++ b/LICENSE @@ -1,23 +1,661 @@ -Copyright (c) 2014, mapillary -All rights reserved. - -Redistribution and use in source and binary forms, with or without -modification, are permitted provided that the following conditions are met: - -* Redistributions of source code must retain the above copyright notice, this - list of conditions and the following disclaimer. - -* Redistributions in binary form must reproduce the above copyright notice, - this list of conditions and the following disclaimer in the documentation - and/or other materials provided with the distribution. - -THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" -AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE -IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE -DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE -FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL -DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR -SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER -CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, -OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE -OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + GNU AFFERO GENERAL PUBLIC LICENSE + Version 3, 19 November 2007 + + Copyright (C) 2007 Free Software Foundation, Inc. + Everyone is permitted to copy and distribute verbatim copies + of this license document, but changing it is not allowed. + + Preamble + + The GNU Affero General Public License is a free, copyleft license for +software and other kinds of works, specifically designed to ensure +cooperation with the community in the case of network server software. + + The licenses for most software and other practical works are designed +to take away your freedom to share and change the works. By contrast, +our General Public Licenses are intended to guarantee your freedom to +share and change all versions of a program--to make sure it remains free +software for all its users. + + 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 +them 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. + + Developers that use our General Public Licenses protect your rights +with two steps: (1) assert copyright on the software, and (2) offer +you this License which gives you legal permission to copy, distribute +and/or modify the software. + + A secondary benefit of defending all users' freedom is that +improvements made in alternate versions of the program, if they +receive widespread use, become available for other developers to +incorporate. Many developers of free software are heartened and +encouraged by the resulting cooperation. However, in the case of +software used on network servers, this result may fail to come about. +The GNU General Public License permits making a modified version and +letting the public access it on a server without ever releasing its +source code to the public. + + The GNU Affero General Public License is designed specifically to +ensure that, in such cases, the modified source code becomes available +to the community. It requires the operator of a network server to +provide the source code of the modified version running there to the +users of that server. Therefore, public use of a modified version, on +a publicly accessible server, gives the public access to the source +code of the modified version. + + An older license, called the Affero General Public License and +published by Affero, was designed to accomplish similar goals. This is +a different license, not a version of the Affero GPL, but Affero has +released a new version of the Affero GPL which permits relicensing under +this license. + + The precise terms and conditions for copying, distribution and +modification follow. + + TERMS AND CONDITIONS + + 0. Definitions. + + "This License" refers to version 3 of the GNU Affero General Public License. + + "Copyright" also means copyright-like laws that apply to other kinds of +works, such as semiconductor masks. + + "The Program" refers to any copyrightable work licensed under this +License. Each licensee is addressed as "you". "Licensees" and +"recipients" may be individuals or organizations. + + To "modify" a work means to copy from or adapt all or part of the work +in a fashion requiring copyright permission, other than the making of an +exact copy. The resulting work is called a "modified version" of the +earlier work or a work "based on" the earlier work. + + A "covered work" means either the unmodified Program or a work based +on the Program. + + To "propagate" a work means to do anything with it that, without +permission, would make you directly or secondarily liable for +infringement under applicable copyright law, except executing it on a +computer or modifying a private copy. Propagation includes copying, +distribution (with or without modification), making available to the +public, and in some countries other activities as well. + + To "convey" a work means any kind of propagation that enables other +parties to make or receive copies. Mere interaction with a user through +a computer network, with no transfer of a copy, is not conveying. + + An interactive user interface displays "Appropriate Legal Notices" +to the extent that it includes a convenient and prominently visible +feature that (1) displays an appropriate copyright notice, and (2) +tells the user that there is no warranty for the work (except to the +extent that warranties are provided), that licensees may convey the +work under this License, and how to view a copy of this License. If +the interface presents a list of user commands or options, such as a +menu, a prominent item in the list meets this criterion. + + 1. Source Code. + + The "source code" for a work means the preferred form of the work +for making modifications to it. "Object code" means any non-source +form of a work. + + A "Standard Interface" means an interface that either is an official +standard defined by a recognized standards body, or, in the case of +interfaces specified for a particular programming language, one that +is widely used among developers working in that language. + + The "System Libraries" of an executable work include anything, other +than the work as a whole, that (a) is included in the normal form of +packaging a Major Component, but which is not part of that Major +Component, and (b) serves only to enable use of the work with that +Major Component, or to implement a Standard Interface for which an +implementation is available to the public in source code form. A +"Major Component", in this context, means a major essential component +(kernel, window system, and so on) of the specific operating system +(if any) on which the executable work runs, or a compiler used to +produce the work, or an object code interpreter used to run it. + + The "Corresponding Source" for a work in object code form means all +the source code needed to generate, install, and (for an executable +work) run the object code and to modify the work, including scripts to +control those activities. However, it does not include the work's +System Libraries, or general-purpose tools or generally available free +programs which are used unmodified in performing those activities but +which are not part of the work. For example, Corresponding Source +includes interface definition files associated with source files for +the work, and the source code for shared libraries and dynamically +linked subprograms that the work is specifically designed to require, +such as by intimate data communication or control flow between those +subprograms and other parts of the work. + + The Corresponding Source need not include anything that users +can regenerate automatically from other parts of the Corresponding +Source. + + The Corresponding Source for a work in source code form is that +same work. + + 2. Basic Permissions. + + All rights granted under this License are granted for the term of +copyright on the Program, and are irrevocable provided the stated +conditions are met. This License explicitly affirms your unlimited +permission to run the unmodified Program. The output from running a +covered work is covered by this License only if the output, given its +content, constitutes a covered work. This License acknowledges your +rights of fair use or other equivalent, as provided by copyright law. + + You may make, run and propagate covered works that you do not +convey, without conditions so long as your license otherwise remains +in force. You may convey covered works to others for the sole purpose +of having them make modifications exclusively for you, or provide you +with facilities for running those works, provided that you comply with +the terms of this License in conveying all material for which you do +not control copyright. Those thus making or running the covered works +for you must do so exclusively on your behalf, under your direction +and control, on terms that prohibit them from making any copies of +your copyrighted material outside their relationship with you. + + Conveying under any other circumstances is permitted solely under +the conditions stated below. Sublicensing is not allowed; section 10 +makes it unnecessary. + + 3. Protecting Users' Legal Rights From Anti-Circumvention Law. + + No covered work shall be deemed part of an effective technological +measure under any applicable law fulfilling obligations under article +11 of the WIPO copyright treaty adopted on 20 December 1996, or +similar laws prohibiting or restricting circumvention of such +measures. + + When you convey a covered work, you waive any legal power to forbid +circumvention of technological measures to the extent such circumvention +is effected by exercising rights under this License with respect to +the covered work, and you disclaim any intention to limit operation or +modification of the work as a means of enforcing, against the work's +users, your or third parties' legal rights to forbid circumvention of +technological measures. + + 4. Conveying Verbatim Copies. + + You may convey 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; +keep intact all notices stating that this License and any +non-permissive terms added in accord with section 7 apply to the code; +keep intact all notices of the absence of any warranty; and give all +recipients a copy of this License along with the Program. + + You may charge any price or no price for each copy that you convey, +and you may offer support or warranty protection for a fee. + + 5. Conveying Modified Source Versions. + + You may convey a work based on the Program, or the modifications to +produce it from the Program, in the form of source code under the +terms of section 4, provided that you also meet all of these conditions: + + a) The work must carry prominent notices stating that you modified + it, and giving a relevant date. + + b) The work must carry prominent notices stating that it is + released under this License and any conditions added under section + 7. This requirement modifies the requirement in section 4 to + "keep intact all notices". + + c) You must license the entire work, as a whole, under this + License to anyone who comes into possession of a copy. This + License will therefore apply, along with any applicable section 7 + additional terms, to the whole of the work, and all its parts, + regardless of how they are packaged. This License gives no + permission to license the work in any other way, but it does not + invalidate such permission if you have separately received it. + + d) If the work has interactive user interfaces, each must display + Appropriate Legal Notices; however, if the Program has interactive + interfaces that do not display Appropriate Legal Notices, your + work need not make them do so. + + A compilation of a covered work with other separate and independent +works, which are not by their nature extensions of the covered work, +and which are not combined with it such as to form a larger program, +in or on a volume of a storage or distribution medium, is called an +"aggregate" if the compilation and its resulting copyright are not +used to limit the access or legal rights of the compilation's users +beyond what the individual works permit. Inclusion of a covered work +in an aggregate does not cause this License to apply to the other +parts of the aggregate. + + 6. Conveying Non-Source Forms. + + You may convey a covered work in object code form under the terms +of sections 4 and 5, provided that you also convey the +machine-readable Corresponding Source under the terms of this License, +in one of these ways: + + a) Convey the object code in, or embodied in, a physical product + (including a physical distribution medium), accompanied by the + Corresponding Source fixed on a durable physical medium + customarily used for software interchange. + + b) Convey the object code in, or embodied in, a physical product + (including a physical distribution medium), accompanied by a + written offer, valid for at least three years and valid for as + long as you offer spare parts or customer support for that product + model, to give anyone who possesses the object code either (1) a + copy of the Corresponding Source for all the software in the + product that is covered by this License, on a durable physical + medium customarily used for software interchange, for a price no + more than your reasonable cost of physically performing this + conveying of source, or (2) access to copy the + Corresponding Source from a network server at no charge. + + c) Convey individual copies of the object code with a copy of the + written offer to provide the Corresponding Source. This + alternative is allowed only occasionally and noncommercially, and + only if you received the object code with such an offer, in accord + with subsection 6b. + + d) Convey the object code by offering access from a designated + place (gratis or for a charge), and offer equivalent access to the + Corresponding Source in the same way through the same place at no + further charge. You need not require recipients to copy the + Corresponding Source along with the object code. If the place to + copy the object code is a network server, the Corresponding Source + may be on a different server (operated by you or a third party) + that supports equivalent copying facilities, provided you maintain + clear directions next to the object code saying where to find the + Corresponding Source. Regardless of what server hosts the + Corresponding Source, you remain obligated to ensure that it is + available for as long as needed to satisfy these requirements. + + e) Convey the object code using peer-to-peer transmission, provided + you inform other peers where the object code and Corresponding + Source of the work are being offered to the general public at no + charge under subsection 6d. + + A separable portion of the object code, whose source code is excluded +from the Corresponding Source as a System Library, need not be +included in conveying the object code work. + + A "User Product" is either (1) a "consumer product", which means any +tangible personal property which is normally used for personal, family, +or household purposes, or (2) anything designed or sold for incorporation +into a dwelling. In determining whether a product is a consumer product, +doubtful cases shall be resolved in favor of coverage. For a particular +product received by a particular user, "normally used" refers to a +typical or common use of that class of product, regardless of the status +of the particular user or of the way in which the particular user +actually uses, or expects or is expected to use, the product. A product +is a consumer product regardless of whether the product has substantial +commercial, industrial or non-consumer uses, unless such uses represent +the only significant mode of use of the product. + + "Installation Information" for a User Product means any methods, +procedures, authorization keys, or other information required to install +and execute modified versions of a covered work in that User Product from +a modified version of its Corresponding Source. The information must +suffice to ensure that the continued functioning of the modified object +code is in no case prevented or interfered with solely because +modification has been made. + + If you convey an object code work under this section in, or with, or +specifically for use in, a User Product, and the conveying occurs as +part of a transaction in which the right of possession and use of the +User Product is transferred to the recipient in perpetuity or for a +fixed term (regardless of how the transaction is characterized), the +Corresponding Source conveyed under this section must be accompanied +by the Installation Information. But this requirement does not apply +if neither you nor any third party retains the ability to install +modified object code on the User Product (for example, the work has +been installed in ROM). + + The requirement to provide Installation Information does not include a +requirement to continue to provide support service, warranty, or updates +for a work that has been modified or installed by the recipient, or for +the User Product in which it has been modified or installed. Access to a +network may be denied when the modification itself materially and +adversely affects the operation of the network or violates the rules and +protocols for communication across the network. + + Corresponding Source conveyed, and Installation Information provided, +in accord with this section must be in a format that is publicly +documented (and with an implementation available to the public in +source code form), and must require no special password or key for +unpacking, reading or copying. + + 7. Additional Terms. + + "Additional permissions" are terms that supplement the terms of this +License by making exceptions from one or more of its conditions. +Additional permissions that are applicable to the entire Program shall +be treated as though they were included in this License, to the extent +that they are valid under applicable law. If additional permissions +apply only to part of the Program, that part may be used separately +under those permissions, but the entire Program remains governed by +this License without regard to the additional permissions. + + When you convey a copy of a covered work, you may at your option +remove any additional permissions from that copy, or from any part of +it. (Additional permissions may be written to require their own +removal in certain cases when you modify the work.) You may place +additional permissions on material, added by you to a covered work, +for which you have or can give appropriate copyright permission. + + Notwithstanding any other provision of this License, for material you +add to a covered work, you may (if authorized by the copyright holders of +that material) supplement the terms of this License with terms: + + a) Disclaiming warranty or limiting liability differently from the + terms of sections 15 and 16 of this License; or + + b) Requiring preservation of specified reasonable legal notices or + author attributions in that material or in the Appropriate Legal + Notices displayed by works containing it; or + + c) Prohibiting misrepresentation of the origin of that material, or + requiring that modified versions of such material be marked in + reasonable ways as different from the original version; or + + d) Limiting the use for publicity purposes of names of licensors or + authors of the material; or + + e) Declining to grant rights under trademark law for use of some + trade names, trademarks, or service marks; or + + f) Requiring indemnification of licensors and authors of that + material by anyone who conveys the material (or modified versions of + it) with contractual assumptions of liability to the recipient, for + any liability that these contractual assumptions directly impose on + those licensors and authors. + + All other non-permissive additional terms are considered "further +restrictions" within the meaning of section 10. If the Program as you +received it, or any part of it, contains a notice stating that it is +governed by this License along with a term that is a further +restriction, you may remove that term. If a license document contains +a further restriction but permits relicensing or conveying under this +License, you may add to a covered work material governed by the terms +of that license document, provided that the further restriction does +not survive such relicensing or conveying. + + If you add terms to a covered work in accord with this section, you +must place, in the relevant source files, a statement of the +additional terms that apply to those files, or a notice indicating +where to find the applicable terms. + + Additional terms, permissive or non-permissive, may be stated in the +form of a separately written license, or stated as exceptions; +the above requirements apply either way. + + 8. Termination. + + You may not propagate or modify a covered work except as expressly +provided under this License. Any attempt otherwise to propagate or +modify it is void, and will automatically terminate your rights under +this License (including any patent licenses granted under the third +paragraph of section 11). + + However, if you cease all violation of this License, then your +license from a particular copyright holder is reinstated (a) +provisionally, unless and until the copyright holder explicitly and +finally terminates your license, and (b) permanently, if the copyright +holder fails to notify you of the violation by some reasonable means +prior to 60 days after the cessation. + + Moreover, your license from a particular copyright holder is +reinstated permanently if the copyright holder notifies you of the +violation by some reasonable means, this is the first time you have +received notice of violation of this License (for any work) from that +copyright holder, and you cure the violation prior to 30 days after +your receipt of the notice. + + Termination of your rights under this section does not terminate the +licenses of parties who have received copies or rights from you under +this License. If your rights have been terminated and not permanently +reinstated, you do not qualify to receive new licenses for the same +material under section 10. + + 9. Acceptance Not Required for Having Copies. + + You are not required to accept this License in order to receive or +run a copy of the Program. Ancillary propagation of a covered work +occurring solely as a consequence of using peer-to-peer transmission +to receive a copy likewise does not require acceptance. However, +nothing other than this License grants you permission to propagate or +modify any covered work. These actions infringe copyright if you do +not accept this License. Therefore, by modifying or propagating a +covered work, you indicate your acceptance of this License to do so. + + 10. Automatic Licensing of Downstream Recipients. + + Each time you convey a covered work, the recipient automatically +receives a license from the original licensors, to run, modify and +propagate that work, subject to this License. You are not responsible +for enforcing compliance by third parties with this License. + + An "entity transaction" is a transaction transferring control of an +organization, or substantially all assets of one, or subdividing an +organization, or merging organizations. If propagation of a covered +work results from an entity transaction, each party to that +transaction who receives a copy of the work also receives whatever +licenses to the work the party's predecessor in interest had or could +give under the previous paragraph, plus a right to possession of the +Corresponding Source of the work from the predecessor in interest, if +the predecessor has it or can get it with reasonable efforts. + + You may not impose any further restrictions on the exercise of the +rights granted or affirmed under this License. For example, you may +not impose a license fee, royalty, or other charge for exercise of +rights granted under this License, and you may not initiate litigation +(including a cross-claim or counterclaim in a lawsuit) alleging that +any patent claim is infringed by making, using, selling, offering for +sale, or importing the Program or any portion of it. + + 11. Patents. + + A "contributor" is a copyright holder who authorizes use under this +License of the Program or a work on which the Program is based. The +work thus licensed is called the contributor's "contributor version". + + A contributor's "essential patent claims" are all patent claims +owned or controlled by the contributor, whether already acquired or +hereafter acquired, that would be infringed by some manner, permitted +by this License, of making, using, or selling its contributor version, +but do not include claims that would be infringed only as a +consequence of further modification of the contributor version. For +purposes of this definition, "control" includes the right to grant +patent sublicenses in a manner consistent with the requirements of +this License. + + Each contributor grants you a non-exclusive, worldwide, royalty-free +patent license under the contributor's essential patent claims, to +make, use, sell, offer for sale, import and otherwise run, modify and +propagate the contents of its contributor version. + + In the following three paragraphs, a "patent license" is any express +agreement or commitment, however denominated, not to enforce a patent +(such as an express permission to practice a patent or covenant not to +sue for patent infringement). To "grant" such a patent license to a +party means to make such an agreement or commitment not to enforce a +patent against the party. + + If you convey a covered work, knowingly relying on a patent license, +and the Corresponding Source of the work is not available for anyone +to copy, free of charge and under the terms of this License, through a +publicly available network server or other readily accessible means, +then you must either (1) cause the Corresponding Source to be so +available, or (2) arrange to deprive yourself of the benefit of the +patent license for this particular work, or (3) arrange, in a manner +consistent with the requirements of this License, to extend the patent +license to downstream recipients. "Knowingly relying" means you have +actual knowledge that, but for the patent license, your conveying the +covered work in a country, or your recipient's use of the covered work +in a country, would infringe one or more identifiable patents in that +country that you have reason to believe are valid. + + If, pursuant to or in connection with a single transaction or +arrangement, you convey, or propagate by procuring conveyance of, a +covered work, and grant a patent license to some of the parties +receiving the covered work authorizing them to use, propagate, modify +or convey a specific copy of the covered work, then the patent license +you grant is automatically extended to all recipients of the covered +work and works based on it. + + A patent license is "discriminatory" if it does not include within +the scope of its coverage, prohibits the exercise of, or is +conditioned on the non-exercise of one or more of the rights that are +specifically granted under this License. You may not convey a covered +work if you are a party to an arrangement with a third party that is +in the business of distributing software, under which you make payment +to the third party based on the extent of your activity of conveying +the work, and under which the third party grants, to any of the +parties who would receive the covered work from you, a discriminatory +patent license (a) in connection with copies of the covered work +conveyed by you (or copies made from those copies), or (b) primarily +for and in connection with specific products or compilations that +contain the covered work, unless you entered into that arrangement, +or that patent license was granted, prior to 28 March 2007. + + Nothing in this License shall be construed as excluding or limiting +any implied license or other defenses to infringement that may +otherwise be available to you under applicable patent law. + + 12. No Surrender of Others' Freedom. + + If 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 convey a +covered work so as to satisfy simultaneously your obligations under this +License and any other pertinent obligations, then as a consequence you may +not convey it at all. For example, if you agree to terms that obligate you +to collect a royalty for further conveying from those to whom you convey +the Program, the only way you could satisfy both those terms and this +License would be to refrain entirely from conveying the Program. + + 13. Remote Network Interaction; Use with the GNU General Public License. + + Notwithstanding any other provision of this License, if you modify the +Program, your modified version must prominently offer all users +interacting with it remotely through a computer network (if your version +supports such interaction) an opportunity to receive the Corresponding +Source of your version by providing access to the Corresponding Source +from a network server at no charge, through some standard or customary +means of facilitating copying of software. This Corresponding Source +shall include the Corresponding Source for any work covered by version 3 +of the GNU General Public License that is incorporated pursuant to the +following paragraph. + + Notwithstanding any other provision of this License, you have +permission to link or combine any covered work with a work licensed +under version 3 of the GNU General Public License into a single +combined work, and to convey the resulting work. The terms of this +License will continue to apply to the part which is the covered work, +but the work with which it is combined will remain governed by version +3 of the GNU General Public License. + + 14. Revised Versions of this License. + + The Free Software Foundation may publish revised and/or new versions of +the GNU Affero 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 that a certain numbered version of the GNU Affero General +Public License "or any later version" applies to it, you have the +option of following the terms and conditions either of that numbered +version or of any later version published by the Free Software +Foundation. If the Program does not specify a version number of the +GNU Affero General Public License, you may choose any version ever published +by the Free Software Foundation. + + If the Program specifies that a proxy can decide which future +versions of the GNU Affero General Public License can be used, that proxy's +public statement of acceptance of a version permanently authorizes you +to choose that version for the Program. + + Later license versions may give you additional or different +permissions. However, no additional obligations are imposed on any +author or copyright holder as a result of your choosing to follow a +later version. + + 15. Disclaimer of Warranty. + + 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. + + 16. Limitation of Liability. + + IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING +WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MODIFIES AND/OR CONVEYS +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. + + 17. Interpretation of Sections 15 and 16. + + If the disclaimer of warranty and limitation of liability provided +above cannot be given local legal effect according to their terms, +reviewing courts shall apply local law that most closely approximates +an absolute waiver of all civil liability in connection with the +Program, unless a warranty or assumption of liability accompanies a +copy of the Program in return for a fee. + + 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 +state 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 Affero General Public License as published + by the Free Software Foundation, either version 3 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 Affero General Public License for more details. + + You should have received a copy of the GNU Affero General Public License + along with this program. If not, see . + +Also add information on how to contact you by electronic and paper mail. + + If your software can interact with users remotely through a computer +network, you should also make sure that it provides a way for users to +get its source. For example, if your program is a web application, its +interface could display a "Source" link that leads users to an archive +of the code. There are many ways you could offer source, and different +solutions will be better for different programs; see section 13 for the +specific requirements. + + You should also get your employer (if you work as a programmer) or school, +if any, to sign a "copyright disclaimer" for the program, if necessary. +For more information on this, and how to apply and follow the GNU AGPL, see +. \ No newline at end of file diff --git a/README.md b/README.md index 769d8a698..d6ccee448 100644 --- a/README.md +++ b/README.md @@ -1,3 +1,11 @@ +# This is OpenDroneMap's OpenSfM fork + +This fork is pretty close to upstream OpenSfM:OpenSfM, but includes changes specific to ODM. + +Please note that all such changes are released under the AGPLv3 license, not BSD. + +--- + ![OpenSfM](doc/images/logo_small.png) @@ -137,6 +145,8 @@ Add `--dense` to include the dense stages, or `--resume` to recover an interrupt * [Benchmarking](doc/benchmark.md) ## ⚖️ License -OpenSfM is BSD-style licensed, as found in the LICENSE file. +Upstream OpenSfM is BSD-style licensed. + +> Changes specific to this OpenDroneMap fork are released under the AGPLv3 license, as found in the LICENSE file. Example data in the README is under [Creative Commons CC-BY 4.0 License](https://creativecommons.org/licenses/by/4.0/) by Wingtra AG, 8045 Zürich, Switzerland. \ No newline at end of file diff --git a/bin/rolling_shutter_plot b/bin/rolling_shutter_plot new file mode 100644 index 000000000..733cfdae5 --- /dev/null +++ b/bin/rolling_shutter_plot @@ -0,0 +1,97 @@ +#!/usr/bin/env python3 + +import argparse +import os.path + +import matplotlib.pyplot as plt +from matplotlib.collections import PatchCollection +from matplotlib.patches import Wedge + +from opensfm import dataset +from opensfm import features +from opensfm import io + + +def plot_features(image, points, orig_points): + h, w, d = image.shape + pixels = features.denormalized_image_coordinates(points, w, h) + pixels_orig = features.denormalized_image_coordinates(orig_points, w, h) + + points[:, 2] = 0.001 + sizes = points[:, 2] * max(w, h) + angles = points[:, 3] + + ax = plt.axes() + ax.imshow(image) + + patches = [] + for p, s, a in zip(pixels, sizes, angles): + patches.append(Wedge(p, s, a + 1, a - 1)) + + orig_points[:, 2] = 0.001 + sizes_orig = orig_points[:, 2] * max(w, h) + angles_orig = orig_points[:, 3] + orig_patches = [] + + for p, s, a in zip(pixels_orig, sizes_orig, angles_orig): + orig_patches.append(Wedge(p, s, a + 1, a - 1)) + + collection = PatchCollection(patches, alpha=0.5, edgecolor='b', facecolor='w') + orig_collection = PatchCollection(orig_patches, alpha=0.5, edgecolor='r', facecolor='w') + + ax.add_collection(orig_collection) + ax.add_collection(collection) + + +if __name__ == "__main__": + parser = argparse.ArgumentParser(description=""" + Plot detected features + NOTE that we require to backup the original features e.g.: + 1. Backup: cp -r path/to/dataset/features path/to/dataset/features.orig + 2. Correct: bin/opensfm correct_rolling_shutter path/to/dataset + 3. Plot: python bin/rolling_shutter_plot path/to/dataset --image image.jpg + """) + parser.add_argument('dataset', + help='path to the dataset to be processed') + parser.add_argument('--image', + help='name of the image to show') + parser.add_argument('--save_figs', + help='save figures instead of showing them', + action='store_true') + args = parser.parse_args() + + data = dataset.DataSet(args.dataset) + + images = [args.image] if args.image else data.images() + for image in images: + # Load the corrected features + features_data = data.load_features(image) + + # Manually load the original features from features.orig + features_orig_path = os.path.join(data.data_path, "features.orig", image + ".features.npz") + features_data_orig = None + if os.path.exists(features_orig_path): + with data.io_handler.open_rb(features_orig_path) as f: + features_data_orig = features.FeaturesData.from_file(f, data.config) + + if not features_data: + continue + if not features_data_orig: + print(f"Skipping {image} because features.orig not found. Did you backup features to features.orig before running correct_rolling_shutter?") + continue + + points = features_data.points + orig_points = features_data_orig.points + + print("plotting {0}/{1} points".format(len(points), len(orig_points))) + plt.figure() + plt.title('Image: ' + image) + fig = plot_features(data.load_image(image), points, orig_points) + + if args.save_figs: + p = os.path.join(args.dataset, 'plot_features') + io.mkdir_p(p) + plt.savefig(os.path.join(p, image + '.jpg'), dpi=1000) + plt.close() + else: + plt.show() diff --git a/opensfm/actions/export_bundler.py b/opensfm/actions/export_bundler.py index 8c9ca5c93..ca135da47 100644 --- a/opensfm/actions/export_bundler.py +++ b/opensfm/actions/export_bundler.py @@ -68,11 +68,7 @@ def export_bundler( if shot_id in shots: shot = shots[shot_id] camera = shot.camera - if shot.camera.projection_type == "brown": - # Will approximate Brown model, not optimal - focal_normalized = camera.focal_x - else: - focal_normalized = camera.focal + focal_normalized = camera.focal scale = max(camera.width, camera.height) focal = focal_normalized * scale k1 = camera.k1 diff --git a/opensfm/actions/export_geocoords.py b/opensfm/actions/export_geocoords.py index 34f6dc695..9be4a313a 100644 --- a/opensfm/actions/export_geocoords.py +++ b/opensfm/actions/export_geocoords.py @@ -9,6 +9,7 @@ from opensfm import io, types from opensfm import geo from opensfm.dataset import DataSet, UndistortedDataSet +from opensfm.reconstruction import bundle_shot_poses logger: logging.Logger = logging.getLogger(__name__) @@ -21,6 +22,8 @@ def run_dataset( reconstruction: bool, dense: bool, output: str, + offset=(0, 0), + mode="affine" ) -> None: """Export reconstructions in geographic coordinates @@ -31,6 +34,8 @@ def run_dataset( reconstruction : export reconstruction.json dense : export dense point cloud (depthmaps/merged.ply) output : path of the output file relative to the dataset + offset : offset to substract from the translation (optional) + mode : Method of georeferencing (optional) [affine|projected] """ @@ -42,6 +47,8 @@ def run_dataset( projection = geo.construct_proj_transformer(proj, inverse=True) t = geo.get_proj_transform_matrix(reference, projection) + t[0, 3] -= offset[0] + t[1, 3] -= offset[1] if transformation: output = output or "geocoords_transformation.txt" @@ -56,8 +63,28 @@ def run_dataset( if reconstruction: reconstructions = data.load_reconstruction() - for r in reconstructions: - geo.transform_reconstruction_with_proj(r, projection) + if mode == "affine": + for r in reconstructions: + geo.transform_reconstruction_with_matrix(r, t) + elif mode == "projected": + camera_priors = data.load_camera_models() + rig_camera_priors = data.load_rig_cameras() + for r in reconstructions: + geo.transform_reconstruction_with_proj(r, projection) + if offset != (0, 0): + for point in r.points.values(): + point.coordinates[0] -= offset[0] + point.coordinates[1] -= offset[1] + for rig_instance in r.rig_instances.values(): + o = rig_instance.pose.get_origin() + rig_instance.pose.set_origin( + [o[0] - offset[0], o[1] - offset[1], o[2]]) + + logger.info("Bundle shot poses") + bundle_shot_poses(r, set(r.shots.keys()), camera_priors, rig_camera_priors, data.config) + else: + raise Exception(f"Invalid mode: {mode}") + output = output or "reconstruction.geocoords.json" data.save_reconstruction(reconstructions, output) diff --git a/opensfm/actions/export_openmvs.py b/opensfm/actions/export_openmvs.py index 8600b18d2..5c78cf59b 100644 --- a/opensfm/actions/export_openmvs.py +++ b/opensfm/actions/export_openmvs.py @@ -45,6 +45,8 @@ def export( ) exporter.add_camera(str(camera.id), K, w, h) + shots_map = {} + for shot in reconstruction.shots.values(): if export_only is not None and shot.id not in export_only: continue @@ -55,6 +57,8 @@ def export( if not os.path.isfile(mask_path): mask_path = "" + shots_map[str(shot.id)] = shot + exporter.add_shot( str(os.path.abspath(image_path)), str(os.path.abspath(mask_path)), @@ -64,6 +68,15 @@ def export( shot.pose.get_origin(), ) + def positive_point_depth(point, shot_id): + if not shot_id in shots_map: + return False + + shot = shots_map[shot_id] + + # Is point in front of the camera? + return shot.pose.transform(point.coordinates)[2] > 0 + for point in reconstruction.points.values(): observations = tracks_manager.get_track_observations(point.id) @@ -72,6 +85,9 @@ def export( else: shots = list(observations) + if shots: + shots = [s for s in shots if positive_point_depth(point, s)] + if shots: coordinates = np.array(point.coordinates, dtype=np.float64) exporter.add_point(coordinates, shots) diff --git a/opensfm/actions/export_visualsfm.py b/opensfm/actions/export_visualsfm.py index 2f7defbdd..119767113 100644 --- a/opensfm/actions/export_visualsfm.py +++ b/opensfm/actions/export_visualsfm.py @@ -54,11 +54,7 @@ def export( shot_size_cache[shot.id] = udata.undistorted_image_size(shot.id) shot_index[shot.id] = i i += 1 - if shot.camera.projection_type == "brown": - # Will approximate Brown model, not optimal - focal_normalized = (shot.camera.focal_x + shot.camera.focal_y) / 2.0 - else: - focal_normalized = shot.camera.focal + focal_normalized = shot.camera.focal words = [ image_path(shot.id, udata), @@ -82,11 +78,11 @@ def export( skipped_points = 0 lines.append("") points = reconstruction.points + shots_keys = list(reconstruction.shots.keys()) lines.append(str(len(points))) points_count_index = len(lines) - 1 for point_id, point in points.items(): - shots = reconstruction.shots coord = point.coordinates color = list(map(int, point.color)) @@ -97,7 +93,7 @@ def export( if export_only is not None and shot_key not in export_only: continue - if shot_key in shots.keys(): + if shot_key in shots_keys: v = obs.point x = (0.5 + v[0]) * shot_size_cache[shot_key][1] y = (0.5 + v[1]) * shot_size_cache[shot_key][0] diff --git a/opensfm/actions/extract_metadata.py b/opensfm/actions/extract_metadata.py index b0476d831..f7e22d792 100644 --- a/opensfm/actions/extract_metadata.py +++ b/opensfm/actions/extract_metadata.py @@ -93,6 +93,9 @@ def _extract_exif(image: str, data: DataSetBase) -> Dict[str, Any]: ): d["model"] = f"unknown_{image}" + if data.config.get("camera_projection_type") != 'AUTO': + d['projection_type'] = data.config['camera_projection_type'].lower() + d["camera"] = exif.camera_id(d) return d diff --git a/opensfm/actions/reconstruct.py b/opensfm/actions/reconstruct.py index 0e96be79d..95385e55f 100644 --- a/opensfm/actions/reconstruct.py +++ b/opensfm/actions/reconstruct.py @@ -1,16 +1,26 @@ # pyre-strict import gc +from typing import Optional from opensfm import io, reconstruction from opensfm.dataset_base import DataSetBase def run_dataset( - data: DataSetBase, algorithm: reconstruction.ReconstructionAlgorithm + data: DataSetBase, algorithm: Optional[reconstruction.ReconstructionAlgorithm] = None ) -> None: """Compute the SfM reconstruction.""" tracks_manager = data.load_tracks_manager() + if algorithm is None: + algo_str = data.config.get("reconstruction_algorithm", "incremental").lower() + try: + algorithm = reconstruction.ReconstructionAlgorithm(algo_str.lower()) + except ValueError: + raise RuntimeError( + f"Unsupported algorithm for reconstruction '{algo_str}'. Allowed values: {list(v.value for v in reconstruction.ReconstructionAlgorithm)}" + ) + if algorithm == reconstruction.ReconstructionAlgorithm.INCREMENTAL: report, reconstructions = reconstruction.incremental_reconstruction( data, tracks_manager @@ -19,6 +29,10 @@ def run_dataset( report, reconstructions = reconstruction.triangulation_reconstruction( data, tracks_manager ) + elif algorithm == reconstruction.ReconstructionAlgorithm.PLANAR: + report, reconstructions = reconstruction.planar_reconstruction( + data, tracks_manager + ) else: raise RuntimeError( f"Unsupported algorithm for reconstruction {algorithm}") diff --git a/opensfm/actions/undistort.py b/opensfm/actions/undistort.py index 13fae2f95..bbf5015c0 100644 --- a/opensfm/actions/undistort.py +++ b/opensfm/actions/undistort.py @@ -1,6 +1,7 @@ # pyre-strict import os -from typing import Optional +from typing import Callable, Optional +import numpy as np from opensfm import dataset, undistort from opensfm.dataset import DataSet @@ -12,6 +13,7 @@ def run_dataset( reconstruction_index: int = 0, tracks: Optional[str] = None, output: str = "undistorted", + imageFilter: Optional[Callable[[str, np.ndarray], np.ndarray]] = None, skip_images: bool = False, ) -> None: """Export reconstruction to NVM_V3 format from VisualSfM @@ -36,5 +38,5 @@ def run_dataset( if reconstructions: r = reconstructions[reconstruction_index] undistort.undistort_reconstruction_with_images( - tracks_manager, r, data, udata, skip_images + tracks_manager, r, data, udata, imageFilter, skip_images ) diff --git a/opensfm/commands/export_geocoords.py b/opensfm/commands/export_geocoords.py index f8b05b9bc..c4f9acc53 100644 --- a/opensfm/commands/export_geocoords.py +++ b/opensfm/commands/export_geocoords.py @@ -20,6 +20,8 @@ def run_impl(self, dataset: DataSet, args: argparse.Namespace) -> None: args.reconstruction, args.dense, args.output, + (args.offset_x, args.offset_y), + args.mode, ) def add_arguments_impl(self, parser: argparse.ArgumentParser) -> None: @@ -51,3 +53,16 @@ def add_arguments_impl(self, parser: argparse.ArgumentParser) -> None: parser.add_argument( "--output", help="Path of the output file relative to the dataset" ) + parser.add_argument( + "--offset-x", type=float, help="Value to add to the final translation X axis", default=0.0 + ) + parser.add_argument( + "--offset-y", type=float, help="Value to add to the final translation Y axis", default=0.0 + ) + parser.add_argument( + "--mode", + help="Georeferencing method", + type=str, + choices=["affine", "projected"], + default="affine", + ) diff --git a/opensfm/commands/reconstruct.py b/opensfm/commands/reconstruct.py index 0df0030f1..de47d4814 100644 --- a/opensfm/commands/reconstruct.py +++ b/opensfm/commands/reconstruct.py @@ -21,5 +21,5 @@ def add_arguments_impl(self, parser: argparse.ArgumentParser) -> None: help="SfM algorithm to use to run reconstrution", type=str, choices=[k.value for k in reconstruction.ReconstructionAlgorithm], - default=reconstruction.ReconstructionAlgorithm.INCREMENTAL.value, + default=None, ) diff --git a/opensfm/commands/undistort.py b/opensfm/commands/undistort.py index 51125380e..16aff412b 100644 --- a/opensfm/commands/undistort.py +++ b/opensfm/commands/undistort.py @@ -18,7 +18,7 @@ def run_impl(self, dataset: DataSet, args: argparse.Namespace) -> None: args.reconstruction_index, args.tracks, args.output, - args.skip_images, + skip_images=args.skip_images, ) def add_arguments_impl(self, parser: argparse.ArgumentParser) -> None: diff --git a/opensfm/config.py b/opensfm/config.py index 6754b5b23..c0f677dba 100644 --- a/opensfm/config.py +++ b/opensfm/config.py @@ -18,7 +18,7 @@ class OpenSfMConfig: # Treat images from unknown camera models as coming from different cameras unknown_camera_models_are_different: bool = False # Default projection type to use when it cannot be inferred from EXIF metadata - default_projection_type: str = "perspective" + default_projection_type: str = "brown" # Default focal length to sensor size ratio to use when it cannot be inferred from EXIF metadata default_focal_prior: float = 0.85 # Enable downloading missing datum grids from the PROJ CDN @@ -53,7 +53,7 @@ class OpenSfMConfig: # Params for SIFT ################################## # Smaller value -> more features - sift_peak_threshold: float = 0.1 + sift_peak_threshold: float = 0.2 # See OpenCV doc sift_edge_threshold: int = 10 # See OpenCV doc @@ -663,6 +663,15 @@ class OpenSfMConfig: # Language for the quality report: "en", "fr", "es", "de", "it" report_language: str = "en" + #################################### + # ODM specific flags + #################################### + # The projection type of the camera : attempt to detect it from metadata (AUTO), + # or set it manually (PERSPECTIVE, BROWN, FISHEYE, SPHERICAL) + camera_projection_type: str = "AUTO" + # SfM algorithm to use to run reconstruction (incremental, triangulation) + reconstruction_algorithm: str = "incremental" + def default_config() -> Dict[str, Any]: """Return default configuration""" diff --git a/opensfm/context.py b/opensfm/context.py index 570cf638c..e39c0282c 100644 --- a/opensfm/context.py +++ b/opensfm/context.py @@ -18,6 +18,7 @@ abspath: str = os.path.dirname(os.path.realpath(__file__)) SENSOR_DATA: str = os.path.join(abspath, "data", "sensor_data.json") +SENSOR_DATA_DB: str = os.path.join(abspath, "data", "sensor_data.sqlite") CAMERA_CALIBRATION: str = os.path.join( abspath, "data", "camera_calibration.yaml") BOW_PATH: str = os.path.join(abspath, "data", "bow") diff --git a/opensfm/data/locale/en.yaml b/opensfm/data/locale/en.yaml index c143c8d4a..d8c5f753f 100644 --- a/opensfm/data/locale/en.yaml +++ b/opensfm/data/locale/en.yaml @@ -24,10 +24,14 @@ dataset: "Dataset" date: "Date" area_covered: "Area Covered" processing_time: "Processing Time" +capture_start: "Capture Start" +capture_end: "Capture End" # Processing summary row labels reconstructed_images: "Reconstructed Images" reconstructed_points: "Reconstructed Points" +reconstructed_points_sparse: "Reconstructed Points (Sparse)" +reconstructed_points_dense: "Reconstructed Points (Dense)" average_track_length: "Average Track Length" camera_calibration_deviation: "Camera Calibration Deviation" front_overlap: "Front Overlap" @@ -38,6 +42,8 @@ geographic_referencing: "Geographic Referencing" gcp_coordinate_system: "GCP Coordinate System" gps_error: "GPS Error" gcp_error: "GCP Error" +alignment: "Alignment" +alignment_errors: "Alignment Errors" # Reconstruction details reprojection_error: "Reprojection Error (norm / px / angular)" diff --git a/opensfm/data/sensor_data.sqlite b/opensfm/data/sensor_data.sqlite new file mode 100644 index 000000000..eae120962 Binary files /dev/null and b/opensfm/data/sensor_data.sqlite differ diff --git a/opensfm/dataset.py b/opensfm/dataset.py index 652597e45..ae88572d1 100644 --- a/opensfm/dataset.py +++ b/opensfm/dataset.py @@ -1,5 +1,6 @@ # pyre-strict import gzip +import hashlib import json import logging import os @@ -418,7 +419,11 @@ def tracks_exists(self, filename: Optional[str] = None) -> bool: def save_tracks_manager( self, tracks_manager: pymap.TracksManager, filename: Optional[str] = None ) -> None: - with self.io_handler.open_wt(self._tracks_manager_file(filename)) as fw: + path = self._tracks_manager_file(filename) + if isinstance(self.io_handler, io.IoFilesystemDefault): + tracks_manager.write_to_file(path) + return + with self.io_handler.open_wt(path) as fw: fw.write(tracks_manager.as_string()) def _reconstruction_file(self, filename: Optional[str]) -> str: @@ -824,6 +829,10 @@ def _undistorted_image_path(self) -> str: def _undistorted_image_file(self, image: str) -> str: """Path of undistorted version of an image.""" + p, ext = os.path.splitext(image) + if ' ' in image: + # NOTE: Downstream tools in ODM (ie OpenMVS) choke on spaces in image paths. + image = p.replace(' ', '') + hashlib.md5(image.encode('utf8')).hexdigest() + ext return os.path.join(self._undistorted_image_path(), image) def load_undistorted_image(self, image: str) -> NDArray: @@ -1630,6 +1639,8 @@ def load_fusion_progress(self) -> Tuple[set, List[int]]: def load_undistorted_tracks_manager(self) -> pymap.TracksManager: filename = os.path.join(self.data_path, "tracks.csv") + if isinstance(self.io_handler, io.IoFilesystemDefault): + return pymap.TracksManager.instanciate_from_file(filename) with self.io_handler.open_rt(filename) as f: return pymap.TracksManager.instanciate_from_string(f.read()) @@ -1637,6 +1648,9 @@ def save_undistorted_tracks_manager( self, tracks_manager: pymap.TracksManager ) -> None: filename = os.path.join(self.data_path, "tracks.csv") + if isinstance(self.io_handler, io.IoFilesystemDefault): + tracks_manager.write_to_file(filename) + return with self.io_handler.open_wt(filename) as fw: fw.write(tracks_manager.as_string()) diff --git a/opensfm/exif.py b/opensfm/exif.py index 958482449..4695651c1 100644 --- a/opensfm/exif.py +++ b/opensfm/exif.py @@ -1,8 +1,8 @@ # pyre-strict - import datetime import logging from codecs import decode, encode +from bs4 import BeautifulSoup from typing import Any, BinaryIO, Callable, Dict, List, Optional, Tuple, Union import exifread @@ -142,11 +142,18 @@ def unescape_string(s: str) -> str: def parse_xmp_string(xmp_str: str) -> Optional[Dict[str, Any]]: - for _ in range(2): + try: + return x2d.parse(xmp_str) + except Exception: + xmp_str = unescape_string(xmp_str) try: return x2d.parse(xmp_str) except Exception: - xmp_str = unescape_string(xmp_str) + xmp_str = str(BeautifulSoup(xmp_str, "xml")) + try: + return x2d.parse(xmp_str) + except Exception: + pass return None @@ -180,6 +187,14 @@ def get_gpano_from_xmp(xmp: List[Dict[str, Any]]) -> Dict[str, Any]: return {} +def get_pix4d_from_xmp(xmp: List[Dict[str, Any]]) -> Dict[str, Any]: + for i in xmp: + for k in i: + if "@Camera:ModelType" in k: + return i + return {} + + class EXIF: def __init__( self, @@ -254,8 +269,19 @@ def extract_model(self) -> str: return self._decode_make_model(model) def extract_projection_type(self) -> str: + projections = ["perspective", "fisheye", "brown", "dual", "equirectangular", "spherical"] + gpano = get_gpano_from_xmp(self.xmp) - return gpano.get("GPano:ProjectionType", self.default_projection_type) + gpano_projection = gpano.get("GPano:ProjectionType", "").lower() + if gpano_projection in projections: + return gpano_projection + + pix4d = get_pix4d_from_xmp(self.xmp) + camera_model = pix4d.get("@Camera:ModelType", "").lower() + if camera_model in projections: + return camera_model + + return self.default_projection_type def extract_focal(self) -> Tuple[float, float]: # DJI CalibratedFocalLength is focal length in pixels — use it if available @@ -696,6 +722,10 @@ def extract_exif(self) -> Dict[str, Any]: def hard_coded_calibration(exif: Dict[str, Any]) -> Optional[Dict[str, Any]]: + # Disabled: ODM relies on EXIF/sensor-DB calibration, not on hard-coded + # per-camera values (ported from 360 branch, commit dbe22158) + return None + focal = exif["focal_ratio"] fmm35 = int(round(focal * 36.0)) make = exif["make"].strip().lower() diff --git a/opensfm/features.py b/opensfm/features.py index 096e13515..aaa551c51 100644 --- a/opensfm/features.py +++ b/opensfm/features.py @@ -365,6 +365,8 @@ def extract_features_sift( sift_edge_threshold = config["sift_edge_threshold"] sift_peak_threshold = float(config["sift_peak_threshold"]) sift_nfeatures = config["sift_nfeatures"] + if sift_nfeatures <= 0: + sift_nfeatures = int(features_count * 3 / 2) sift_octave_layers = config["sift_octave_layers"] sift_sigma = float(config["sift_sigma"]) while True: @@ -393,6 +395,8 @@ def extract_features_sift( detector = cv2.FeatureDetector_create("SIFT") descriptor = cv2.DescriptorExtractor_create("SIFT") detector.setDouble("edgeThreshold", sift_edge_threshold) + detector.setDouble("contrastThreshold", sift_peak_threshold) + detector.setDouble("nFeatures", sift_nfeatures) points = detector.detect(image) logger.debug("Found {0} points in {1}s".format(len(points), time.time() - t)) @@ -415,6 +419,22 @@ def extract_features_sift( return points, desc +def extract_features_popsift( + image: np.ndarray, config: Dict[str, Any], features_count: int +) -> Tuple[np.ndarray, np.ndarray]: + from opensfm import pypopsift + + sift_edge_threshold = float(config["sift_edge_threshold"]) + sift_peak_threshold = float(config["sift_peak_threshold"]) + + points, desc = pypopsift.popsift(image, peak_threshold=sift_peak_threshold, + edge_threshold=sift_edge_threshold, + target_num_features=features_count, + use_root=bool(config["feature_root"])) + + return points, desc + + def extract_features_surf( image: NDArray, config: Dict[str, Any], features_count: int ) -> Tuple[NDArray, NDArray]: @@ -634,11 +654,13 @@ def extract_features( points, desc = extract_features_hahog(image_gray, config, features_count) elif feature_type == "ORB": points, desc = extract_features_orb(image_gray, config, features_count) + elif feature_type == 'SIFT_GPU': + points, desc = extract_features_popsift(image_gray, config, features_count) elif feature_type == 'DSPSIFT': points, desc = extract_features_dspsift(image_gray, config, features_count) else: raise ValueError( - "Unknown feature type (must be SURF, SIFT, AKAZE, HAHOG, ORB or DSPSIFT)" + "Unknown feature type (must be SURF, SIFT, AKAZE, HAHOG, ORB, SIFT_GPU or DSPSIFT)" ) xs = points[:, 0].round().astype(int) diff --git a/opensfm/features_processing.py b/opensfm/features_processing.py index 7927a68fd..2147d0c48 100644 --- a/opensfm/features_processing.py +++ b/opensfm/features_processing.py @@ -188,7 +188,7 @@ def read_images( expected: int, force: bool, ) -> None: - full_queue_timeout = 600 + full_queue_timeout = 3600 for image in images: logger.info( f"Reading data for image {image} (queue-size={queue.qsize()})") diff --git a/opensfm/geo.py b/opensfm/geo.py index 06ff8836d..764c22929 100644 --- a/opensfm/geo.py +++ b/opensfm/geo.py @@ -539,10 +539,10 @@ def transform_reconstruction_with_matrix( A, b = transformation[:3, :3], transformation[:3, 3] A1 = np.linalg.inv(A) - for shot in reconstruction.shots.values(): - R = shot.pose.get_rotation_matrix() - shot.pose.set_rotation_matrix(np.dot(R, A1)) - shot.pose.set_origin(np.dot(A, shot.pose.get_origin()) + b) + for rig_instance in reconstruction.rig_instances.values(): + R = rig_instance.pose.get_rotation_matrix() + rig_instance.pose.set_rotation_matrix(np.dot(R, A1)) + rig_instance.pose.set_origin(np.dot(A, rig_instance.pose.get_origin()) + b) for point in reconstruction.points.values(): point.coordinates = list(np.dot(A, point.coordinates) + b) diff --git a/opensfm/io.py b/opensfm/io.py index ca12a61f5..b9d04fdf5 100644 --- a/opensfm/io.py +++ b/opensfm/io.py @@ -23,9 +23,17 @@ import pyproj from numpy.typing import NDArray from opensfm import context, features, geo, pygeometry, pymap, types +import rasterio +import rawpy +from rasterio.plot import reshape_as_image +import warnings +warnings.filterwarnings("ignore", category=rasterio.errors.NotGeoreferencedWarning) + from PIL import Image logger: logging.Logger = logging.getLogger(__name__) +logging.getLogger("rasterio").setLevel(logging.WARNING) +logging.getLogger("PIL").setLevel(logging.WARNING) JSONType = Any # pyre-ignore[33] @@ -1320,8 +1328,14 @@ def open_rt(path: str) -> TextIO: def imread( path: str, grayscale: bool = False, unchanged: bool = False, anydepth: bool = False ) -> NDArray: - with open(path, "rb") as fb: - return imread_from_fileobject(fb, grayscale, unchanged, anydepth) + _, ext = os.path.splitext(path) + if ext.lower() == ".tiff" or ext.lower() == ".tif": + return imread_rasterio(path, grayscale, unchanged, anydepth) + elif ext.lower() in [".dng", ".raw", ".nef"]: + return imread_rawpy(path, grayscale, unchanged, anydepth) + else: + with open(path, "rb") as fb: + return imread_from_fileobject(fb, grayscale, unchanged, anydepth) def imread_from_fileobject( @@ -1370,12 +1384,89 @@ def imread_from_fileobject( if len(image.shape) == 3: # Turn BGR to RGB (or BGRA to RGBA) image[:, :, :3] = image[:, :, [2, 1, 0]] + elif len(image.shape) == 2: + image = image[..., np.newaxis] # Make sure we always have a band dimension + return image + + +def imread_rasterio(path, grayscale=False, unchanged=False, anydepth=False): + """Load image as an array ignoring EXIF orientation.""" + if grayscale: + raise IOError("Grayscale not implemented") + + with rasterio.open(path, "r") as f: + image = reshape_as_image(f.read()) + + return _imread_postprocess(image, grayscale, unchanged, anydepth) + + +def imread_rawpy(path, grayscale=False, unchanged=False, anydepth=False): + if grayscale: + raise IOError("Grayscale not implemented") + + with rawpy.imread(path) as r: + image = r.postprocess(output_bps=16, use_camera_wb=True, use_auto_wb=False) + + return _imread_postprocess(image, grayscale, unchanged, anydepth) + + +def _imread_postprocess(image, grayscale=False, unchanged=False, anydepth=False): + if image is None: + raise IOError("Unable to load image") + + if not anydepth or not unchanged: + # Normalize to 8bit + min_value = float(image.min()) + value_range = float(image.max()) - min_value + + image = image.astype(np.float32) + image -= min_value + image *= 255.0 / value_range + np.around(image, out=image) + image[image > 255] = 255 + image[image < 0] = 0 + image = image.astype(np.uint8) + + if not unchanged: + # Convert to RGB + if image.shape[2] == 1: + image = np.repeat(image[:, :, :], 3, axis=2) + elif image.shape[2] > 3: + image = image[:,:,:3] + return image def imwrite(path: str, image: NDArray) -> None: - with open(path, "wb") as fwb: - return imwrite_from_fileobject(fwb, image, path) + _, ext = os.path.splitext(path) + if ext.lower() == ".tiff" or ext.lower() == ".tif": + return imwrite_rasterio(path, image) + else: + with open(path, "wb") as fwb: + return imwrite_from_fileobject(fwb, image, path) + + +def imwrite_rasterio(path: str, image: NDArray) -> None: + if len(image.shape) == 2: + image = np.repeat(image[:, :, np.newaxis], 1, axis=2) + + profile = { + 'width': image.shape[1], + 'height': image.shape[0], + 'count': image.shape[2], + 'blockxsize': 512, + 'blockysize': 512, + 'compress': 'lzw', + 'predictor': 3 if image.dtype == np.float32 else 2, + 'driver': 'GTiff', + 'dtype': image.dtype, + 'interleave': 'pixel', + 'tiled': False + } + + with rasterio.open(path, "w", **profile) as f: + for b in range(0, image.shape[2]): + f.write(image[:, :, b], b + 1) def imwrite_from_fileobject( @@ -1547,18 +1638,33 @@ def imread( unchanged: bool = False, anydepth: bool = False, ) -> NDArray: - with cls.open_rb(path) as fb: - return imread_from_fileobject(fb, grayscale, unchanged, anydepth) + _, ext = os.path.splitext(path) + if ext.lower() == ".tiff" or ext.lower() == ".tif": + return imread_rasterio(path, grayscale, unchanged, anydepth) + elif ext.lower() in [".dng", ".raw", ".nef"]: + return imread_rawpy(path, grayscale, unchanged, anydepth) + else: + with cls.open_rb(path) as fb: + return imread_from_fileobject(fb, grayscale, unchanged, anydepth) @classmethod def imwrite(cls, path: str, image: NDArray) -> None: - with cls.open_wb(path) as fwb: - imwrite_from_fileobject(fwb, image, path) + _, ext = os.path.splitext(path) + if ext.lower() == ".tiff" or ext.lower() == ".tif": + imwrite_rasterio(path, image) + else: + with cls.open_wb(path) as fwb: + imwrite_from_fileobject(fwb, image, path) @classmethod def image_size(cls, path: str) -> Tuple[int, int]: - with cls.open_rb(path) as fb: - return image_size_from_fileobject(fb) + try: + with cls.open_rb(path) as fb: + return image_size_from_fileobject(fb) + except Exception: + # Fallback to rasterio (RGB 32bit floats fail with PIL) + with rasterio.open(path, "r") as r: + return r.height, r.width @classmethod def timestamp(cls, path: str) -> float: diff --git a/opensfm/large/lru_cache.py b/opensfm/large/lru_cache.py new file mode 100644 index 000000000..ac9dc544c --- /dev/null +++ b/opensfm/large/lru_cache.py @@ -0,0 +1,272 @@ +#!/usr/bin/env python2 +# -*- coding: utf-8 -*- + +""" +Memory-aware LRU Cache function decorator +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + +A modification of the builtin ``functools.lru_cache`` decorator that takes an +additional keyword argument, ``use_memory_up_to``. The cache is considered full +if there are fewer than ``use_memory_up_to`` bytes of memory available. + +If ``use_memory_up_to`` is set, then ``maxsize`` has no effect. + +Uses the ``vmem`` module to get the available memory. +""" + +import vmem +try: + from functools import RLock, namedtuple +except ImportError: + from threading import RLock + from collections import namedtuple +from functools import update_wrapper + +_CacheInfo = namedtuple("CacheInfo", ["hits", "misses", "maxsize", "currsize"]) + +class _HashedSeq(list): + """ This class guarantees that hash() will be called no more than once + per element. This is important because the lru_cache() will hash + the key multiple times on a cache miss. + + """ + + __slots__ = 'hashvalue' + + def __init__(self, tup, hash=hash): + self[:] = tup + self.hashvalue = hash(tup) + + def __hash__(self): + return self.hashvalue + +def _make_key(args, kwds, typed, + kwd_mark = (object(),), + fasttypes = {int, str, frozenset, type(None)}, + sorted=sorted, tuple=tuple, type=type, len=len): + """Make a cache key from optionally typed positional and keyword arguments + + The key is constructed in a way that is flat as possible rather than + as a nested structure that would take more memory. + + If there is only a single argument and its data type is known to cache + its hash value, then that argument is returned without a wrapper. This + saves space and improves lookup speed. + + """ + key = args + if kwds: + sorted_items = sorted(kwds.items()) + key += kwd_mark + for item in sorted_items: + key += item + if typed: + key += tuple(type(v) for v in args) + if kwds: + key += tuple(type(v) for k, v in sorted_items) + elif len(key) == 1 and type(key[0]) in fasttypes: + return key[0] + return _HashedSeq(key) + +def lru_cache(maxsize=128, typed=False, use_memory_up_to=False): + """Least-recently-used cache decorator. + + *use_memory_up_to* is an integer representing the number of bytes of memory + that must be available on the system in order for a value to be cached. If + it is set, *maxsize* has no effect. + + If *maxsize* is set to None, the LRU features are disabled and the cache + can grow without bound. + + If *typed* is True, arguments of different types will be cached separately. + For example, f(3.0) and f(3) will be treated as distinct calls with + distinct results. + + Arguments to the cached function must be hashable. + + View the cache statistics named tuple (hits, misses, maxsize, currsize) + with f.cache_info(). Clear the cache and statistics with f.cache_clear(). + Access the underlying function with f.__wrapped__. + + See: http://en.wikipedia.org/wiki/Cache_algorithms#Least_Recently_Used + + """ + if use_memory_up_to: + maxsize=None + + # Users should only access the lru_cache through its public API: + # cache_info, cache_clear, and f.__wrapped__ + # The internals of the lru_cache are encapsulated for thread safety and + # to allow the implementation to change (including a possible C version). + + # Constants shared by all lru cache instances: + sentinel = object() # unique object used to signal cache misses + make_key = _make_key # build a key from the function arguments + PREV, NEXT, KEY, RESULT = 0, 1, 2, 3 # names for the link fields + + def decorating_function(user_function): + cache = {} + + class nonloc: + hits = 0 + misses = 0 + root = [] # root of the circular doubly linked list + full = False + + cache_get = cache.get # bound method to lookup a key or return None + lock = RLock() # because linkedlist updates aren't threadsafe + nonloc.root[:] = [nonloc.root, nonloc.root, None, None] # initialize by pointing to self + + if use_memory_up_to: + + def wrapper(*args, **kwds): + # Size limited caching that tracks accesses by recency + key = make_key(args, kwds, typed) + with lock: + link = cache_get(key) + if link is not None: + # Move the link to the front of the circular queue + link_prev, link_next, _key, result = link + link_prev[NEXT] = link_next + link_next[PREV] = link_prev + last = nonloc.root[PREV] + last[NEXT] = nonloc.root[PREV] = link + link[PREV] = last + link[NEXT] = nonloc.root + nonloc.hits += 1 + return result + result = user_function(*args, **kwds) + with lock: + if key in cache: + # Getting here means that this same key was added to the + # cache while the lock was released. Since the link + # update is already done, we need only return the + # computed result and update the count of misses. + pass + elif nonloc.full: + # Use the old root to store the new key and result. + oldroot = nonloc.root + oldroot[KEY] = key + oldroot[RESULT] = result + # Empty the oldest link and make it the new root. + # Keep a reference to the old key and old result to + # prevent their ref counts from going to zero during the + # update. That will prevent potentially arbitrary object + # clean-up code (i.e. __del__) from running while we're + # still adjusting the links. + nonloc.root = oldroot[NEXT] + oldkey = nonloc.root[KEY] + oldresult = nonloc.root[RESULT] + nonloc.root[KEY] = nonloc.root[RESULT] = None + # Now update the cache dictionary. + del cache[oldkey] + # Save the potentially reentrant cache[key] assignment + # for last, after the root and links have been put in + # a consistent state. + cache[key] = oldroot + else: + # Put result in a new link at the front of the queue. + last = nonloc.root[PREV] + link = [last, nonloc.root, key, result] + last[NEXT] = nonloc.root[PREV] = cache[key] = link + nonloc.full = (vmem.virtual_memory().available + < use_memory_up_to) + nonloc.misses += 1 + return result + + elif maxsize == 0: + + def wrapper(*args, **kwds): + # No caching -- just a statistics update after a successful call + result = user_function(*args, **kwds) + nonloc.misses += 1 + return result + + elif maxsize is None: + + def wrapper(*args, **kwds): + # Simple caching without ordering or size limit + key = make_key(args, kwds, typed) + result = cache_get(key, sentinel) + if result is not sentinel: + nonloc.hits += 1 + return result + result = user_function(*args, **kwds) + cache[key] = result + nonloc.misses += 1 + return result + + else: + + def wrapper(*args, **kwds): + # Size limited caching that tracks accesses by recency + key = make_key(args, kwds, typed) + with lock: + link = cache_get(key) + if link is not None: + # Move the link to the front of the circular queue + link_prev, link_next, _key, result = link + link_prev[NEXT] = link_next + link_next[PREV] = link_prev + last = nonloc.root[PREV] + last[NEXT] = nonloc.root[PREV] = link + link[PREV] = last + link[NEXT] = nonloc.root + nonloc.hits += 1 + return result + result = user_function(*args, **kwds) + with lock: + if key in cache: + # Getting here means that this same key was added to the + # cache while the lock was released. Since the link + # update is already done, we need only return the + # computed result and update the count of misses. + pass + elif nonloc.full: + # Use the old root to store the new key and result. + oldroot = nonloc.root + oldroot[KEY] = key + oldroot[RESULT] = result + # Empty the oldest link and make it the new root. + # Keep a reference to the old key and old result to + # prevent their ref counts from going to zero during the + # update. That will prevent potentially arbitrary object + # clean-up code (i.e. __del__) from running while we're + # still adjusting the links. + nonloc.root = oldroot[NEXT] + oldkey = nonloc.root[KEY] + oldresult = nonloc.root[RESULT] + nonloc.root[KEY] = nonloc.root[RESULT] = None + # Now update the cache dictionary. + del cache[oldkey] + # Save the potentially reentrant cache[key] assignment + # for last, after the root and links have been put in + # a consistent state. + cache[key] = oldroot + else: + # Put result in a new link at the front of the queue. + last = nonloc.root[PREV] + link = [last, nonloc.root, key, result] + last[NEXT] = nonloc.root[PREV] = cache[key] = link + nonloc.full = (len(cache) >= maxsize) + nonloc.misses += 1 + return result + + def cache_info(): + """Report cache statistics""" + with lock: + return _CacheInfo(nonloc.hits, misses, maxsize, len(cache)) + + def cache_clear(): + """Clear the cache and cache statistics""" + with lock: + cache.clear() + nonloc.root[:] = [nonloc.root, nonloc.root, None, None] + nonloc.hits = nonloc.misses = 0 + nonloc.full = False + + wrapper.cache_info = cache_info + wrapper.cache_clear = cache_clear + return update_wrapper(wrapper, user_function) + + return decorating_function \ No newline at end of file diff --git a/opensfm/large/tools.py b/opensfm/large/tools.py index 57eba9500..f44e09223 100644 --- a/opensfm/large/tools.py +++ b/opensfm/large/tools.py @@ -25,6 +25,8 @@ types, ) from opensfm.large.metadataset import MetaDataSet +import vmem +from opensfm.large.lru_cache import lru_cache logger: logging.Logger = logging.getLogger(__name__) @@ -205,7 +207,7 @@ def add_camera_constraints_hard( ) -@lru_cache(25) +@lru_cache(use_memory_up_to=vmem.virtual_memory().available * 0.9) def load_reconstruction( path: str, index: int ) -> Tuple[str, Tuple[types.Reconstruction, pymap.TracksManager]]: diff --git a/opensfm/log.py b/opensfm/log.py index 6ee33c3dc..4337f512c 100644 --- a/opensfm/log.py +++ b/opensfm/log.py @@ -3,20 +3,16 @@ import os from typing import Optional +import vmem + def setup() -> None: logging.basicConfig( - format="%(asctime)s %(levelname)s: %(message)s", level=logging.DEBUG + format="%(asctime)s %(levelname)s: %(message)s", level=logging.DEBUG, force=True ) def memory_available() -> Optional[int]: """Available memory in MB. - - Only works on linux and returns None otherwise. """ - lines = os.popen("free -t -m").readlines() - if not lines: - return None - available_mem = int(lines[1].split()[6]) - return available_mem + return vmem.virtual_memory().available / 1024 / 1024 diff --git a/opensfm/masking.py b/opensfm/masking.py index 8fc98f43f..61449ddea 100644 --- a/opensfm/masking.py +++ b/opensfm/masking.py @@ -66,6 +66,9 @@ def load_features_mask( if points is None or len(points) == 0: return np.array([], dtype=bool) + # ODM fork hack: we never mask features (they can add good keypoints) + return np.ones((points.shape[0],), dtype=bool) + if mask_image is None: mask_image = _load_combined_mask(data, image) if mask_image is None: diff --git a/opensfm/reconstruction.py b/opensfm/reconstruction.py index f920911b8..de3d49e7e 100644 --- a/opensfm/reconstruction.py +++ b/opensfm/reconstruction.py @@ -6,11 +6,12 @@ import logging import math from abc import ABC, abstractmethod -from collections import defaultdict +from collections import defaultdict, deque from itertools import combinations from timeit import default_timer as timer from typing import Any, Dict, List, Optional, Set, Tuple, Union - +import random +import networkx as nx import cv2 import numpy as np from numpy.typing import NDArray @@ -23,6 +24,7 @@ pysfm, reconstruction_helpers as helpers, rig, + tracking, types, ) from opensfm.align import align_reconstruction, apply_similarity @@ -33,9 +35,12 @@ logger: logging.Logger = logging.getLogger(__name__) +CommonTracksIndex = Dict[Tuple[str, str], Tuple[int, int, int]] + class ReconstructionAlgorithm(str, enum.Enum): INCREMENTAL = "incremental" TRIANGULATION = "triangulation" + PLANAR = "planar" def log_bundle_stats(bundle_type: str, bundle_report: Dict[str, Any]) -> None: @@ -244,8 +249,7 @@ def two_view_reconstruction_plane_based( b1, b2, R.T, -R.T.dot(t), threshold) motion_inliers.append(inliers) - # pyre-fixme[6]: For 1st argument expected `Union[_SupportsArray[dtype[typing.Any... - best = np.argmax(map(len, motion_inliers)) + best = np.argmax(list(map(len, motion_inliers))) R, t, n, d = motions[best] inliers = motion_inliers[best] return cv2.Rodrigues(R)[0].ravel(), t, inliers @@ -1275,3 +1279,423 @@ def lap_times(self) -> List[Tuple[str, float]]: def total_time(self) -> float: return self.laps[-1][2] - self.laps[0][2] + +def triangulate_planar(tracks_manager: pymap.TracksManager, reconstruction: types.Reconstruction, track: str, threshold: float) -> None: + """Triangulate track using a main plane and add point to reconstruction.""" + os_list, bs, ids = [], [], [] + + plane_center = np.array([0, 0, 1]) + plane_normal = np.array([0, 0, 1]) + + observations = tracks_manager.get_track_observations(track) + for shot_id, obs in observations.items(): + if shot_id in reconstruction.shots: + shot = reconstruction.shots[shot_id] + os_list.append(shot.pose.get_origin()) + b = shot.camera.pixel_bearing(obs.point) + r = shot.pose.get_rotation_matrix().T + bs.append(r.dot(b)) + ids.append(shot_id) + + if len(os_list) >= 2: + Xs = [] + + for i in range(len(os_list)): + o = os_list[i] + b = bs[i] + X = np.zeros(3) + + # https://math.stackexchange.com/questions/100439/determine-where-a-vector-will-intersect-a-plane + # [t == (a1*n1 + a2*n2 + a3*n3 - n1*o1 - n2*o2 - n3*o3)/(d1*n1 + d2*n2 + d3*n3) + d = b.dot(plane_normal) + if d == 0: + continue + + t = (plane_center.dot(plane_normal) - plane_normal.dot(o)) / d + + X[0] = o[0] + b[0] * t + X[1] = o[1] + b[1] * t + X[2] = o[2] + b[2] * t + Xs.append(X) + + if len(Xs) < 2: + return + + avgX = np.average(Xs, axis=0) + maxDist = np.max([np.linalg.norm(avgX - X) for X in Xs]) + if maxDist < threshold: + point = reconstruction.create_point(track, avgX) + point.color = np.array(list(observations.values())[0].color, dtype=int) + +def retriangulate_planar( + tracks_manager: pymap.TracksManager, reconstruction: types.Reconstruction, threshold: float = 0.02 +) -> Dict[str, Any]: + """Retrianguate all points""" + chrono = Chronometer() + report = {} + report['num_points_before'] = len(reconstruction.points) + + reconstruction.points = {} + + all_shots_ids = set(tracks_manager.get_shot_ids()) + + tracks = set() + for image in reconstruction.shots.keys(): + if image in all_shots_ids: + tracks.update(tracks_manager.get_shot_observations(image).keys()) + + for track in tracks: + triangulate_planar(tracks_manager, reconstruction, track, threshold) + + reconstruction.add_correspondences_from_tracks_manager(tracks_manager) + + report['num_points_after'] = len(reconstruction.points) + chrono.lap('retriangulate') + report['wall_time'] = chrono.total_time() + return report + +def find_planar_homography( + common_tracks_index: CommonTracksIndex, + common_tracks_data: NDArray[Any], + pair: Tuple[str, str], + graph: nx.DiGraph, + error_threshold: float, + ransac_error_threshold: float, +) -> Tuple[Optional[NDArray[Any]], Optional[Set[str]]]: + # Based on TRASAC: https://people.eecs.berkeley.edu/~yima/psfile/Planar-CVPR12.pdf + # Note: This implementation uses cv2.estimateAffinePartial2D (4-DoF) instead of + # a true projective homography (8-DoF). This is an intentional optimization for ODM + # as this works better for nadir-facing drone flights. + + num_trials = 15 + max_iters = 1000 + trials = [] + + num_tracks, rec_start, rec_end = common_tracks_index[pair] + ct_frame = common_tracks_data[rec_start:rec_end].reshape((num_tracks, 8)) + pair_tracks, p1, p2, hp1, hp2 = ct_frame[:,:1].T[0], ct_frame[:,1:3], ct_frame[:,4:6], ct_frame[:,1:4], ct_frame[:,4:7] + + adjacent_pairs = [] + for p in pair: + for n in graph[p]: + for ap in ((n,p), (p,n)): + if ap in common_tracks_index and ap != pair: + adjacent_pairs.append((ap, common_tracks_index[ap][0])) + + adjacent_pairs = sorted(adjacent_pairs, key=lambda e: e[1], reverse=True) + adjacent_pairs = adjacent_pairs[:2] + + i = 0 + while len(trials) < num_trials and i < max_iters: + i += 1 + + r = { + 'H': None, + 'inliers': set(), + 'pair_inliers': set(), + } + + # Sample 4 + sample_ids = random.sample(range(0, num_tracks), 4) + + # Compute homography + track_points1 = p1.take(sample_ids, axis=0) + track_points2 = p2.take(sample_ids, axis=0) + + H, _ = cv2.estimateAffinePartial2D(track_points1, track_points2) + if H is None: + # Degenerate sample (e.g. collinear points), try another one + continue + r['H'] = np.vstack([H, [0,0,1]]) + + # Classify each track + for j in range(num_tracks): + err = np.linalg.norm(hp2[j] - r['H'].dot(hp1[j])) + if err < error_threshold: + r['pair_inliers'].add(pair_tracks[j]) + + if len(r['pair_inliers']) <= 4: + continue + + r['inliers'] = set(r['pair_inliers']) + + for adj_pair, _ in adjacent_pairs: + adj_num_tracks, adj_rec_start, adj_rec_end = common_tracks_index[adj_pair] + adj_ct_frame = common_tracks_data[adj_rec_start:adj_rec_end].reshape((adj_num_tracks, 8)) + adj_pair_tracks, new_p1, new_p2, new_hp1, new_hp2 = adj_ct_frame[:,:1].T[0], adj_ct_frame[:,1:3], adj_ct_frame[:,4:6], adj_ct_frame[:,1:4], adj_ct_frame[:,4:7] + + inliers_common_tracks = set(adj_pair_tracks).intersection(r['pair_inliers']) + if len(inliers_common_tracks) <= 4: + continue + + inliers_p1 = [] + inliers_p2 = [] + for x in range(len(adj_pair_tracks)): + if adj_pair_tracks[x] in inliers_common_tracks: + inliers_p1.append(new_p1[x]) + inliers_p2.append(new_p2[x]) + inliers_p1 = np.reshape(inliers_p1, newshape=(len(inliers_p1), 2)) + inliers_p2 = np.reshape(inliers_p2, newshape=(len(inliers_p2), 2)) + + H, _ = cv2.estimateAffinePartial2D(inliers_p1, inliers_p2) + + if H is None: + continue + + H = np.vstack([H, [0,0,1]]) + + for j in range(adj_num_tracks): + err = np.linalg.norm(new_hp2[j] - H.dot(new_hp1[j])) + if err < error_threshold: + r['inliers'].add(adj_pair_tracks[j]) + + trials.append(r) + + max_inliers = -1 + best_t = None + for trial in trials: + num_inliers = len(trial['inliers']) + if num_inliers > max_inliers: + best_t = trial + max_inliers = num_inliers + + if best_t is None: + return None, None + + return best_t['H'], best_t['pair_inliers'] + +def Rt_from_H( + H: NDArray[Any], + K: NDArray[Any], + K1: NDArray[Any] +) -> Tuple[NDArray[Any], NDArray[Any]]: + """ + NOTE: This is not a generalized homography decomposition. But assumes + a nadir drone flight (constant altitude, zero pitch/roll). It extracts R directly + from the 2D affine matrix and forces the camera altitude to never change (t[2] = 0). + """ + H = H.copy() + + _, s, _ = np.linalg.svd(K1.dot(H).dot(K)) + H /= s[1] + + R = H.copy() + R[0][2] = R[1][2] = 0 + R[2][2] = 1.0 + + t = H[:,2] + t[0] /= K[0][0] + t[1] /= K[1][1] + t[2] = 0 + + return R, t + + +def _compute_planar_homography( + args: Tuple[Tuple[str, str], CommonTracksIndex, NDArray[Any], nx.DiGraph] +) -> Tuple[Tuple[str, str], Optional[NDArray[Any]], Optional[NDArray[Any]]]: + error_threshold = 0.002 + ransac_error_threshold = 0.004 + + pair, common_tracks_index, common_tracks_data, graph, = args + num_tracks, _, _ = common_tracks_index[pair] + + H, plane_inliers = find_planar_homography(common_tracks_index, common_tracks_data, pair, graph, error_threshold, ransac_error_threshold) + if H is None: + logger.warning("Could not compute homography for %s" % str(pair)) + return (pair, None, None) + + num_outliers = num_tracks - len(plane_inliers) + logger.info("%s <=> %s inliers: %s outliers: %s" % (pair[0], pair[1], len(plane_inliers), num_outliers)) + + return (pair, H, np.linalg.inv(H)) + +def planar_reconstruction( + data: DataSetBase, tracks_manager: pymap.TracksManager +) -> Tuple[Dict[str, Any], List[types.Reconstruction]]: + """Run the entire planar reconstruction pipeline.""" + + logger.info("Starting planar reconstruction") + report = {} + chrono = Chronometer() + images = tracks_manager.get_shot_ids() + if len(images) == 0: + logger.warning("No images") + exit(1) + + min_inliers = data.config["five_point_algo_min_inliers"] + + data.init_reference(images) + bundle_gcp = data.config["bundle_use_gcp"] + if data.config["align_method"] == "auto": + logger.info("Switching align_method to orientation_prior") + data.config["align_method"] = "orientation_prior" + #data.config["bundle_max_iterations"] = 10 + if bundle_gcp: + data.config["bundle_use_gcp"] = False + processes = data.config["processes"] + multithread = processes > 1 + + gcp = data.load_ground_control_points() + rig_assignments = rig.rig_assignments_per_image(data.load_rig_assignments()) + common_tracks_data, common_tracks_index = tracking.np_all_common_tracks_with_features(tracks_manager, min_common=min_inliers) + + camera_priors = data.load_camera_models() + rig_camera_priors = data.load_rig_cameras() + + camera_id = data.load_exif(images[0])["camera"] + camera = camera_priors[camera_id] + + K = camera.get_K() + K1 = np.linalg.inv(K) + + graph = nx.DiGraph() + + for pair in common_tracks_index: + graph.add_node(pair[0]) + graph.add_node(pair[1]) + num_tracks, _, _ = common_tracks_index[pair] + graph.add_edge(*pair, weight=num_tracks) + graph.add_edge(*(tuple(reversed(pair))), weight=num_tracks) + + # Find initial image by looking at the most centralized nodes + centrality = sorted(nx.degree_centrality(graph).items(), key=lambda c: (c[1], c[0]), reverse=True) + + # Compute homography graph + central_image = centrality[0][0] + + nodes = deque([central_image]) + parallel_args = [] + + while len(nodes) > 0: + node = nodes.popleft() + edges = sorted(graph[node].items(), key=lambda e: (e[1]['weight'], e[0]), reverse=True) + + for n, attrs in edges: + # Has this node's homography been computed ? + # (is there an edge with an homography connected to this node?) + connected = False + for e in graph[n]: + if 'H' in graph[n][e]: + connected = True + break + + if connected: + continue + + # Compute H + pair = (node, n) + if not pair in common_tracks_index: + pair = tuple(reversed(pair)) + + if not pair in common_tracks_index: + # Should never happen? + logger.warning("%s not in common tracks" % str(pair)) + continue + + # Mark edges with candidate pairs to be computed + # in parallel + if multithread: + graph[pair[0]][pair[1]]['H'] = True + graph[pair[1]][pair[0]]['H'] = True + parallel_args.append((pair, common_tracks_index, common_tracks_data, graph)) + else: + _, H, H1 = _compute_planar_homography((pair, common_tracks_index, common_tracks_data, graph)) + if H is not None: + graph[pair[0]][pair[1]]['H'] = H + graph[pair[1]][pair[0]]['H'] = H1 + + nodes.append(n) + + if multithread: + for pair, H, H1 in context.parallel_map(_compute_planar_homography, parallel_args, processes, backend="loky"): + if H is not None: + graph[pair[0]][pair[1]]['H'] = H + graph[pair[1]][pair[0]]['H'] = H1 + else: + del graph[pair[0]][pair[1]]['H'] + del graph[pair[1]][pair[0]]['H'] + + # Remove edges that have no homographies + edges = list(nx.edges(graph)) + for u,v in edges: + if 'H' not in graph[u][v]: + graph.remove_edge(u, v) + + # Remove isolated nodes + nodes = list(nx.nodes(graph)) + for n in nodes: + if len(graph[n]) == 0: + graph.remove_node(n) + logger.info("Image %s could not be added (not enough overlap with other images)" % n) + + #nx.write_gpickle(graph, "/datasets/brighton2/graph.gpickle") + # graph = nx.read_gpickle("/datasets/brighton2/graph.gpickle") + # print("loaded existing graph") + # centrality = sorted(nx.degree_centrality(graph).items(), key=lambda c: (c[1], c[0]), reverse=True) + # central_image = centrality[0][0] + + # The most central image may have been pruned above, recompute if needed + if central_image not in graph: + centrality = sorted(nx.degree_centrality(graph).items(), key=lambda c: (c[1], c[0]), reverse=True) + if len(centrality) == 0: + logger.warning("No images could be added (not enough overlap between images)") + report["decision"] = "Failed: no overlapping images" + return report, [] + central_image = centrality[0][0] + + rec = types.Reconstruction() + rec.reference = data.load_reference() + rec.cameras = camera_priors + rec.rig_cameras = rig_camera_priors + + # Compute homography chain for each node in the graph + # starting from the center image and add them to the reconstruction + logger.info("Creating reconstruction poses from homography graph") + for im, graph_path in nx.single_source_shortest_path(graph, central_image).items(): + H = np.identity(3) + for n1,n2 in zip(graph_path, graph_path[1:]): + H = graph[n1][n2]['H'].dot(H) + + R, t = Rt_from_H(H, K, K1) + + shot = rec.create_shot( + im, + camera_id, + pygeometry.Pose(R, t) + ) + + shot.metadata = helpers.get_image_metadata(data, im) + shot.scale = 1.0 + + # Triangulate and bundle + retriangulate_planar(tracks_manager, rec, 0.01) + logger.info("Triangulated %s points on plane" % len(rec.points)) + logger.info("Bundle shot poses") + bundle_shot_poses(rec, set(rec.shots.keys()), camera_priors, rig_camera_priors, data.config) + + retriangulate(tracks_manager, rec, data.config) + logger.info("Reconstructed %s points" % len(rec.points)) + + logger.info("Bundle adjustment") + + align_reconstruction(rec, gcp, data.config) + bundle(rec, camera_priors, rig_camera_priors, gcp, 0, data.config) + remove_outliers(rec, data.config) + retriangulate(tracks_manager, rec, data.config) + logger.info("Reconstructed %s points" % len(rec.points)) + + if bundle_gcp: + data.config["bundle_use_gcp"] = True + + align_reconstruction(rec, gcp, data.config) + bundle(rec, camera_priors, rig_camera_priors, gcp, 0, data.config) + remove_outliers(rec, data.config) + + align_reconstruction(rec, gcp, data.config) + paint_reconstruction(data, tracks_manager, rec) + + report["decision"] = "Success" + report["memory_usage"] = context.current_memory_usage() + + return report, [rec] diff --git a/opensfm/report.py b/opensfm/report.py index b6ad08158..9f492cf1e 100644 --- a/opensfm/report.py +++ b/opensfm/report.py @@ -118,7 +118,7 @@ def _quality_color_lower_is_better( class Report: - def __init__(self, data: DataSet, title: Optional[str] = None, accent_color: Optional[Tuple[int, int, int]] = None) -> None: + def __init__(self, data: DataSet, stats: Optional[Dict[str, Any]] = None, title: Optional[str] = None, accent_color: Optional[Tuple[int, int, int]] = None) -> None: self.output_path: str = os.path.join(data.data_path, "stats") self.dataset_name: str = os.path.basename(data.data_path) self.io_handler: io.IoFilesystemBase = data.io_handler @@ -136,7 +136,7 @@ def __init__(self, data: DataSet, title: Optional[str] = None, accent_color: Opt self.pdf.set_auto_page_break(auto=True, margin=MARGIN) self.pdf.add_page() - self.stats: Dict[str, Any] = data.load_stats() + self.stats: Dict[str, Any] = stats if stats is not None else data.load_stats() def save_report(self, filename: str) -> None: bytestring = self.pdf.output() @@ -406,25 +406,22 @@ def make_title(self) -> None: self.pdf.set_xy(title_x, MARGIN) self.pdf.set_font("Helvetica", "B", FONT_TITLE) self.pdf.set_text_color(*self.accent_color) - self.pdf.cell(0, 10, "OpenSfM", align="L") + self.pdf.cell(0, 10, "ODM", align="L") self.pdf.set_xy(title_x, MARGIN + 9) self.pdf.set_font("Helvetica", "", FONT_H1) self.pdf.set_text_color(*COLOR_PANEL) self.pdf.cell(0, 8, self.locale.t("quality_report"), align="L") - # Version number (right aligned) + # Version number (right aligned), read from ODM's VERSION file + version_file = os.path.join(os.path.dirname( + __file__), "..", "..", "..", "..", "..", "VERSION") + version = "" try: - out, _ = subprocess.Popen( - ["git", "describe", "--tags"], - stdout=subprocess.PIPE, - stderr=subprocess.PIPE, - ).communicate() - version = out.strip().decode() - except BaseException as e: - logger.warning( - f"Exception thrown while extracting 'git' version, {e}") - version = "" + with open(version_file, 'r') as f: + version = f.read().strip() + except Exception as e: + logger.warning("Invalid version file" + version_file + ": " + str(e)) version = "unknown" if version == "" else version @@ -439,8 +436,13 @@ def make_title(self) -> None: def make_dataset_summary(self) -> None: self._make_section(self.locale.t("dataset_summary")) + # NOTE: ODM's odm_report stage provides a humanized total time; + # fall back to OpenSfM's own step timings when running outside ODM. + odm_time_human = self.stats.get( + "odm_processing_statistics", {}).get("total_time_human") + rows = [ - [self.locale.t("dataset"), self.dataset_name], + #[self.locale.t("dataset"), self.dataset_name], [self.locale.t("date"), self.stats["processing_statistics"]["date"]], [ @@ -450,9 +452,15 @@ def make_dataset_summary(self) -> None: ], [ self.locale.t("processing_time"), - self.locale.format_time( + odm_time_human + if odm_time_human is not None + else self.locale.format_time( self.stats["processing_statistics"]["steps_times"]["Total Time"]), ], + [self.locale.t("capture_start"), + self.stats["processing_statistics"]["start_date"]], + [self.locale.t("capture_end"), + self.stats["processing_statistics"]["end_date"]], ] self._make_table(None, rows, True) self.pdf.set_xy(MARGIN, self.pdf.get_y() + SECTION_GAP) @@ -486,6 +494,9 @@ def make_processing_summary(self) -> None: if self._has_meaningful_gcp(): geo_string.append(self.locale.t("gcp")) + if "align" in self.stats: + geo_string = [self.locale.t("alignment")] + ratio_shots = rec_shots / init_shots * 100 if init_shots > 0 else -1 ratio_points = rec_points / init_points * 100 if init_points > 0 else -1 avg_track = self.stats["reconstruction_statistics"]["average_track_length"] @@ -507,7 +518,7 @@ def make_processing_summary(self) -> None: shots_thresholds, col_sizes, ) self._make_graded_row( - self.locale.t("reconstructed_points"), + self.locale.t("reconstructed_points_sparse"), ratio_points, f"{rec_points} / {init_points} ({ratio_points:.1f}%)", points_thresholds, col_sizes, @@ -564,14 +575,31 @@ def make_processing_summary(self) -> None: self.locale.t("ground_sampling_distance"), self.locale.format_gsd(gsd) if gsd > 0 else "N/A", ], - [ - self.locale.t("reconstructed_components"), - f"{self.stats['reconstruction_statistics']['components']}", - ], + # [ + # self.locale.t("reconstructed_components"), + # f"{self.stats['reconstruction_statistics']['components']}", + # ], [self.locale.t("geographic_referencing"), " + ".join(geo_string) if geo_string else self.locale.t("none")], ] + # Dense points (if available) + if self.stats.get("point_cloud_statistics"): + if self.stats["point_cloud_statistics"].get("dense"): + rows.insert(0, [ + self.locale.t("reconstructed_points_dense"), + f"{self.stats['point_cloud_statistics']['stats']['statistic'][0]['count']:,} points", + ]) + + # Alignment RMSE (coarse + fine) instead of GPS/GCP errors + if "align" in self.stats: + rows.append([ + self.locale.t("alignment_errors"), + self.locale.format_distance_label( + self.stats["align"]["coarse"]["rmse_3d"] + + self.stats["align"]["fine"]["rmse_3d"]), + ]) + # CRS of the georeferenced products (centralized decision: GCP CRS if # projected, else UTM). Fall back to the GCP CRS for older stats files. output_crs = self.stats.get("output_coordinate_system") or self.stats.get( @@ -585,49 +613,51 @@ def make_processing_summary(self) -> None: self._make_table(None, rows, True) self.pdf.set_xy(MARGIN, self.pdf.get_y() + TABLE_GAP) - # GPS error graded row (based on multiples of input std dev) - if self.stats["reconstruction_statistics"]["has_gps"] and "average_error" in self.stats.get("gps_errors", {}): - gps_avg_error = self.stats["gps_errors"]["average_error"] - avg_gps_std = self.stats["gps_errors"].get("average_gps_std") - if avg_gps_std: - # Use the average of XYZ std as reference - ref_std = (avg_gps_std["x"] + - avg_gps_std["y"] + avg_gps_std["z"]) / 3.0 - if ref_std > 1e-9: - # Thresholds: good <= 2*std, avg <= 3*std, bad > 4*std (lower is better) - gps_thresholds = (2.0 * ref_std, 3.0 * - ref_std, 4.0 * ref_std) - self._make_graded_row_lower( - self.locale.t("gps_error"), - gps_avg_error, - self.locale.format_distance_label(gps_avg_error), - gps_thresholds, col_sizes, - ) - - # GCP error graded row (based on GSD, same Z thresholds: good <= 3*GSD, avg <= 4*GSD, bad > 5*GSD) - if self._has_meaningful_gcp() and gsd > 0: - gcp_only = self.stats["gcp_errors"].get("gcp_only", {}) - gcp_avg_error = gcp_only.get( - "average_error", self.stats["gcp_errors"]["average_error"]) - gcp_thresholds = (3.0 * gsd, 4.0 * gsd, 5.0 * gsd) - self._make_graded_row_lower( - self.locale.t("gcp_error"), - gcp_avg_error, - self.locale.format_distance_label(gcp_avg_error, precision=3), - gcp_thresholds, col_sizes, - ) + # GPS/GCP/CP error graded rows (skipped in alignment mode) + if "align" not in self.stats: + # GPS error graded row (based on multiples of input std dev) + if self.stats["reconstruction_statistics"]["has_gps"] and "average_error" in self.stats.get("gps_errors", {}): + gps_avg_error = self.stats["gps_errors"]["average_error"] + avg_gps_std = self.stats["gps_errors"].get("average_gps_std") + if avg_gps_std: + # Use the average of XYZ std as reference + ref_std = (avg_gps_std["x"] + + avg_gps_std["y"] + avg_gps_std["z"]) / 3.0 + if ref_std > 1e-9: + # Thresholds: good <= 2*std, avg <= 3*std, bad > 4*std (lower is better) + gps_thresholds = (2.0 * ref_std, 3.0 * + ref_std, 4.0 * ref_std) + self._make_graded_row_lower( + self.locale.t("gps_error"), + gps_avg_error, + self.locale.format_distance_label(gps_avg_error), + gps_thresholds, col_sizes, + ) + + # GCP error graded row (based on GSD, same Z thresholds: good <= 3*GSD, avg <= 4*GSD, bad > 5*GSD) + if self._has_meaningful_gcp() and gsd > 0: + gcp_only = self.stats["gcp_errors"].get("gcp_only", {}) + gcp_avg_error = gcp_only.get( + "average_error", self.stats["gcp_errors"]["average_error"]) + gcp_thresholds = (3.0 * gsd, 4.0 * gsd, 5.0 * gsd) + self._make_graded_row_lower( + self.locale.t("gcp_error"), + gcp_avg_error, + self.locale.format_distance_label(gcp_avg_error, precision=3), + gcp_thresholds, col_sizes, + ) - # CP error graded row - if self._has_meaningful_cp() and gsd > 0: - cp_only = self.stats["gcp_errors"]["cp_only"] - cp_avg_error = cp_only["average_error"] - cp_thresholds = (3.0 * gsd, 4.0 * gsd, 5.0 * gsd) - self._make_graded_row_lower( - self.locale.t("gcp_error"), - cp_avg_error, - self.locale.format_distance_label(cp_avg_error, precision=3), - cp_thresholds, col_sizes, - ) + # CP error graded row + if self._has_meaningful_cp() and gsd > 0: + cp_only = self.stats["gcp_errors"]["cp_only"] + cp_avg_error = cp_only["average_error"] + cp_thresholds = (3.0 * gsd, 4.0 * gsd, 5.0 * gsd) + self._make_graded_row_lower( + self.locale.t("gcp_error"), + cp_avg_error, + self.locale.format_distance_label(cp_avg_error, precision=3), + cp_thresholds, col_sizes, + ) self.pdf.set_xy(MARGIN, self.pdf.get_y() + TABLE_GAP) @@ -857,26 +887,26 @@ def make_gps_details(self) -> None: self._make_table(["", self.locale.t("absolute"), self.locale.t("relative")], rows, False) self.pdf.set_xy(MARGIN, self.pdf.get_y() + TABLE_GAP) - rows = [] - columns_names = [ - self.locale.t("gps_bias"), - self.locale.t("scale"), - self.locale.t("translation"), - self.locale.t("rotation"), - ] - for camera, params in self.stats["camera_errors"].items(): - bias = params["bias"] - s, t, R = bias["scale"], bias["translation"], bias["rotation"] - rows.append( - [ - camera, - f"{s:.2f}", - f"{t[0]:.2f} {t[1]:.2f} {t[2]:.2f}", - f"{R[0]:.2f} {R[1]:.2f} {R[2]:.2f}", - ] - ) - if rows: - self._make_table(columns_names, rows) + # rows = [] + # columns_names = [ + # self.locale.t("gps_bias"), + # self.locale.t("scale"), + # self.locale.t("translation"), + # self.locale.t("rotation"), + # ] + # for camera, params in self.stats["camera_errors"].items(): + # bias = params["bias"] + # s, t, R = bias["scale"], bias["translation"], bias["rotation"] + # rows.append( + # [ + # camera, + # f"{s:.2f}", + # f"{t[0]:.2f} {t[1]:.2f} {t[2]:.2f}", + # f"{R[0]:.2f} {R[1]:.2f} {R[2]:.2f}", + # ] + # ) + # if rows: + # self._make_table(columns_names, rows) self.pdf.set_xy(MARGIN, self.pdf.get_y() + TABLE_GAP) @@ -1004,6 +1034,121 @@ def make_gcp_details(self) -> None: self.pdf.set_xy(MARGIN, self.pdf.get_y() + TABLE_GAP) + def _read_gcp_stats_file(self, filename): + file_path = os.path.join(self.output_path, "ground_control_points.json") + + with self.io_handler.open_rt(file_path) as fin: + return io.json_load(fin) + + def make_gcp_error_details(self) -> None: + gcp_stats_path = os.path.join(self.output_path, "ground_control_points.json") + if not self.io_handler.isfile(gcp_stats_path): + return + + self._make_section("Ground Control Point Error") + + gcp_stats = self._read_gcp_stats_file(gcp_stats_path) + + rows = [] + column_names = ["ID", "Error X (m)", "Error Y (m)", "Error Z (m)"] + + for gcp in gcp_stats: + gcp_id = ''.join(c if ord(c) < 256 else '?' for c in gcp["id"]) # latin-1 only due to fpdf2 requirements + row = [gcp_id] + row.append(f"{gcp['error'][0]:.3f}") + row.append(f"{gcp['error'][1]:.3f}") + row.append(f"{gcp['error'][2]:.3f}") + + rows.append(row) + + self._make_table(column_names, rows) + self.pdf.set_xy(MARGIN, self.pdf.get_y() + MARGIN / 2) + + def make_align_details(self) -> None: + self._make_section("Alignment Errors Details") + + # Alignment + rows = [] + columns_names = ["", "DSM (Coarse)", "ICP (Fine)", "RMS Error"] + for comp in ["x", "y", "z", "3d"]: + row = [comp.upper() + " Error (meters)"] + row.append(f"{self.stats['align']['coarse']['rmse_' + comp]:.3f}") + row.append(f"{self.stats['align']['fine']['rmse_' + comp]:.3f}") + row.append(f"{(self.stats['align']['coarse']['rmse_' + comp] + self.stats['align']['fine']['rmse_' + comp]):.3f}") + rows.append(row) + + self._make_table(columns_names, rows) + self.pdf.set_xy(MARGIN, self.pdf.get_y() + MARGIN / 2) + + dsm_feature_matches = os.path.join(self.output_path, "codem", "dsm_feature_matches.png") + if self.io_handler.isfile(dsm_feature_matches): + self._make_centered_image(dsm_feature_matches, 80) + + def make_survey_data(self) -> None: + self._make_section("Survey Data") + + self._make_centered_image( + os.path.join(self.output_path, "overlap.png"), 90 + ) + self._make_centered_image( + os.path.join(self.output_path, "overlap_diagram_legend.png"), 3 + ) + + self.pdf.set_xy(MARGIN, self.pdf.get_y() + MARGIN / 2) + + def _add_image_label(self, text: str) -> None: + self.pdf.set_font_size(FONT_SMALL) + self.pdf.text(self.pdf.get_x() + CONTENT_WIDTH / 2 - self.pdf.get_string_width(text) / 2, self.pdf.get_y() - 5, text) + + def make_preview(self) -> bool: + ortho = os.path.join(self.output_path, "ortho.png") + dsm = os.path.join(self.output_path, "dsm.png") + dtm = os.path.join(self.output_path, "dtm.png") + count = 0 + + if self.io_handler.isfile(ortho) or self.io_handler.isfile(dsm): + self._make_section("Previews") + + if self.io_handler.isfile(ortho): + self._make_centered_image(ortho, 110) + self._add_image_label("Orthophoto") + count += 1 + + if self.io_handler.isfile(dsm) and self.stats.get('dsm_statistics'): + self._make_centered_image(dsm, 110) + self._add_image_label("Digital Surface Model") + + self._make_centered_image( + os.path.join(self.output_path, "dsm_gradient.png"), 4 + ) + self.pdf.set_font_size(FONT_SMALL) + min_text = "{:,.2f}m".format(self.stats['dsm_statistics']['min']) + max_text = "{:,.2f}m".format(self.stats['dsm_statistics']['max']) + self.pdf.text(self.pdf.get_x() + 40, self.pdf.get_y() - 5, min_text) + self.pdf.text(self.pdf.get_x() + 40 + 110.5 - self.pdf.get_string_width(max_text), self.pdf.get_y() - 5, max_text) + count += 1 + + if self.io_handler.isfile(dtm) and self.stats.get('dtm_statistics'): + if count >= 2: + self.add_page_break() + + self._make_centered_image(dtm, 110) + self._add_image_label("Digital Terrain Model") + + self._make_centered_image( + os.path.join(self.output_path, "dsm_gradient.png"), 4 + ) + self.pdf.set_font_size(FONT_SMALL) + min_text = "{:,.2f}m".format(self.stats['dtm_statistics']['min']) + max_text = "{:,.2f}m".format(self.stats['dtm_statistics']['max']) + self.pdf.text(self.pdf.get_x() + 40, self.pdf.get_y() - 5, min_text) + self.pdf.text(self.pdf.get_x() + 40 + 110.5 - self.pdf.get_string_width(max_text), self.pdf.get_y() - 5, max_text) + + self.pdf.set_xy(MARGIN, self.pdf.get_y() + MARGIN) + + return True + return False + def make_orientation_details(self) -> None: if "opk_errors" not in self.stats: return @@ -1287,7 +1432,7 @@ def generate_report(self) -> None: self.make_tracks_details() self.make_camera_models_details() - self.make_rig_cameras_details() + #self.make_rig_cameras_details() self.add_page_break() report_images = self.io_handler.ls(self.output_path) @@ -1299,7 +1444,16 @@ def generate_report(self) -> None: self.add_page_break() self.make_overlap_summary() - self.make_gps_details() - self.make_gcp_details() + if "align" not in self.stats: + self.make_gps_details() + self.make_gcp_details() + self.make_gcp_error_details() + else: + self.make_align_details() self.make_orientation_details() self.make_processing_time_details() + + if self.io_handler.isfile(os.path.join(self.output_path, "overlap.png")): + self.make_survey_data() + + self.make_preview() diff --git a/opensfm/sensors.py b/opensfm/sensors.py index b91f115df..6ce10d281 100644 --- a/opensfm/sensors.py +++ b/opensfm/sensors.py @@ -3,11 +3,25 @@ from typing import Any, Dict, List import yaml -from opensfm import context, io +import os +import sqlite3 +from opensfm import context, io @lru_cache(1) def sensor_data() -> Dict[str, Any]: + if os.path.isfile(context.SENSOR_DATA_DB): + try: + conn = sqlite3.connect(context.SENSOR_DATA_DB) + cur = conn.cursor() + cur.execute("SELECT id,focal FROM sensors") + rows = cur.fetchall() + conn.close() + return {r[0]: r[1] for r in rows} + except Exception as e: + print("Cannot query %s: %s" % (context.SENSOR_DATA_DB, str(e))) + # Fallback to sensor_data.json + with io.open_rt(context.SENSOR_DATA) as f: data = io.json_load(f) diff --git a/opensfm/src/lib/CMakeLists.txt b/opensfm/src/lib/CMakeLists.txt index d643a08e3..b537beebd 100644 --- a/opensfm/src/lib/CMakeLists.txt +++ b/opensfm/src/lib/CMakeLists.txt @@ -43,9 +43,16 @@ endif() ####### Find Dependencies ####### # Note: Common dependencies are in ../cmake/OpenSfMCommon.cmake +if (ADDITIONAL_INCLUDE_DIRS) + include_directories(${ADDITIONAL_INCLUDE_DIRS}) +endif() +if (YET_ADDITIONAL_INCLUDE_DIRS) + include_directories(${YET_ADDITIONAL_INCLUDE_DIRS}) +endif() + find_package(LAPACK) find_package(SuiteSparse) -find_package(Ceres) +find_package(Ceres REQUIRED) find_package(absl REQUIRED) find_library(TCMALLOC_LIBRARY tcmalloc) diff --git a/opensfm/src/lib/bundle/src/bundle_adjuster.cc b/opensfm/src/lib/bundle/src/bundle_adjuster.cc index 0795c5ca9..29d0fcf0c 100644 --- a/opensfm/src/lib/bundle/src/bundle_adjuster.cc +++ b/opensfm/src/lib/bundle/src/bundle_adjuster.cc @@ -11,6 +11,7 @@ #include #include +#include #include #include #include @@ -674,8 +675,14 @@ void BundleAdjuster::Run() { for (auto& [_, cam] : cameras_) { auto& data = cam.GetValueData(); solver.AddParameterBlock(data.data(), data.size()); - ApplyParameterLocking(solver, data.data(), data.size(), - cam.GetParametersToOptimize()); + + // Keep aspect ratio (index 6) constant (BROWN only) + auto params = cam.GetParametersToOptimize(); + if (!params.empty() && + cam.GetValue().GetProjectionType() == geometry::ProjectionType::BROWN) { + params.erase(std::remove(params.begin(), params.end(), 6), params.end()); + } + ApplyParameterLocking(solver, data.data(), data.size(), params); // Add a barrier for constraining transition of dual to stay in [0, 1] const auto camera = cam.GetValue(); diff --git a/opensfm/src/lib/features/src/akaze_bind.cc b/opensfm/src/lib/features/src/akaze_bind.cc index 4131d23d6..bd86bd378 100644 --- a/opensfm/src/lib/features/src/akaze_bind.cc +++ b/opensfm/src/lib/features/src/akaze_bind.cc @@ -29,7 +29,7 @@ py::tuple akaze(foundation::pyarray_uint8 image, AKAZEOptions options) { cv::Mat desc; evolution.Compute_Descriptors(kpts, desc); - evolution.Show_Computation_Times(); + // evolution.Show_Computation_Times(); // Convert to numpy. cv::Mat keys(kpts.size(), 4, CV_32F); diff --git a/opensfm/src/lib/map/src/tracks_manager.cc b/opensfm/src/lib/map/src/tracks_manager.cc index af320e595..1bc70fef4 100644 --- a/opensfm/src/lib/map/src/tracks_manager.cc +++ b/opensfm/src/lib/map/src/tracks_manager.cc @@ -130,24 +130,79 @@ static map::TracksManager ParseTracksBuffer(const char* data, size_t size, return manager; } +struct TrackLengths{ + uint16_t imageLen; + uint16_t trackIdLen; +}; + +struct TrackRecord{ + int featureID; + float x; + float y; + float scale; + uint8_t r; + uint8_t g; + uint8_t b; + // int segm; + // int inst; +}; + +map::TracksManager InstanciateFromFilenameBinaryV2(const std::string &filename) { + std::ifstream fs(filename, std::ios_base::binary); + + map::TracksManager manager; + + // Skip version line + char buffer[131072]; + std::string version; + std::getline(fs, version); + + TrackLengths tl; + TrackRecord tr; + + while(!fs.eof()){ + fs.read(reinterpret_cast(&tl), sizeof(TrackLengths)); + fs.read(buffer, tl.imageLen + tl.trackIdLen); + std::string image(buffer, tl.imageLen); + std::string trackID(buffer + tl.imageLen, tl.trackIdLen); + + fs.read(reinterpret_cast(&tr), sizeof(TrackRecord)); + auto observation = InstanciateObservation(tr.x, tr.y, tr.scale, tr.featureID, tr.r, tr.g, tr.b, -1, -1); + manager.AddObservation(image, trackID, observation); + } + + return manager; +} + template void WriteToStreamCurrentVersion(S& ostream, const map::TracksManager& manager) { ostream << manager.TRACKS_HEADER << "_v" << manager.TRACKS_VERSION << std::endl; + + TrackLengths tl; + TrackRecord tr; + const auto shotsIDs = manager.GetShotIds(); for (const auto& shotID : shotsIDs) { const auto observations = manager.GetShotObservations(shotID); for (const auto& observation : observations) { - ostream << shotID << "\t" << observation.first << "\t" - << observation.second.feature_id << "\t" - << observation.second.point(0) << "\t" - << observation.second.point(1) << "\t" << observation.second.scale - << "\t" << observation.second.color(0) << "\t" - << observation.second.color(1) << "\t" - << observation.second.color(2) << "\t" - << observation.second.segmentation_id << "\t" - << observation.second.instance_id << std::endl; + + tl.imageLen = shotID.length(); + tl.trackIdLen = observation.first.length(); + tr.featureID = observation.second.feature_id; + tr.x = observation.second.point(0); + tr.y = observation.second.point(1); + tr.scale = observation.second.scale; + tr.r = static_cast(observation.second.color(0)); + tr.g = static_cast(observation.second.color(1)); + tr.b = static_cast(observation.second.color(2)); + // tr.segm = observation.second.segmentation_id; + // tr.inst = observation.second.instance_id; + + ostream.write(reinterpret_cast(&tl), sizeof(tl)); + ostream << shotID << observation.first; + ostream.write(reinterpret_cast(&tr), sizeof(tr)); } } } @@ -687,11 +742,14 @@ TracksManager TracksManager::InstanciateFromFile(const std::string& filename) { const char* data = mmap.data(); const size_t size = mmap.length(); const auto [version, offset] = DetectVersion(data, size); + if (version == 102) { + return InstanciateFromFilenameBinaryV2(filename); + } return ParseTracksBuffer(data + offset, size - offset, version); } void TracksManager::WriteToFile(const std::string& filename) const { - std::ofstream ostream(filename); + std::ofstream ostream(filename, std::ios_base::binary); if (ostream.is_open()) { WriteToStreamCurrentVersion(ostream, *this); } else { @@ -700,17 +758,19 @@ void TracksManager::WriteToFile(const std::string& filename) const { } TracksManager TracksManager::InstanciateFromString(const std::string& str) { - const auto [version, offset] = DetectVersion(str.data(), str.size()); - return ParseTracksBuffer(str.data() + offset, str.size() - offset, version); +// const auto [version, offset] = DetectVersion(str.data(), str.size()); +// return ParseTracksBuffer(str.data() + offset, str.size() - offset, version); + throw std::runtime_error("Instantiation from string not supported: " + str); } std::string TracksManager::AsString() const { - std::stringstream sstream; - WriteToStreamCurrentVersion(sstream, *this); - return sstream.str(); +// std::stringstream sstream; +// WriteToStreamCurrentVersion(sstream, *this); +// return sstream.str(); + throw std::runtime_error("Serialization from string not supported"); } std::string TracksManager::TRACKS_HEADER = "OPENSFM_TRACKS_VERSION"; -int TracksManager::TRACKS_VERSION = 2; +int TracksManager::TRACKS_VERSION = 102; } // namespace map diff --git a/opensfm/src/lib/map/test/tracks_manager_test.cc b/opensfm/src/lib/map/test/tracks_manager_test.cc index db891b1b4..f6b46d641 100644 --- a/opensfm/src/lib/map/test/tracks_manager_test.cc +++ b/opensfm/src/lib/map/test/tracks_manager_test.cc @@ -200,19 +200,18 @@ TEST_F(TracksManagerTest, HasIOFileConsistency) { ::testing::WhenSorted(::testing::ElementsAre("1", "2", "3"))); EXPECT_THAT(manager_new.GetTrackIds(), ::testing::WhenSorted(::testing::ElementsAre("1"))); - EXPECT_EQ(track, manager_new.GetTrackObservations("1")); -} -TEST_F(TracksManagerTest, HasIOStringConsistency) { - const auto serialized = manager.AsString(); - const map::TracksManager manager_new = - map::TracksManager::InstanciateFromString(serialized); + std::unordered_map track_from_file = manager_new.GetTrackObservations("1"); - EXPECT_THAT(manager_new.GetShotIds(), - ::testing::WhenSorted(::testing::ElementsAre("1", "2", "3"))); - EXPECT_THAT(manager_new.GetTrackIds(), - ::testing::WhenSorted(::testing::ElementsAre("1"))); - EXPECT_EQ(track, manager_new.GetTrackObservations("1")); + // Test that point, scale, color, and feature ID are correctly recovered from the file. + // Segmentation and instance IDs are not saved in the ODM file format. For this reason, + // track_from_file will not exactly match track (it is missing these fields). + for(auto const& it : track) { + EXPECT_EQ(it.second.point, track_from_file[it.first].point); + EXPECT_EQ(it.second.scale, track_from_file[it.first].scale); + EXPECT_EQ(it.second.color, track_from_file[it.first].color); + EXPECT_EQ(it.second.feature_id, track_from_file[it.first].feature_id); + } } } // namespace diff --git a/opensfm/src/lib/third_party/openmvs/Interface.h b/opensfm/src/lib/third_party/openmvs/Interface.h index 7fe6b4e4b..9a2436096 100644 --- a/opensfm/src/lib/third_party/openmvs/Interface.h +++ b/opensfm/src/lib/third_party/openmvs/Interface.h @@ -416,6 +416,23 @@ struct Interface inline Pos3d GetTranslation() const { return R*(-C); } inline void SetTranslation(const Pos3d& T) { C = R.t()*(-T); } + // combine poses + inline Pose operator * (const Pose& P) const { + return Pose(R*P.R, P.R.t()*C+P.C); + } + inline Pose& operator *= (const Pose& P) { + R = R*P.R; C = P.R.t()*C+P.C; return *this; + } + + // project point: world to local coordinates + inline Pos3d operator * (const Pos3d& X) const { + return R * (X - C); + } + // back-project point: local to world coordinates + inline Pos3d operator / (const Pos3d& X) const { + return R.t() * X + C; + } + template void serialize(Archive& ar, const unsigned int /*version*/) { ar & R; @@ -631,6 +648,29 @@ struct Interface return platforms[image.platformID].GetPose(image.cameraID, image.poseID); } + // apply similarity transform + void Transform(const Mat33d& rotation, const Pos3d& translation, const double scale) { + for (Platform& platform : platforms) { + for (Platform::Pose& pose : platform.poses) { + pose.R = pose.R * rotation.t(); + pose.C = rotation * pose.C * scale + translation; + } + } + for (Vertex& vertex : vertices) { + vertex.X = rotation * Pos3d(vertex.X) * scale + translation; + } + for (Normal& normal : verticesNormal) { + normal.n = rotation * Pos3d(normal.n); + } + for (Line& line : lines) { + line.pt1 = rotation * Pos3d(line.pt1) * scale + translation; + line.pt2 = rotation * Pos3d(line.pt2) * scale + translation; + } + for (Normal& normal : linesNormal) { + normal.n = rotation * Pos3d(normal.n); + } + } + template void serialize(Archive& ar, const unsigned int version) { ar & platforms; diff --git a/opensfm/src/lib/third_party/vlfeat/CMakeLists.txt b/opensfm/src/lib/third_party/vlfeat/CMakeLists.txt index c981fed10..8bfaeb8e3 100644 --- a/opensfm/src/lib/third_party/vlfeat/CMakeLists.txt +++ b/opensfm/src/lib/third_party/vlfeat/CMakeLists.txt @@ -18,7 +18,8 @@ else() set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -Wno-tautological-constant-out-of-range-compare") endif() -if( ${CMAKE_SYSTEM_PROCESSOR} STREQUAL "aarch64" ) +# vlfeat SSE2 disabled on arm +if( (${CMAKE_SYSTEM_PROCESSOR} STREQUAL "aarch64") OR (${CMAKE_SYSTEM_PROCESSOR} STREQUAL "arm64") ) add_definitions( -DVL_DISABLE_SSE2 ) endif() diff --git a/opensfm/stats.py b/opensfm/stats.py index 1b06833cf..94c2a50fb 100644 --- a/opensfm/stats.py +++ b/opensfm/stats.py @@ -1,5 +1,6 @@ # pyre-strict import datetime +import json import logging import math import os @@ -23,7 +24,7 @@ import numpy as np from numpy.typing import NDArray -from opensfm import feature_loader, geo, geometry, report, io, multiview, pygeometry, pymap, types +from opensfm import feature_loader, features, geo, geometry, report, io, multiview, pygeometry, pymap, types from opensfm.dataset import DataSet, DataSetBase logger: logging.Logger = logging.getLogger(__name__) @@ -199,6 +200,7 @@ def gcp_errors( all_errors = [] gcp_details: List[Dict[str, Any]] = [] + gcp_stats = [] for gcp in gcps: if not gcp.lla: continue @@ -253,6 +255,37 @@ def gcp_errors( "sigma": sigma_xyz, }) + # Export GCP stats for ODM (ground_control_points.json) + observations = [] + for obs in gcp.observations: + if obs.shot_id not in rec.shots: + continue + shot = rec.shots[obs.shot_id] + + reprojected = shot.project(gcp_enu) + annotated = obs.projection + + r_pixel = features.denormalized_image_coordinates(np.array([[reprojected[0], reprojected[1]]]), shot.camera.width, shot.camera.height)[0] + r_pixel[0] /= shot.camera.width + r_pixel[1] /= shot.camera.height + + a_pixel = features.denormalized_image_coordinates(np.array([[annotated[0], annotated[1]]]), shot.camera.width, shot.camera.height)[0] + a_pixel[0] /= shot.camera.width + a_pixel[1] /= shot.camera.height + + observations.append({ + 'shot_id': obs.shot_id, + 'annotated': list(a_pixel), + 'reprojected': list(r_pixel) + }) + + gcp_stats.append({ + 'id': gcp.id, + 'coordinates': list(gcp_enu), + 'observations': observations, + 'error': list(error) + }) + # Separate GCP-only and CP-only errors gcp_only_errors = [ e for e, d in zip(all_errors, [dd for dd in gcp_details if dd["error"] is not None]) @@ -281,6 +314,10 @@ def gcp_errors( crs = data.load_gcp_coordinate_system() if crs: stats["coordinate_system"] = crs + + with data.io_handler.open_wt(os.path.join(data.data_path, "stats", "ground_control_points.json")) as f: + f.write(json.dumps(gcp_stats, indent=4)) + return stats @@ -747,6 +784,14 @@ def processing_statistics( except FileNotFoundError: stats["date"] = "unknown" + start_ct, end_ct = start_end_capture_time(reconstructions) + if start_ct is not None and end_ct is not None: + stats["start_date"] = datetime.datetime.utcfromtimestamp(start_ct).strftime("%d/%m/%Y at %H:%M:%S") + stats["end_date"] = datetime.datetime.utcfromtimestamp(end_ct).strftime("%d/%m/%Y at %H:%M:%S") + else: + stats["start_date"] = "unknown" + stats["end_date"] = "unknown" + default_max = 1e30 min_x, min_y, max_x, max_y = default_max, default_max, 0, 0 for rec in reconstructions: @@ -1258,9 +1303,11 @@ def save_topview( ) plt.yticks( [im_size_y, im_size_y / 2, 0], - [0, f"{int(size_y / 2):.0f}", f"{size_y:.0f} meters"], + [f"{size_y:.0f} meters", f"{int(size_y / 2):.0f}", 0], fontsize="small", ) + plt.gca().invert_yaxis() + with io_handler.open_wb(os.path.join(output_path, "topview.png")) as fwb: plt.savefig( fwb, @@ -2017,3 +2064,22 @@ def decimate_points( for point_id in random_ids: rec.remove_point(point_id) + + +def start_end_capture_time(reconstructions): + end_ct = float('-inf') + start_ct = float('inf') + + for reconstruction in reconstructions: + for shot in reconstruction.shots.values(): + if shot.metadata.capture_time.has_value: + v = shot.metadata.capture_time.value + if v > end_ct: + end_ct = v + if v < start_ct: + start_ct = v + + if end_ct != float('-inf') and start_ct != float('inf'): + return (start_ct, end_ct) + else: + return (None, None) diff --git a/opensfm/synthetic_data/synthetic_scene.py b/opensfm/synthetic_data/synthetic_scene.py index 2bb9f25ae..6518b75df 100644 --- a/opensfm/synthetic_data/synthetic_scene.py +++ b/opensfm/synthetic_data/synthetic_scene.py @@ -548,3 +548,4 @@ def compare( if len(rig_instance_assignment) > 0 else 1.0, } + diff --git a/opensfm/test/test_commands.py b/opensfm/test/test_commands.py index d38df4c57..53e66bdfe 100644 --- a/opensfm/test/test_commands.py +++ b/opensfm/test/test_commands.py @@ -6,7 +6,7 @@ import pytest -from opensfm import commands, dataset +from opensfm import commands, dataset, reconstruction from opensfm.test import data_generation, utils @@ -44,6 +44,7 @@ def test_run_all(tmpdir: Any) -> None: commands.export_colmap, commands.compute_statistics, commands.export_report, + commands.export_geocoords, ] output_rec_path = join(data.data_path, "rec_prior.json") @@ -53,6 +54,17 @@ def test_run_all(tmpdir: Any) -> None: join(data.data_path, "reconstruction.json"), "--output", output_rec_path, + ], + commands.export_geocoords: [ + "--proj", + "+proj=utm +zone=33 +datum=WGS84 +units=m +no_defs", + "--offset-x", + "100", + "--offset-y", + "200", + "--mode", + "projected", + "--reconstruction" ] } @@ -63,6 +75,7 @@ def test_run_all(tmpdir: Any) -> None: check_reconstruction(data) check_prior(data, output_rec_path) + check_export_geocoords(data) def check_reconstruction(data: dataset.DataSet) -> None: @@ -78,3 +91,47 @@ def check_prior(data: dataset.DataSet, output_rec_path: str) -> None: utils.assert_shots_equal(shot, prior_rec[0].shots[shot_id]) assert len(prior_rec[0].points) > 1000 + + +def check_export_geocoords(data: dataset.DataSet) -> None: + import os + geocoords_path = join(data.data_path, "reconstruction.geocoords.json") + assert os.path.exists(geocoords_path) + + # Load and check if points/shots exist in the exported file + geocoords_rec = data.load_reconstruction("reconstruction.geocoords.json") + assert len(geocoords_rec) > 0 + assert len(geocoords_rec[0].shots) == 3 + assert len(geocoords_rec[0].points) > 1000 + + +def test_command_reconstruct_config(tmpdir: Any, monkeypatch: Any) -> None: + """ + Mocks Reconstruct command algorithm setting from config.py::reconstruction_algorithm + """ + data = data_generation.create_berlin_test_folder(tmpdir) + # Set the config to use triangulation + data_generation.save_config({"reconstruction_algorithm": "triangulation"}, data.data_path) + + called_algorithm = None + + def mock_triangulation(dataset: Any, tracks: Any) -> Any: + nonlocal called_algorithm + called_algorithm = "triangulation" + return {}, [] + + def mock_incremental(dataset: Any, tracks: Any) -> Any: + nonlocal called_algorithm + called_algorithm = "incremental" + return {}, [] + + # Mock the actual algorithms so it doesn't take forever to run + monkeypatch.setattr(reconstruction, "triangulation_reconstruction", mock_triangulation) + monkeypatch.setattr(reconstruction, "incremental_reconstruction", mock_incremental) + monkeypatch.setattr(dataset.DataSet, "load_tracks_manager", lambda self, *args, **kwargs: None) + + # Run the command with no --algorithm flag, it should fall back to config + command = commands.reconstruct.Command() + run_command(command, [data.data_path]) + + assert called_algorithm == "triangulation" diff --git a/opensfm/test/test_extract_metadata.py b/opensfm/test/test_extract_metadata.py index 801c7aeec..6db2b6eeb 100644 --- a/opensfm/test/test_extract_metadata.py +++ b/opensfm/test/test_extract_metadata.py @@ -56,6 +56,14 @@ def save_camera_models(self, camera_models: Dict[str, Any]) -> None: def init_reference(self) -> None: self.init_reference_calls += 1 + def open_image_file(self, image: str) -> Any: + class DummyFile: + def __enter__(self) -> Any: return self + def __exit__(self, *args: Any) -> None: pass + return DummyFile() + + def image_size(self, image: str) -> Any: + return (1000, 1000) def test_run_dataset_force_reextracts_and_initializes_reference( monkeypatch: Any, @@ -157,3 +165,24 @@ def fake_run_dataset(dataset: Any, force: bool = False) -> None: command.run_impl(dataset, args) assert called == {"dataset": dataset, "force": True} + + +def test_extract_exif_overrides_projection_type(monkeypatch: Any) -> None: + """ + Tests flag camera_projection_type (Overrides if != AUTO[Default]) + """ + data = FakeData(["image.jpg"], {}, {}) + data.config["camera_projection_type"] = "BROWN" + data.config["use_exif_size"] = False + data.config["default_projection_type"] = "perspective" + data.config["unknown_camera_models_are_different"] = False + + monkeypatch.setattr( + extract_metadata_action.exif, + "extract_exif_from_file", + lambda *a, **k: {"model": "m", "camera": "c", "projection_type": "perspective"} + ) + monkeypatch.setattr(extract_metadata_action.exif, "camera_id", lambda d: "c") + + d = extract_metadata_action._extract_exif("image.jpg", data) + assert d["projection_type"] == "brown" diff --git a/opensfm/test/test_features_sift_gpu.py b/opensfm/test/test_features_sift_gpu.py new file mode 100644 index 000000000..39ebfe79e --- /dev/null +++ b/opensfm/test/test_features_sift_gpu.py @@ -0,0 +1,43 @@ +# pyre-strict +"""Test that feature_type SIFT_GPU dispatches to pypopsift (mocked).""" + +import sys +import types +from typing import Any, Tuple + +import numpy as np +import pytest +from numpy.typing import NDArray +from opensfm import features + + +def test_sift_gpu_dispatch(monkeypatch: pytest.MonkeyPatch) -> None: + sift_gpu_called: bool = False + + def fake_popsift( + image: NDArray[Any], **kwargs: Any + ) -> Tuple[NDArray[Any], NDArray[Any]]: + nonlocal sift_gpu_called + sift_gpu_called = True + return ( + np.array([[4.0, 4.0, 1.0, 0.0]]), # one keypoint: x, y, size, angle + np.array([[1, 2, 3]], dtype=np.uint8), # its descriptor + ) + + module = types.ModuleType("opensfm.pypopsift") + module.popsift = fake_popsift # type: ignore[attr-defined] + monkeypatch.setitem(sys.modules, "opensfm.pypopsift", module) + + config = { + "feature_type": "sift_gpu", # lowercase on purpose, must be normalized + "feature_root": False, + "feature_process_size": 2048, + "feature_min_frames": 100, + "sift_edge_threshold": 10, + "sift_peak_threshold": 0.02, + } + image = np.zeros((9, 9, 3), dtype=np.uint8) + + features.extract_features(image, config, False) + + assert sift_gpu_called diff --git a/opensfm/test/test_stats.py b/opensfm/test/test_stats.py index c49635cc0..de00c3ca7 100644 --- a/opensfm/test/test_stats.py +++ b/opensfm/test/test_stats.py @@ -17,7 +17,7 @@ def test_processing_statistics_normal( processing_statistics = stats.processing_statistics(dataset, [reference]) assert list(processing_statistics.keys()) == [ - "steps_times", "date", "area"] + "steps_times", "date", "start_date", "end_date", "area"] assert processing_statistics["steps_times"] == { "Feature Extraction": -1, "Features Matching": -1, @@ -42,7 +42,7 @@ def test_processing_statistics_null( processing_statistics = stats.processing_statistics(dataset, [null_scene]) assert list(processing_statistics.keys()) == [ - "steps_times", "date", "area"] + "steps_times", "date", "start_date", "end_date", "area"] assert processing_statistics["steps_times"] == { "Feature Extraction": -1, "Features Matching": -1, diff --git a/opensfm/test/test_undistort.py b/opensfm/test/test_undistort.py index 4dfae75cc..56c1e79cb 100644 --- a/opensfm/test/test_undistort.py +++ b/opensfm/test/test_undistort.py @@ -3,6 +3,7 @@ import numpy as np from opensfm import pygeometry, pymap, types, undistort +from unittest.mock import MagicMock def test_perspective_views_of_a_panorama() -> None: @@ -193,3 +194,29 @@ def test_validity_from_remap_no_remap() -> None: result = undistort._validity_from_remap({}, 20, 20, [shot], 1000) assert "s1" in result assert np.all(result["s1"] == 255) + +def test_image_filter_is_applied() -> None: + """Check that imageFilter is called during undistortion.""" + rec = types.Reconstruction() + shot = _make_undistorted_shot(rec, "s1") + + data = MagicMock() + data.config = {"undistorted_image_max_size": 1000} + original_image = np.ones((10, 10, 3), dtype=np.uint8) + data.load_image.return_value = original_image + data.load_mask.return_value = None + data.load_segmentation.return_value = None + data.load_detection.return_value = None + + udata = MagicMock() + + filter_called = False + def mock_filter(shot_id, img): + nonlocal filter_called + filter_called = True + assert shot_id == "s1" + return img * 2 + + arguments = (shot, [shot], data, udata, mock_filter, {}) + undistort.undistort_image_and_masks(arguments) + assert filter_called diff --git a/opensfm/tracking.py b/opensfm/tracking.py index e83dc6e76..08fbc0822 100644 --- a/opensfm/tracking.py +++ b/opensfm/tracking.py @@ -292,6 +292,41 @@ def process_pair(pair): return common_tracks +def np_all_common_tracks_with_features( + tracks_manager: pymap.TracksManager, + min_common: int = 50, +): + """Same as all_common_tracks, but as a numpy structure + dictionary for indexing""" + all_pairs = tracks_manager.get_all_pairs_connectivity().items() + common_tracks_data = [] + common_tracks_index = {} + rec_start = 0 + for pair, size in all_pairs: + if size < min_common: + continue + + tuples = tracks_manager.get_all_common_observations(pair[0], pair[1]) + num_tracks = len(tuples) + rec_end = rec_start + num_tracks * 8 + + for j in range(num_tracks): + pair_track, p1, p2 = tuples[j] + + common_tracks_data.append(float(pair_track)) + common_tracks_data.append(p1.point[0]) + common_tracks_data.append(p1.point[1]) + common_tracks_data.append(1.0) # Homogeneous + common_tracks_data.append(p2.point[0]) + common_tracks_data.append(p2.point[1]) + common_tracks_data.append(1.0) # Homogeneous + common_tracks_data.append(0.0) # Align to 8bit boundary + + common_tracks_index[pair] = (num_tracks, rec_start, rec_end) + rec_start = rec_end + + return np.reshape(common_tracks_data, newshape=(len(common_tracks_data), )), common_tracks_index + + def _good_track(track: List[Tuple[str, int]], min_length: int) -> bool: if len(track) < min_length: return False diff --git a/opensfm/undistort.py b/opensfm/undistort.py index dcdb1650e..12e759069 100644 --- a/opensfm/undistort.py +++ b/opensfm/undistort.py @@ -1,7 +1,8 @@ # pyre-strict import itertools import logging -from typing import Dict, Iterator, List, Optional, Tuple +import os +from typing import Callable, Dict, Iterator, List, Optional, Tuple import cv2 import numpy as np @@ -21,6 +22,9 @@ logger: logging.Logger = logging.getLogger(__name__) +ImageFilterType = Callable[[str, NDArray], NDArray] +CameraMappingsType = Dict[str, Tuple[NDArray, NDArray]] + def undistort_reconstruction( tracks_manager: Optional[pymap.TracksManager], @@ -86,20 +90,92 @@ def undistort_reconstruction( return undistorted_shots +def dump_camera_mapping_cache(camera_mappings_dict: Dict[str, Tuple[NDArray, NDArray]], dest_file: str) -> None: + ids = [] + outs = {} + + def cast_to_max(arr: NDArray) -> NDArray: + m = np.max(arr) + if m <= 255: + return arr.astype(np.uint8) + elif m <= 65535: + return arr.astype(np.uint16) + else: + return arr.astype(np.uint32) + + # Dump the camera mappings information for each camera + # compressing the numbers by rounding to nearest int * mul + # and then by attempting to offset the values + # the original value of the map can then be computed with: + # map1[px_y, px_x] <--> ((compressed_map1[px_y,px_x] + offset[0]) / 10.0) + px_x + # map2[px_y, px_x] <--> ((compressed_map2[px_y,px_x] + offset[1]) / 10.0) + px_y + # (note our values are rounded to 1/10th of a pixel) + idx = 0 + for key, (map1, map2) in camera_mappings_dict.items(): + ids.append(key) + map1 = np.round(map1).astype(np.float32) + map2 = np.round(map2).astype(np.float32) + mul = 10.0 # keep precision up to 1/10th of a pixel + + i = np.arange(map1.shape[0]).reshape(-1, 1) + j = np.arange(map1.shape[1]) + map1 -= j + map2 -= i + map1 = np.round(map1 * mul) + map2 = np.round(map2 * mul) + + offset = np.array([np.min(map1), np.min(map2)]) + map1 -= offset[0] + map2 -= offset[1] + + map1 = cast_to_max(map1) + map2 = cast_to_max(map2) + + outs['%s_x' % idx] = map1 + outs['%s_y' % idx] = map2 + outs['%s_offset' % idx] = offset + outs['%s_mul' % idx] = np.array([mul]) + + idx += 1 + + np.savez_compressed(dest_file, ids=ids, **outs) + + def undistort_reconstruction_with_images( tracks_manager: Optional[pymap.TracksManager], reconstruction: types.Reconstruction, data: DataSetBase, udata: UndistortedDataSet, + imageFilter: Optional[ImageFilterType] = None, skip_images: bool = False, ) -> Dict[pymap.Shot, List[pymap.Shot]]: undistorted_shots = undistort_reconstruction( tracks_manager, reconstruction, data, udata ) if not skip_images: + # Precompute camera mappings for the .npz dump + camera_mappings_dict = {} + for shot_id, subshots in undistorted_shots.items(): + shot = reconstruction.shots[shot_id] + cam_id = shot.camera.id + if cam_id not in camera_mappings_dict and shot.camera.projection_type in _REMAP_CAMERA_TYPES: + [undistorted_shot] = subshots + height, width = shot.camera.height, shot.camera.width + map1, map2 = pygeometry.compute_camera_mapping( + shot.camera, undistorted_shot.camera, width, height + ) + camera_mappings_dict[cam_id] = (map1, map2) + + if camera_mappings_dict: + dump_camera_mapping_cache(camera_mappings_dict, os.path.join(data.data_path, "camera_mappings.npz")) + arguments = [] for shot_id, subshots in undistorted_shots.items(): - arguments.append((reconstruction.shots[shot_id], subshots, data, udata)) + arguments.append((reconstruction.shots[shot_id], subshots, data, udata, imageFilter, camera_mappings_dict)) + + # Sort arguments by camera ID so that joblib batches images from the same camera together. + # This maximizes cache hits in the worker processes and minimizes memory usage. + arguments.sort(key=lambda x: x[0].camera.id) processes = data.config["processes"] @@ -126,9 +202,9 @@ def undistort_reconstruction_with_images( def undistort_image_and_masks( - arguments: Tuple[pymap.Shot, List[pymap.Shot], DataSetBase, UndistortedDataSet], + arguments: Tuple[pymap.Shot, List[pymap.Shot], DataSetBase, UndistortedDataSet, Optional[ImageFilterType], CameraMappingsType], ) -> None: - shot, undistorted_shots, data, udata = arguments + shot, undistorted_shots, data, udata, imageFilter, camera_mappings_dict = arguments log.setup() logger.debug("Undistorting image {}".format(shot.id)) max_size = data.config["undistorted_image_max_size"] @@ -136,8 +212,15 @@ def undistort_image_and_masks( # Undistort image image = data.load_image(shot.id, unchanged=True, anydepth=True) if image is not None: + if imageFilter is not None: + image = imageFilter(shot.id, image) + height, width = image.shape[:2] - remap = _compute_camera_remap(shot, undistorted_shots, height, width) + remap = None + if shot.camera.id in camera_mappings_dict: + remap = {undistorted_shots[0].id: camera_mappings_dict[shot.camera.id]} + else: + remap = _compute_camera_remap(shot, undistorted_shots, height, width) undistorted = undistort_image( shot, undistorted_shots, image, cv2.INTER_AREA, max_size, @@ -178,7 +261,7 @@ def undistort_image( original: Optional[NDArray], interpolation: int, max_size: int, - remap: Optional[Dict[str, Tuple[NDArray, NDArray]]] = None, + remap: Optional[CameraMappingsType] = None, ) -> Dict[str, NDArray]: """Undistort an image into a set of undistorted ones. @@ -254,7 +337,7 @@ def _compute_camera_remap( undistorted_shots: List[pymap.Shot], height: int, width: int, -) -> Dict[str, Tuple[NDArray, NDArray]]: +) -> CameraMappingsType: """Compute per-undistorted-view remap coordinate maps. For standard camera types returns ``{shot_id: (map1, map2)}``. diff --git a/pyproject.toml b/pyproject.toml index 157916417..2d345d142 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -14,12 +14,14 @@ authors = [ ] dependencies = [ + "beautifulsoup4>=4.9.1", "cloudpickle>=0.4.0", "exifread>=2.1.2", "flask>=2.3.2", "fpdf2>=2.4.6", "igraph>=1.0.0", "joblib>=1.0.0", + "lxml>=4.5.1", "lz4>=4.4.5", "matplotlib", "networkx>=2.5", @@ -29,7 +31,10 @@ dependencies = [ "pyproj>=1.9.5.1", "python-dateutil>=2.7", "pyyaml>=5.4", + "rasterio>=1.3", + "rawpy>=0.18", "scipy>=1.10.0", + "vmem>=1.0.2", "xmltodict>=0.10.2", "opencv-python>=4.8.0", "rerun-sdk==0.27.2", diff --git a/viewer/legacy/reconstruction.html b/viewer/legacy/reconstruction.html index 229c5d897..cbac7b4de 100644 --- a/viewer/legacy/reconstruction.html +++ b/viewer/legacy/reconstruction.html @@ -317,6 +317,7 @@ var reconstruction_groups = []; var point_clouds = []; var camera_lines = []; + var camera_planes = []; var gps_lines = []; var imagePlanes = []; var imagePlaneCameras = []; @@ -337,8 +338,10 @@ }; var options = { - cameraSize: 0.9, - pointSize: 0.7, + cameraSize: 0.07, + pointSize: 0.04, + // cameraSize: 0.9, + // pointSize: 0.7, imagePlaneSize: 50, showThumbnail: true, showImagePlane: false, @@ -667,6 +670,41 @@ return lines; } + function createPlane(id, shot){ + var geometry = new THREE.PlaneGeometry( 2, 2 ); + var plane = new THREE.Mesh( geometry, new THREE.MeshBasicMaterial( {color: reconstruction_color(id), opacity: 0.5, transparent: true, side: THREE.DoubleSide} ) ); + var oc = opticalCenter(shot); + oc.setZ(shot["plane"]["depth"]); + var la = new THREE.Vector3().addVectors(oc, new THREE.Vector3(shot["plane"]["normal"][0], shot["plane"]["normal"][1], shot["plane"]["normal"][2])); + plane.position.copy(oc); + plane.lookAt(la); + return plane; + } + + function initCameraPlanes(reconstruction) { + var planes = [] + var i = 0; + for (var shot_id in reconstruction.shots) { + if (reconstruction.shots.hasOwnProperty(shot_id)) { + var shot = reconstruction.shots[shot_id]; + if (shot["plane"]){ + planes.push(createPlane(i++, shot)); + } + } + } + // Main plane + if (reconstruction.plane){ + var geometry = new THREE.PlaneGeometry( 120, 120 ); + var plane = new THREE.Mesh( geometry, new THREE.MeshBasicMaterial( {color: reconstruction_color(99), opacity: 0.3, transparent: true, side: THREE.DoubleSide} ) ); + var center = new THREE.Vector3(0, 0, reconstruction["plane"]["depth"]); + var la = new THREE.Vector3().addVectors(center, new THREE.Vector3(reconstruction["plane"]["normal"][0], reconstruction["plane"]["normal"][1], reconstruction["plane"]["normal"][2])); + plane.position.copy(center); + plane.lookAt(la); + planes.push(plane); + } + return planes; + } + function updateCameraLines() { for (var i = 0; i < camera_lines.length; ++i) { var linegeo = cameraLineGeo(camera_lines[i].reconstruction, camera_lines[i].shot_id); @@ -873,7 +911,7 @@ return rotation; } - function reconstruction_color(id){ + function reconstruction_color(id, opacity){ var s = 1 if(id !== 0){ s = Math.sin(id) * 10000; @@ -881,6 +919,7 @@ } var color = new THREE.Color(); color.setHex( s * 0xffffff ); + if (opacity !== undefined) color.opacity = opacity; return color; } @@ -943,6 +982,13 @@ camera_lines.push(lines[i]); } + // Planes + var planes = initCameraPlanes(reconstruction); + for (var i = 0; i < planes.length; ++i) { + group.add(planes[i]); + camera_planes.push(planes[i]); + } + // GPS positions for (var shot_id in reconstruction.shots) { if (reconstruction.shots.hasOwnProperty(shot_id)) { diff --git a/viewer/src/ui/OpensfmViewer.js b/viewer/src/ui/OpensfmViewer.js index 8c0e529c7..eb6964d22 100644 --- a/viewer/src/ui/OpensfmViewer.js +++ b/viewer/src/ui/OpensfmViewer.js @@ -38,12 +38,12 @@ export class OpensfmViewer extends EventEmitter { const cvm = CameraVisualizationMode.Homogeneous; const opm = OriginalPositionMode.Hidden; const spatialConfiguration = { - cameraSize: 0.5, + cameraSize: 0.1,//0.5, cameraVisualizationMode: cvm, cellGridDepth: 3, cellsVisible: false, originalPositionMode: opm, - pointSize: 0.2, + pointSize: 0.01,//0.2, pointsVisible: true, }; @@ -86,7 +86,7 @@ export class OpensfmViewer extends EventEmitter { const thumbnailVisible = false; const controllerOptions = { - axesVisible: true, + axesVisible: false, cameraControlMode, commandsVisible, gridVisible: true, @@ -145,7 +145,7 @@ export class OpensfmViewer extends EventEmitter { mode: cameraControlMode, }); this._customRenderer = new CustomRenderer(this._viewer); - this._customRenderer.add(this._axesRenderer); + //this._customRenderer.add(this._axesRenderer); this._customRenderer.add(this._earthRenderer); this._viewer.addCustomRenderer(this._customRenderer); }