diff --git a/LICENSE b/clojure/LICENSE similarity index 100% rename from LICENSE rename to clojure/LICENSE diff --git a/README.md b/clojure/README.md similarity index 100% rename from README.md rename to clojure/README.md diff --git a/project.clj b/clojure/project.clj similarity index 100% rename from project.clj rename to clojure/project.clj diff --git a/resources/public/index.html b/clojure/resources/public/index.html similarity index 100% rename from resources/public/index.html rename to clojure/resources/public/index.html diff --git a/resources/public/js/main.js b/clojure/resources/public/js/main.js similarity index 100% rename from resources/public/js/main.js rename to clojure/resources/public/js/main.js diff --git a/resources/public/js/optimized/cljs/core.cljs b/clojure/resources/public/js/optimized/cljs/core.cljs similarity index 100% rename from resources/public/js/optimized/cljs/core.cljs rename to clojure/resources/public/js/optimized/cljs/core.cljs diff --git a/resources/public/js/optimized/cljs/core.js b/clojure/resources/public/js/optimized/cljs/core.js similarity index 100% rename from resources/public/js/optimized/cljs/core.js rename to clojure/resources/public/js/optimized/cljs/core.js diff --git a/resources/public/js/optimized/cljs/core/constants.js b/clojure/resources/public/js/optimized/cljs/core/constants.js similarity index 100% rename from resources/public/js/optimized/cljs/core/constants.js rename to clojure/resources/public/js/optimized/cljs/core/constants.js diff --git a/resources/public/js/optimized/cljsc_opts.edn b/clojure/resources/public/js/optimized/cljsc_opts.edn similarity index 100% rename from resources/public/js/optimized/cljsc_opts.edn rename to clojure/resources/public/js/optimized/cljsc_opts.edn diff --git a/resources/public/js/optimized/cljsjs/p5/production/p5.min.inc.js b/clojure/resources/public/js/optimized/cljsjs/p5/production/p5.min.inc.js similarity index 100% rename from resources/public/js/optimized/cljsjs/p5/production/p5.min.inc.js rename to clojure/resources/public/js/optimized/cljsjs/p5/production/p5.min.inc.js diff --git a/resources/public/js/optimized/clojure/string.cljs b/clojure/resources/public/js/optimized/clojure/string.cljs similarity index 100% rename from resources/public/js/optimized/clojure/string.cljs rename to clojure/resources/public/js/optimized/clojure/string.cljs diff --git a/resources/public/js/optimized/clojure/string.js b/clojure/resources/public/js/optimized/clojure/string.js similarity index 100% rename from resources/public/js/optimized/clojure/string.js rename to clojure/resources/public/js/optimized/clojure/string.js diff --git a/resources/public/js/optimized/curlnoise/core.js b/clojure/resources/public/js/optimized/curlnoise/core.js similarity index 100% rename from resources/public/js/optimized/curlnoise/core.js rename to clojure/resources/public/js/optimized/curlnoise/core.js diff --git a/resources/public/js/optimized/goog/array/array.js b/clojure/resources/public/js/optimized/goog/array/array.js similarity index 100% rename from resources/public/js/optimized/goog/array/array.js rename to clojure/resources/public/js/optimized/goog/array/array.js diff --git a/resources/public/js/optimized/goog/asserts/asserts.js b/clojure/resources/public/js/optimized/goog/asserts/asserts.js similarity index 100% rename from resources/public/js/optimized/goog/asserts/asserts.js rename to clojure/resources/public/js/optimized/goog/asserts/asserts.js diff --git a/resources/public/js/optimized/goog/base.js b/clojure/resources/public/js/optimized/goog/base.js similarity index 100% rename from resources/public/js/optimized/goog/base.js rename to clojure/resources/public/js/optimized/goog/base.js diff --git a/resources/public/js/optimized/goog/debug/entrypointregistry.js b/clojure/resources/public/js/optimized/goog/debug/entrypointregistry.js similarity index 100% rename from resources/public/js/optimized/goog/debug/entrypointregistry.js rename to clojure/resources/public/js/optimized/goog/debug/entrypointregistry.js diff --git a/resources/public/js/optimized/goog/debug/error.js b/clojure/resources/public/js/optimized/goog/debug/error.js similarity index 100% rename from resources/public/js/optimized/goog/debug/error.js rename to clojure/resources/public/js/optimized/goog/debug/error.js diff --git a/resources/public/js/optimized/goog/disposable/disposable.js b/clojure/resources/public/js/optimized/goog/disposable/disposable.js similarity index 100% rename from resources/public/js/optimized/goog/disposable/disposable.js rename to clojure/resources/public/js/optimized/goog/disposable/disposable.js diff --git a/resources/public/js/optimized/goog/disposable/idisposable.js b/clojure/resources/public/js/optimized/goog/disposable/idisposable.js similarity index 100% rename from resources/public/js/optimized/goog/disposable/idisposable.js rename to clojure/resources/public/js/optimized/goog/disposable/idisposable.js diff --git a/resources/public/js/optimized/goog/dom/asserts.js b/clojure/resources/public/js/optimized/goog/dom/asserts.js similarity index 100% rename from resources/public/js/optimized/goog/dom/asserts.js rename to clojure/resources/public/js/optimized/goog/dom/asserts.js diff --git a/resources/public/js/optimized/goog/dom/browserfeature.js b/clojure/resources/public/js/optimized/goog/dom/browserfeature.js similarity index 100% rename from resources/public/js/optimized/goog/dom/browserfeature.js rename to clojure/resources/public/js/optimized/goog/dom/browserfeature.js diff --git a/resources/public/js/optimized/goog/dom/dom.js b/clojure/resources/public/js/optimized/goog/dom/dom.js similarity index 100% rename from resources/public/js/optimized/goog/dom/dom.js rename to clojure/resources/public/js/optimized/goog/dom/dom.js diff --git a/resources/public/js/optimized/goog/dom/htmlelement.js b/clojure/resources/public/js/optimized/goog/dom/htmlelement.js similarity index 100% rename from resources/public/js/optimized/goog/dom/htmlelement.js rename to clojure/resources/public/js/optimized/goog/dom/htmlelement.js diff --git a/resources/public/js/optimized/goog/dom/nodetype.js b/clojure/resources/public/js/optimized/goog/dom/nodetype.js similarity index 100% rename from resources/public/js/optimized/goog/dom/nodetype.js rename to clojure/resources/public/js/optimized/goog/dom/nodetype.js diff --git a/resources/public/js/optimized/goog/dom/safe.js b/clojure/resources/public/js/optimized/goog/dom/safe.js similarity index 100% rename from resources/public/js/optimized/goog/dom/safe.js rename to clojure/resources/public/js/optimized/goog/dom/safe.js diff --git a/resources/public/js/optimized/goog/dom/tagname.js b/clojure/resources/public/js/optimized/goog/dom/tagname.js similarity index 100% rename from resources/public/js/optimized/goog/dom/tagname.js rename to clojure/resources/public/js/optimized/goog/dom/tagname.js diff --git a/resources/public/js/optimized/goog/dom/tags.js b/clojure/resources/public/js/optimized/goog/dom/tags.js similarity index 100% rename from resources/public/js/optimized/goog/dom/tags.js rename to clojure/resources/public/js/optimized/goog/dom/tags.js diff --git a/resources/public/js/optimized/goog/dom/vendor.js b/clojure/resources/public/js/optimized/goog/dom/vendor.js similarity index 100% rename from resources/public/js/optimized/goog/dom/vendor.js rename to clojure/resources/public/js/optimized/goog/dom/vendor.js diff --git a/resources/public/js/optimized/goog/events/browserevent.js b/clojure/resources/public/js/optimized/goog/events/browserevent.js similarity index 100% rename from resources/public/js/optimized/goog/events/browserevent.js rename to clojure/resources/public/js/optimized/goog/events/browserevent.js diff --git a/resources/public/js/optimized/goog/events/browserfeature.js b/clojure/resources/public/js/optimized/goog/events/browserfeature.js similarity index 100% rename from resources/public/js/optimized/goog/events/browserfeature.js rename to clojure/resources/public/js/optimized/goog/events/browserfeature.js diff --git a/resources/public/js/optimized/goog/events/event.js b/clojure/resources/public/js/optimized/goog/events/event.js similarity index 100% rename from resources/public/js/optimized/goog/events/event.js rename to clojure/resources/public/js/optimized/goog/events/event.js diff --git a/resources/public/js/optimized/goog/events/eventid.js b/clojure/resources/public/js/optimized/goog/events/eventid.js similarity index 100% rename from resources/public/js/optimized/goog/events/eventid.js rename to clojure/resources/public/js/optimized/goog/events/eventid.js diff --git a/resources/public/js/optimized/goog/events/events.js b/clojure/resources/public/js/optimized/goog/events/events.js similarity index 100% rename from resources/public/js/optimized/goog/events/events.js rename to clojure/resources/public/js/optimized/goog/events/events.js diff --git a/resources/public/js/optimized/goog/events/eventtype.js b/clojure/resources/public/js/optimized/goog/events/eventtype.js similarity index 100% rename from resources/public/js/optimized/goog/events/eventtype.js rename to clojure/resources/public/js/optimized/goog/events/eventtype.js diff --git a/resources/public/js/optimized/goog/events/listenable.js b/clojure/resources/public/js/optimized/goog/events/listenable.js similarity index 100% rename from resources/public/js/optimized/goog/events/listenable.js rename to clojure/resources/public/js/optimized/goog/events/listenable.js diff --git a/resources/public/js/optimized/goog/events/listener.js b/clojure/resources/public/js/optimized/goog/events/listener.js similarity index 100% rename from resources/public/js/optimized/goog/events/listener.js rename to clojure/resources/public/js/optimized/goog/events/listener.js diff --git a/resources/public/js/optimized/goog/events/listenermap.js b/clojure/resources/public/js/optimized/goog/events/listenermap.js similarity index 100% rename from resources/public/js/optimized/goog/events/listenermap.js rename to clojure/resources/public/js/optimized/goog/events/listenermap.js diff --git a/resources/public/js/optimized/goog/fs/url.js b/clojure/resources/public/js/optimized/goog/fs/url.js similarity index 100% rename from resources/public/js/optimized/goog/fs/url.js rename to clojure/resources/public/js/optimized/goog/fs/url.js diff --git a/resources/public/js/optimized/goog/functions/functions.js b/clojure/resources/public/js/optimized/goog/functions/functions.js similarity index 100% rename from resources/public/js/optimized/goog/functions/functions.js rename to clojure/resources/public/js/optimized/goog/functions/functions.js diff --git a/resources/public/js/optimized/goog/html/legacyconversions.js b/clojure/resources/public/js/optimized/goog/html/legacyconversions.js similarity index 100% rename from resources/public/js/optimized/goog/html/legacyconversions.js rename to clojure/resources/public/js/optimized/goog/html/legacyconversions.js diff --git a/resources/public/js/optimized/goog/html/safehtml.js b/clojure/resources/public/js/optimized/goog/html/safehtml.js similarity index 100% rename from resources/public/js/optimized/goog/html/safehtml.js rename to clojure/resources/public/js/optimized/goog/html/safehtml.js diff --git a/resources/public/js/optimized/goog/html/safescript.js b/clojure/resources/public/js/optimized/goog/html/safescript.js similarity index 100% rename from resources/public/js/optimized/goog/html/safescript.js rename to clojure/resources/public/js/optimized/goog/html/safescript.js diff --git a/resources/public/js/optimized/goog/html/safestyle.js b/clojure/resources/public/js/optimized/goog/html/safestyle.js similarity index 100% rename from resources/public/js/optimized/goog/html/safestyle.js rename to clojure/resources/public/js/optimized/goog/html/safestyle.js diff --git a/resources/public/js/optimized/goog/html/safestylesheet.js b/clojure/resources/public/js/optimized/goog/html/safestylesheet.js similarity index 100% rename from resources/public/js/optimized/goog/html/safestylesheet.js rename to clojure/resources/public/js/optimized/goog/html/safestylesheet.js diff --git a/resources/public/js/optimized/goog/html/safeurl.js b/clojure/resources/public/js/optimized/goog/html/safeurl.js similarity index 100% rename from resources/public/js/optimized/goog/html/safeurl.js rename to clojure/resources/public/js/optimized/goog/html/safeurl.js diff --git a/resources/public/js/optimized/goog/html/trustedresourceurl.js b/clojure/resources/public/js/optimized/goog/html/trustedresourceurl.js similarity index 100% rename from resources/public/js/optimized/goog/html/trustedresourceurl.js rename to clojure/resources/public/js/optimized/goog/html/trustedresourceurl.js diff --git a/resources/public/js/optimized/goog/html/uncheckedconversions.js b/clojure/resources/public/js/optimized/goog/html/uncheckedconversions.js similarity index 100% rename from resources/public/js/optimized/goog/html/uncheckedconversions.js rename to clojure/resources/public/js/optimized/goog/html/uncheckedconversions.js diff --git a/resources/public/js/optimized/goog/i18n/bidi.js b/clojure/resources/public/js/optimized/goog/i18n/bidi.js similarity index 100% rename from resources/public/js/optimized/goog/i18n/bidi.js rename to clojure/resources/public/js/optimized/goog/i18n/bidi.js diff --git a/resources/public/js/optimized/goog/iter/iter.js b/clojure/resources/public/js/optimized/goog/iter/iter.js similarity index 100% rename from resources/public/js/optimized/goog/iter/iter.js rename to clojure/resources/public/js/optimized/goog/iter/iter.js diff --git a/resources/public/js/optimized/goog/labs/useragent/browser.js b/clojure/resources/public/js/optimized/goog/labs/useragent/browser.js similarity index 100% rename from resources/public/js/optimized/goog/labs/useragent/browser.js rename to clojure/resources/public/js/optimized/goog/labs/useragent/browser.js diff --git a/resources/public/js/optimized/goog/labs/useragent/engine.js b/clojure/resources/public/js/optimized/goog/labs/useragent/engine.js similarity index 100% rename from resources/public/js/optimized/goog/labs/useragent/engine.js rename to clojure/resources/public/js/optimized/goog/labs/useragent/engine.js diff --git a/resources/public/js/optimized/goog/labs/useragent/platform.js b/clojure/resources/public/js/optimized/goog/labs/useragent/platform.js similarity index 100% rename from resources/public/js/optimized/goog/labs/useragent/platform.js rename to clojure/resources/public/js/optimized/goog/labs/useragent/platform.js diff --git a/resources/public/js/optimized/goog/labs/useragent/util.js b/clojure/resources/public/js/optimized/goog/labs/useragent/util.js similarity index 100% rename from resources/public/js/optimized/goog/labs/useragent/util.js rename to clojure/resources/public/js/optimized/goog/labs/useragent/util.js diff --git a/resources/public/js/optimized/goog/math/box.js b/clojure/resources/public/js/optimized/goog/math/box.js similarity index 100% rename from resources/public/js/optimized/goog/math/box.js rename to clojure/resources/public/js/optimized/goog/math/box.js diff --git a/resources/public/js/optimized/goog/math/coordinate.js b/clojure/resources/public/js/optimized/goog/math/coordinate.js similarity index 100% rename from resources/public/js/optimized/goog/math/coordinate.js rename to clojure/resources/public/js/optimized/goog/math/coordinate.js diff --git a/resources/public/js/optimized/goog/math/integer.js b/clojure/resources/public/js/optimized/goog/math/integer.js similarity index 100% rename from resources/public/js/optimized/goog/math/integer.js rename to clojure/resources/public/js/optimized/goog/math/integer.js diff --git a/resources/public/js/optimized/goog/math/irect.js b/clojure/resources/public/js/optimized/goog/math/irect.js similarity index 100% rename from resources/public/js/optimized/goog/math/irect.js rename to clojure/resources/public/js/optimized/goog/math/irect.js diff --git a/resources/public/js/optimized/goog/math/long.js b/clojure/resources/public/js/optimized/goog/math/long.js similarity index 100% rename from resources/public/js/optimized/goog/math/long.js rename to clojure/resources/public/js/optimized/goog/math/long.js diff --git a/resources/public/js/optimized/goog/math/math.js b/clojure/resources/public/js/optimized/goog/math/math.js similarity index 100% rename from resources/public/js/optimized/goog/math/math.js rename to clojure/resources/public/js/optimized/goog/math/math.js diff --git a/resources/public/js/optimized/goog/math/rect.js b/clojure/resources/public/js/optimized/goog/math/rect.js similarity index 100% rename from resources/public/js/optimized/goog/math/rect.js rename to clojure/resources/public/js/optimized/goog/math/rect.js diff --git a/resources/public/js/optimized/goog/math/size.js b/clojure/resources/public/js/optimized/goog/math/size.js similarity index 100% rename from resources/public/js/optimized/goog/math/size.js rename to clojure/resources/public/js/optimized/goog/math/size.js diff --git a/resources/public/js/optimized/goog/object/object.js b/clojure/resources/public/js/optimized/goog/object/object.js similarity index 100% rename from resources/public/js/optimized/goog/object/object.js rename to clojure/resources/public/js/optimized/goog/object/object.js diff --git a/resources/public/js/optimized/goog/reflect/reflect.js b/clojure/resources/public/js/optimized/goog/reflect/reflect.js similarity index 100% rename from resources/public/js/optimized/goog/reflect/reflect.js rename to clojure/resources/public/js/optimized/goog/reflect/reflect.js diff --git a/resources/public/js/optimized/goog/string/const.js b/clojure/resources/public/js/optimized/goog/string/const.js similarity index 100% rename from resources/public/js/optimized/goog/string/const.js rename to clojure/resources/public/js/optimized/goog/string/const.js diff --git a/resources/public/js/optimized/goog/string/string.js b/clojure/resources/public/js/optimized/goog/string/string.js similarity index 100% rename from resources/public/js/optimized/goog/string/string.js rename to clojure/resources/public/js/optimized/goog/string/string.js diff --git a/resources/public/js/optimized/goog/string/stringbuffer.js b/clojure/resources/public/js/optimized/goog/string/stringbuffer.js similarity index 100% rename from resources/public/js/optimized/goog/string/stringbuffer.js rename to clojure/resources/public/js/optimized/goog/string/stringbuffer.js diff --git a/resources/public/js/optimized/goog/string/typedstring.js b/clojure/resources/public/js/optimized/goog/string/typedstring.js similarity index 100% rename from resources/public/js/optimized/goog/string/typedstring.js rename to clojure/resources/public/js/optimized/goog/string/typedstring.js diff --git a/resources/public/js/optimized/goog/structs/map.js b/clojure/resources/public/js/optimized/goog/structs/map.js similarity index 100% rename from resources/public/js/optimized/goog/structs/map.js rename to clojure/resources/public/js/optimized/goog/structs/map.js diff --git a/resources/public/js/optimized/goog/structs/structs.js b/clojure/resources/public/js/optimized/goog/structs/structs.js similarity index 100% rename from resources/public/js/optimized/goog/structs/structs.js rename to clojure/resources/public/js/optimized/goog/structs/structs.js diff --git a/resources/public/js/optimized/goog/style/style.js b/clojure/resources/public/js/optimized/goog/style/style.js similarity index 100% rename from resources/public/js/optimized/goog/style/style.js rename to clojure/resources/public/js/optimized/goog/style/style.js diff --git a/resources/public/js/optimized/goog/uri/uri.js b/clojure/resources/public/js/optimized/goog/uri/uri.js similarity index 100% rename from resources/public/js/optimized/goog/uri/uri.js rename to clojure/resources/public/js/optimized/goog/uri/uri.js diff --git a/resources/public/js/optimized/goog/uri/utils.js b/clojure/resources/public/js/optimized/goog/uri/utils.js similarity index 100% rename from resources/public/js/optimized/goog/uri/utils.js rename to clojure/resources/public/js/optimized/goog/uri/utils.js diff --git a/resources/public/js/optimized/goog/useragent/useragent.js b/clojure/resources/public/js/optimized/goog/useragent/useragent.js similarity index 100% rename from resources/public/js/optimized/goog/useragent/useragent.js rename to clojure/resources/public/js/optimized/goog/useragent/useragent.js diff --git a/resources/public/js/optimized/process/env.cljs b/clojure/resources/public/js/optimized/process/env.cljs similarity index 100% rename from resources/public/js/optimized/process/env.cljs rename to clojure/resources/public/js/optimized/process/env.cljs diff --git a/resources/public/js/optimized/process/env.js b/clojure/resources/public/js/optimized/process/env.js similarity index 100% rename from resources/public/js/optimized/process/env.js rename to clojure/resources/public/js/optimized/process/env.js diff --git a/resources/public/js/optimized/quil/core.cljc b/clojure/resources/public/js/optimized/quil/core.cljc similarity index 100% rename from resources/public/js/optimized/quil/core.cljc rename to clojure/resources/public/js/optimized/quil/core.cljc diff --git a/resources/public/js/optimized/quil/core.js b/clojure/resources/public/js/optimized/quil/core.js similarity index 100% rename from resources/public/js/optimized/quil/core.js rename to clojure/resources/public/js/optimized/quil/core.js diff --git a/resources/public/js/optimized/quil/middleware.cljc b/clojure/resources/public/js/optimized/quil/middleware.cljc similarity index 100% rename from resources/public/js/optimized/quil/middleware.cljc rename to clojure/resources/public/js/optimized/quil/middleware.cljc diff --git a/resources/public/js/optimized/quil/middleware.js b/clojure/resources/public/js/optimized/quil/middleware.js similarity index 100% rename from resources/public/js/optimized/quil/middleware.js rename to clojure/resources/public/js/optimized/quil/middleware.js diff --git a/resources/public/js/optimized/quil/middlewares/deprecated_options.cljc b/clojure/resources/public/js/optimized/quil/middlewares/deprecated_options.cljc similarity index 100% rename from resources/public/js/optimized/quil/middlewares/deprecated_options.cljc rename to clojure/resources/public/js/optimized/quil/middlewares/deprecated_options.cljc diff --git a/resources/public/js/optimized/quil/middlewares/deprecated_options.js b/clojure/resources/public/js/optimized/quil/middlewares/deprecated_options.js similarity index 100% rename from resources/public/js/optimized/quil/middlewares/deprecated_options.js rename to clojure/resources/public/js/optimized/quil/middlewares/deprecated_options.js diff --git a/resources/public/js/optimized/quil/middlewares/fun_mode.cljc b/clojure/resources/public/js/optimized/quil/middlewares/fun_mode.cljc similarity index 100% rename from resources/public/js/optimized/quil/middlewares/fun_mode.cljc rename to clojure/resources/public/js/optimized/quil/middlewares/fun_mode.cljc diff --git a/resources/public/js/optimized/quil/middlewares/fun_mode.js b/clojure/resources/public/js/optimized/quil/middlewares/fun_mode.js similarity index 100% rename from resources/public/js/optimized/quil/middlewares/fun_mode.js rename to clojure/resources/public/js/optimized/quil/middlewares/fun_mode.js diff --git a/resources/public/js/optimized/quil/middlewares/navigation_2d.cljc b/clojure/resources/public/js/optimized/quil/middlewares/navigation_2d.cljc similarity index 100% rename from resources/public/js/optimized/quil/middlewares/navigation_2d.cljc rename to clojure/resources/public/js/optimized/quil/middlewares/navigation_2d.cljc diff --git a/resources/public/js/optimized/quil/middlewares/navigation_2d.js b/clojure/resources/public/js/optimized/quil/middlewares/navigation_2d.js similarity index 100% rename from resources/public/js/optimized/quil/middlewares/navigation_2d.js rename to clojure/resources/public/js/optimized/quil/middlewares/navigation_2d.js diff --git a/resources/public/js/optimized/quil/middlewares/navigation_3d.cljc b/clojure/resources/public/js/optimized/quil/middlewares/navigation_3d.cljc similarity index 100% rename from resources/public/js/optimized/quil/middlewares/navigation_3d.cljc rename to clojure/resources/public/js/optimized/quil/middlewares/navigation_3d.cljc diff --git a/resources/public/js/optimized/quil/middlewares/navigation_3d.js b/clojure/resources/public/js/optimized/quil/middlewares/navigation_3d.js similarity index 100% rename from resources/public/js/optimized/quil/middlewares/navigation_3d.js rename to clojure/resources/public/js/optimized/quil/middlewares/navigation_3d.js diff --git a/resources/public/js/optimized/quil/sketch.cljs b/clojure/resources/public/js/optimized/quil/sketch.cljs similarity index 100% rename from resources/public/js/optimized/quil/sketch.cljs rename to clojure/resources/public/js/optimized/quil/sketch.cljs diff --git a/resources/public/js/optimized/quil/sketch.js b/clojure/resources/public/js/optimized/quil/sketch.js similarity index 100% rename from resources/public/js/optimized/quil/sketch.js rename to clojure/resources/public/js/optimized/quil/sketch.js diff --git a/resources/public/js/optimized/quil/util.cljc b/clojure/resources/public/js/optimized/quil/util.cljc similarity index 100% rename from resources/public/js/optimized/quil/util.cljc rename to clojure/resources/public/js/optimized/quil/util.cljc diff --git a/resources/public/js/optimized/quil/util.js b/clojure/resources/public/js/optimized/quil/util.js similarity index 100% rename from resources/public/js/optimized/quil/util.js rename to clojure/resources/public/js/optimized/quil/util.js diff --git a/src/curlnoise/core.cljc b/clojure/src/curlnoise/core.cljc similarity index 100% rename from src/curlnoise/core.cljc rename to clojure/src/curlnoise/core.cljc diff --git a/python/curl_flow.py b/python/curl_flow.py new file mode 100644 index 0000000..750d9b8 --- /dev/null +++ b/python/curl_flow.py @@ -0,0 +1,117 @@ +# -*- coding: utf-8 -*- + +# Scratch code which (very slowly) animates 3D flow via curl noise. + +import functools + +import numpy as np +import vispy +import vispy.scene +from vispy.scene import visuals + +import opensimplex + +def gradient(fn, eps=1e-3): + eps_inv = 1.0 / eps + # Numerical gradient by finite differences + def _grad(x, y, z): + p = fn(x, y, z) + p_dx = fn(x+eps, y, z) + p_dy = fn(x, y+eps, z) + p_dz = fn(x, y, z+eps) + return (p_dx - p)*eps_inv, (p_dy - p)*eps_inv, (p_dz - p)*eps_inv + return _grad + +def curl_3d(grads): + # 'grads' should be of shape (..., 3, 3); + # 2nd-to-last dimension is gradient of potential x, y, z. + # Last dimension is gradient of that w.r.t. x, y, and z. + cx = grads[..., 2, 1] - grads[..., 1, 2] + cy = grads[..., 0, 2] - grads[..., 2, 0] + cz = grads[..., 1, 0] - grads[..., 0, 1] + return np.stack((cx, cy, cz), axis=-1) + +def generate(grid): + x_spx = opensimplex.OpenSimplex(seed=0) + y_spx = opensimplex.OpenSimplex(seed=12345) + z_spx = opensimplex.OpenSimplex(seed=45678) + grad_x = gradient(x_spx.noise3d) + grad_y = gradient(y_spx.noise3d) + grad_z = gradient(z_spx.noise3d) + # grad_x = gradient(functools.partial(x_spx.noise4d, t)) + # grad_y = gradient(functools.partial(y_spx.noise4d, t)) + # grad_z = gradient(functools.partial(z_spx.noise4d, t)) + # grad_x = gradient(lambda x,y,z: x_spx.noise3d(x + t, y, z)) + # grad_y = gradient(lambda x,y,z: y_spx.noise3d(x + t, y, z)) + # grad_z = gradient(lambda x,y,z: z_spx.noise3d(x + t, y, z)) + grads = np.array([ + (grad_x(x,y,z), grad_y(x,y,z), grad_z(x,y,z)) + for (x,y,z) in grid + ]) + curl = curl_3d(grads) + return curl + +class Data(object): + def __init__(self, view): + self.view = view + self.visual = vispy.scene.visuals.Line( + color=(1,1,1,0.75), + connect='segments', + ) + self.view.add(self.visual) + self.view.camera = 'turntable' # or try 'arcball' + #view.camera = 'arcball' + s = 0.1 + count = 9 + xs = zs = np.linspace(-s, s, count) + ys = np.array([0]) + self.points = np.array([i.flatten() for i in np.meshgrid(xs,ys,zs)]).T + self.points_old = np.zeros((0,3)) + self.update() + def update(self, ev=None): + t = 0 if ev is None else ev.elapsed + # Get velocity for current points: + curl = generate(self.points) + # + a1 = self.points + a2 = self.points + curl*0.005 + # So I guess if I give argument 'pos', it needs to be of shape + # (2*N, 3) and consist of alternating starting & ending vertices? + # The docs don't say a thing about this. + lines = np.hstack((a1, a2)).reshape(self.points.shape[0]*2, -1) + self.points_old = np.vstack((self.points_old, lines)) + maxpoints = self.points.shape[0] * 200 + extra = self.points_old.shape[0] - maxpoints + if extra > 0: + self.points_old = self.points_old[extra:] + self.visual.set_data( + pos=self.points_old, + #arrows=np.hstack((a1, a2)), + ) + self.points = a2 + # self.scatter = visuals.Markers() + # self.scatter.set_data(self.d, edge_color=None, face_color=(1, 0.5, 1, .5), size=4) + # view.add(self.scatter) + #m = np.array([[np.cos(t), np.sin(t*1.01), np.cos(t*1.02)]]) + #d2 = self.d*m + #self.scatter.set_data(d2, edge_color=None, face_color=(1, 0.5, 1, .5), size=4) + +def main(): + # + # Make a canvas and add simple view + # + canvas = vispy.scene.SceneCanvas(keys='interactive', show=True) + view = canvas.central_widget.add_view() + import sys + # add a colored 3D axis for orientation + axis = visuals.XYZAxis(parent=view.scene) + timer = vispy.app.Timer() + da = Data(view) + # Problem is how slow update() is: + timer.connect(da.update) + timer.start(0.05) + if sys.flags.interactive != 1: + vispy.app.run() + +if __name__ == '__main__': + main() diff --git a/python/curlnoise_fibers.py b/python/curlnoise_fibers.py new file mode 100644 index 0000000..afd586d --- /dev/null +++ b/python/curlnoise_fibers.py @@ -0,0 +1,151 @@ +#!/usr/bin/env python3 +# -*- coding: utf-8 -*- + +# (setq python-shell-interpreter (shell-command-to-string "nix-shell --command \"which python3 | tr -d '\n'\"")) + +import functools +import math + +import numpy as np +import vispy +import vispy.scene +from vispy.scene import visuals + +import opensimplex + +def gradient(fn, eps=1e-3): + eps_inv = 1.0 / eps + # Numerical gradient by finite differences + def _grad(x, y, z): + p = fn(x, y, z) + p_dx = fn(x+eps, y, z) + p_dy = fn(x, y+eps, z) + p_dz = fn(x, y, z+eps) + return (p_dx - p)*eps_inv, (p_dy - p)*eps_inv, (p_dz - p)*eps_inv + return _grad + +def curl_3d(grads): + # 'grads' should be of shape (..., 3, 3). Each 3x3 matrix should + # be the Jacobian of the potential. + cx = grads[..., 2, 1] - grads[..., 1, 2] + cy = grads[..., 0, 2] - grads[..., 2, 0] + cz = grads[..., 1, 0] - grads[..., 0, 1] + return np.stack((cx, cy, cz), axis=-1) + +class Potential(object): + def __init__(self): + self.x_spx = opensimplex.OpenSimplex(seed=0) + self.y_spx = opensimplex.OpenSimplex(seed=12345) + self.z_spx = opensimplex.OpenSimplex(seed=45678) + def eval(self, x, y, z): + y2 = y# + 0.5*math.sin(3*x) + 0.5*math.sin(2.5*z) + x2 = x + z2 = z + f1 = np.array([ + self.x_spx.noise3d(x2, y2, z2), + self.y_spx.noise3d(x2, y2, z2), + self.z_spx.noise3d(x2, y2, z2), + ]) + f2 = np.array([z*0.5, 0, 0]) + return f1# + f2 + def grad(self, x, y, z, eps=1e-3): + # Returns gradients in matrix form. + # + # Element (i, j) is gradient of i'th component with respect to + # j'th component, where components are (x,y,z) - i.e. row 0 is + # (d/dx Px, d/dy Px, d/dz Px). I'm fairly sure this is just + # the Jacobian. + p = self.eval(x, y, z) + p_dx = self.eval(x+eps, y, z) + p_dy = self.eval(x, y+eps, z) + p_dz = self.eval(x, y, z+eps) + return (np.stack((p_dx, p_dy, p_dz)) - p).T / eps + +def generate(grid): + p = Potential() + grads = np.array([p.grad(*pt) for pt in grid]) + curl = curl_3d(grads) + return curl + +class Data(object): + def __init__(self, view): + self.use_tubes = False + self.view = view + s = 0.1 + self.s = s + count = 8 + xs = zs = np.linspace(-s, s, count) + ys = np.array([0]) + self.points = np.array([i.flatten() for i in np.meshgrid(xs,ys,zs)]).T + self.points_old = None + if self.use_tubes: + self.visual = vispy.scene.visuals.Line( + color=(1,1,1,0.75), + connect='segments', + ) + self.view.add(self.visual) + self.view.camera = 'turntable' + self.update() + else: + p = self.points + points = [p] + for i in range(200): + print(i) + curl = generate(p) + p2 = p + curl*0.1*s + p = p2 + points.append(p) + points = np.stack(points, axis=1) / s + # points = (count*count, N, 3) where first dimension chooses which + # trajectory, and second dimension proceeds along time/iterations + # of that trajectory. + for traj in points: + tube = vispy.scene.visuals.Tube(points=traj, radius=0.4*s) + self.view.add(tube) + self.view.camera = 'turntable' + def update(self, ev=None): + if not self.use_tubes: + return + t = 0 if ev is None else ev.elapsed + # Get velocity for current points: + curl = generate(self.points) + # + a1 = self.points + a2 = self.points + curl*0.1*self.s + lines = np.hstack((a1, a2)).reshape(self.points.shape[0]*2, -1) / self.s + self.points_old = np.vstack((self.points_old, lines)) + maxpoints = self.points.shape[0] * 200 + extra = self.points_old.shape[0] - maxpoints + if extra > 0: + self.points_old = self.points_old[extra:] + self.visual.set_data( + pos=self.points_old, + #arrows=np.hstack((a1, a2)), + ) + self.points = a2 + # self.scatter = visuals.Markers() + # self.scatter.set_data(self.d, edge_color=None, face_color=(1, 0.5, 1, .5), size=4) + # view.add(self.scatter) + #m = np.array([[np.cos(t), np.sin(t*1.01), np.cos(t*1.02)]]) + #d2 = self.d*m + #self.scatter.set_data(d2, edge_color=None, face_color=(1, 0.5, 1, .5), size=4) + +def main(): + # + # Make a canvas and add simple view + # + canvas = vispy.scene.SceneCanvas(keys='interactive', show=True) + view = canvas.central_widget.add_view() + import sys + # add a colored 3D axis for orientation + axis = visuals.XYZAxis(parent=view.scene) + timer = vispy.app.Timer() + da = Data(view) + # Problem is how slow update() is: + timer.connect(da.update) + timer.start(0.05) + if sys.flags.interactive != 1: + vispy.app.run() + +if __name__ == '__main__': + main() diff --git a/python/mac_m1_notes.txt b/python/mac_m1_notes.txt new file mode 100644 index 0000000..8c3e6a6 --- /dev/null +++ b/python/mac_m1_notes.txt @@ -0,0 +1,5 @@ +How I coaxed this into working on ARM64 & Apple Silicon: +- Yes, you need Rosetta. +- Install Conda. +- Run: + conda install -c conda-force python=3.9 vispy pytest pyqt diff --git a/python/opensimplex.py b/python/opensimplex.py new file mode 100644 index 0000000..09cbf92 --- /dev/null +++ b/python/opensimplex.py @@ -0,0 +1,1922 @@ +# https://github.com/lmas/opensimplex + +# Based on: https://gist.github.com/KdotJPG/b1270127455a94ac5d19 + +import sys +from ctypes import c_int64 +from math import floor as _floor + +if sys.version_info[0] < 3: + def floor(num): + return int(_floor(num)) +else: + floor = _floor + + +STRETCH_CONSTANT_2D = -0.211324865405187 # (1/Math.sqrt(2+1)-1)/2 +SQUISH_CONSTANT_2D = 0.366025403784439 # (Math.sqrt(2+1)-1)/2 +STRETCH_CONSTANT_3D = -1.0 / 6 # (1/Math.sqrt(3+1)-1)/3 +SQUISH_CONSTANT_3D = 1.0 / 3 # (Math.sqrt(3+1)-1)/3 +STRETCH_CONSTANT_4D = -0.138196601125011 # (1/Math.sqrt(4+1)-1)/4 +SQUISH_CONSTANT_4D = 0.309016994374947 # (Math.sqrt(4+1)-1)/4 + +NORM_CONSTANT_2D = 47 +NORM_CONSTANT_3D = 103 +NORM_CONSTANT_4D = 30 + +DEFAULT_SEED = 0 + + +# Gradients for 2D. They approximate the directions to the +# vertices of an octagon from the center. +GRADIENTS_2D = ( + 5, 2, 2, 5, + -5, 2, -2, 5, + 5, -2, 2, -5, + -5, -2, -2, -5, +) + +# Gradients for 3D. They approximate the directions to the +# vertices of a rhombicuboctahedron from the center, skewed so +# that the triangular and square facets can be inscribed inside +# circles of the same radius. +GRADIENTS_3D = ( + -11, 4, 4, -4, 11, 4, -4, 4, 11, + 11, 4, 4, 4, 11, 4, 4, 4, 11, + -11, -4, 4, -4, -11, 4, -4, -4, 11, + 11, -4, 4, 4, -11, 4, 4, -4, 11, + -11, 4, -4, -4, 11, -4, -4, 4, -11, + 11, 4, -4, 4, 11, -4, 4, 4, -11, + -11, -4, -4, -4, -11, -4, -4, -4, -11, + 11, -4, -4, 4, -11, -4, 4, -4, -11, +) + +# Gradients for 4D. They approximate the directions to the +# vertices of a disprismatotesseractihexadecachoron from the center, +# skewed so that the tetrahedral and cubic facets can be inscribed inside +# spheres of the same radius. +GRADIENTS_4D = ( + 3, 1, 1, 1, 1, 3, 1, 1, 1, 1, 3, 1, 1, 1, 1, 3, + -3, 1, 1, 1, -1, 3, 1, 1, -1, 1, 3, 1, -1, 1, 1, 3, + 3, -1, 1, 1, 1, -3, 1, 1, 1, -1, 3, 1, 1, -1, 1, 3, + -3, -1, 1, 1, -1, -3, 1, 1, -1, -1, 3, 1, -1, -1, 1, 3, + 3, 1, -1, 1, 1, 3, -1, 1, 1, 1, -3, 1, 1, 1, -1, 3, + -3, 1, -1, 1, -1, 3, -1, 1, -1, 1, -3, 1, -1, 1, -1, 3, + 3, -1, -1, 1, 1, -3, -1, 1, 1, -1, -3, 1, 1, -1, -1, 3, + -3, -1, -1, 1, -1, -3, -1, 1, -1, -1, -3, 1, -1, -1, -1, 3, + 3, 1, 1, -1, 1, 3, 1, -1, 1, 1, 3, -1, 1, 1, 1, -3, + -3, 1, 1, -1, -1, 3, 1, -1, -1, 1, 3, -1, -1, 1, 1, -3, + 3, -1, 1, -1, 1, -3, 1, -1, 1, -1, 3, -1, 1, -1, 1, -3, + -3, -1, 1, -1, -1, -3, 1, -1, -1, -1, 3, -1, -1, -1, 1, -3, + 3, 1, -1, -1, 1, 3, -1, -1, 1, 1, -3, -1, 1, 1, -1, -3, + -3, 1, -1, -1, -1, 3, -1, -1, -1, 1, -3, -1, -1, 1, -1, -3, + 3, -1, -1, -1, 1, -3, -1, -1, 1, -1, -3, -1, 1, -1, -1, -3, + -3, -1, -1, -1, -1, -3, -1, -1, -1, -1, -3, -1, -1, -1, -1, -3, +) + + +def overflow(x): + # Since normal python ints and longs can be quite humongous we have to use + # this hack to make them be able to overflow + return c_int64(x).value + + +class OpenSimplex(object): + """ + OpenSimplex n-dimensional gradient noise functions. + """ + + def __init__(self, seed=DEFAULT_SEED): + """ + Initiate the class using a permutation array generated from a 64-bit seed number. + """ + # Generates a proper permutation (i.e. doesn't merely perform N + # successive pair swaps on a base array) + perm = self._perm = [0] * 256 # Have to zero fill so we can properly loop over it later + perm_grad_index_3D = self._perm_grad_index_3D = [0] * 256 + source = [i for i in range(0, 256)] + seed = overflow(seed * 6364136223846793005 + 1442695040888963407) + seed = overflow(seed * 6364136223846793005 + 1442695040888963407) + seed = overflow(seed * 6364136223846793005 + 1442695040888963407) + for i in range(255, -1, -1): + seed = overflow(seed * 6364136223846793005 + 1442695040888963407) + r = int((seed + 31) % (i + 1)) + if r < 0: + r += i + 1 + perm[i] = source[r] + perm_grad_index_3D[i] = int((perm[i] % (len(GRADIENTS_3D) / 3)) * 3) + source[r] = source[i] + + def _extrapolate2d(self, xsb, ysb, dx, dy): + perm = self._perm + index = perm[(perm[xsb & 0xFF] + ysb) & 0xFF] & 0x0E + + g1, g2 = GRADIENTS_2D[index:index + 2] + return g1 * dx + g2 * dy + + def _extrapolate3d(self, xsb, ysb, zsb, dx, dy, dz): + perm = self._perm + index = self._perm_grad_index_3D[ + (perm[(perm[xsb & 0xFF] + ysb) & 0xFF] + zsb) & 0xFF + ] + + g1, g2, g3 = GRADIENTS_3D[index:index + 3] + return g1 * dx + g2 * dy + g3 * dz + + def _extrapolate4d(self, xsb, ysb, zsb, wsb, dx, dy, dz, dw): + perm = self._perm + index = perm[( + perm[( + perm[(perm[xsb & 0xFF] + ysb) & 0xFF] + zsb + ) & 0xFF] + wsb + ) & 0xFF] & 0xFC + + g1, g2, g3, g4 = GRADIENTS_4D[index:index + 4] + return g1 * dx + g2 * dy + g3 * dz + g4 * dw + + def noise2d(self, x, y): + """ + Generate 2D OpenSimplex noise from X,Y coordinates. + """ + # Place input coordinates onto grid. + stretch_offset = (x + y) * STRETCH_CONSTANT_2D + xs = x + stretch_offset + ys = y + stretch_offset + + # Floor to get grid coordinates of rhombus (stretched square) super-cell origin. + xsb = floor(xs) + ysb = floor(ys) + + # Skew out to get actual coordinates of rhombus origin. We'll need these later. + squish_offset = (xsb + ysb) * SQUISH_CONSTANT_2D + xb = xsb + squish_offset + yb = ysb + squish_offset + + # Compute grid coordinates relative to rhombus origin. + xins = xs - xsb + yins = ys - ysb + + # Sum those together to get a value that determines which region we're in. + in_sum = xins + yins + + # Positions relative to origin point. + dx0 = x - xb + dy0 = y - yb + + value = 0 + + # Contribution (1,0) + dx1 = dx0 - 1 - SQUISH_CONSTANT_2D + dy1 = dy0 - 0 - SQUISH_CONSTANT_2D + attn1 = 2 - dx1 * dx1 - dy1 * dy1 + extrapolate = self._extrapolate2d + if attn1 > 0: + attn1 *= attn1 + value += attn1 * attn1 * extrapolate(xsb + 1, ysb + 0, dx1, dy1) + + # Contribution (0,1) + dx2 = dx0 - 0 - SQUISH_CONSTANT_2D + dy2 = dy0 - 1 - SQUISH_CONSTANT_2D + attn2 = 2 - dx2 * dx2 - dy2 * dy2 + if attn2 > 0: + attn2 *= attn2 + value += attn2 * attn2 * extrapolate(xsb + 0, ysb + 1, dx2, dy2) + + if in_sum <= 1: # We're inside the triangle (2-Simplex) at (0,0) + zins = 1 - in_sum + if zins > xins or zins > yins: # (0,0) is one of the closest two triangular vertices + if xins > yins: + xsv_ext = xsb + 1 + ysv_ext = ysb - 1 + dx_ext = dx0 - 1 + dy_ext = dy0 + 1 + else: + xsv_ext = xsb - 1 + ysv_ext = ysb + 1 + dx_ext = dx0 + 1 + dy_ext = dy0 - 1 + else: # (1,0) and (0,1) are the closest two vertices. + xsv_ext = xsb + 1 + ysv_ext = ysb + 1 + dx_ext = dx0 - 1 - 2 * SQUISH_CONSTANT_2D + dy_ext = dy0 - 1 - 2 * SQUISH_CONSTANT_2D + else: # We're inside the triangle (2-Simplex) at (1,1) + zins = 2 - in_sum + if zins < xins or zins < yins: # (0,0) is one of the closest two triangular vertices + if xins > yins: + xsv_ext = xsb + 2 + ysv_ext = ysb + 0 + dx_ext = dx0 - 2 - 2 * SQUISH_CONSTANT_2D + dy_ext = dy0 + 0 - 2 * SQUISH_CONSTANT_2D + else: + xsv_ext = xsb + 0 + ysv_ext = ysb + 2 + dx_ext = dx0 + 0 - 2 * SQUISH_CONSTANT_2D + dy_ext = dy0 - 2 - 2 * SQUISH_CONSTANT_2D + else: # (1,0) and (0,1) are the closest two vertices. + dx_ext = dx0 + dy_ext = dy0 + xsv_ext = xsb + ysv_ext = ysb + xsb += 1 + ysb += 1 + dx0 = dx0 - 1 - 2 * SQUISH_CONSTANT_2D + dy0 = dy0 - 1 - 2 * SQUISH_CONSTANT_2D + + # Contribution (0,0) or (1,1) + attn0 = 2 - dx0 * dx0 - dy0 * dy0 + if attn0 > 0: + attn0 *= attn0 + value += attn0 * attn0 * extrapolate(xsb, ysb, dx0, dy0) + + # Extra Vertex + attn_ext = 2 - dx_ext * dx_ext - dy_ext * dy_ext + if attn_ext > 0: + attn_ext *= attn_ext + value += attn_ext * attn_ext * extrapolate(xsv_ext, ysv_ext, dx_ext, dy_ext) + + return value / NORM_CONSTANT_2D + + def noise3d(self, x, y, z): + """ + Generate 3D OpenSimplex noise from X,Y,Z coordinates. + """ + # Place input coordinates on simplectic honeycomb. + stretch_offset = (x + y + z) * STRETCH_CONSTANT_3D + xs = x + stretch_offset + ys = y + stretch_offset + zs = z + stretch_offset + + # Floor to get simplectic honeycomb coordinates of rhombohedron (stretched cube) super-cell origin. + xsb = floor(xs) + ysb = floor(ys) + zsb = floor(zs) + + # Skew out to get actual coordinates of rhombohedron origin. We'll need these later. + squish_offset = (xsb + ysb + zsb) * SQUISH_CONSTANT_3D + xb = xsb + squish_offset + yb = ysb + squish_offset + zb = zsb + squish_offset + + # Compute simplectic honeycomb coordinates relative to rhombohedral origin. + xins = xs - xsb + yins = ys - ysb + zins = zs - zsb + + # Sum those together to get a value that determines which region we're in. + in_sum = xins + yins + zins + + # Positions relative to origin point. + dx0 = x - xb + dy0 = y - yb + dz0 = z - zb + + value = 0 + extrapolate = self._extrapolate3d + if in_sum <= 1: # We're inside the tetrahedron (3-Simplex) at (0,0,0) + + # Determine which two of (0,0,1), (0,1,0), (1,0,0) are closest. + a_point = 0x01 + a_score = xins + b_point = 0x02 + b_score = yins + if a_score >= b_score and zins > b_score: + b_score = zins + b_point = 0x04 + elif a_score < b_score and zins > a_score: + a_score = zins + a_point = 0x04 + + # Now we determine the two lattice points not part of the tetrahedron that may contribute. + # This depends on the closest two tetrahedral vertices, including (0,0,0) + wins = 1 - in_sum + if wins > a_score or wins > b_score: # (0,0,0) is one of the closest two tetrahedral vertices. + c = b_point if (b_score > a_score) else a_point # Our other closest vertex is the closest out of a and b. + + if (c & 0x01) == 0: + xsv_ext0 = xsb - 1 + xsv_ext1 = xsb + dx_ext0 = dx0 + 1 + dx_ext1 = dx0 + else: + xsv_ext0 = xsv_ext1 = xsb + 1 + dx_ext0 = dx_ext1 = dx0 - 1 + + if (c & 0x02) == 0: + ysv_ext0 = ysv_ext1 = ysb + dy_ext0 = dy_ext1 = dy0 + if (c & 0x01) == 0: + ysv_ext1 -= 1 + dy_ext1 += 1 + else: + ysv_ext0 -= 1 + dy_ext0 += 1 + else: + ysv_ext0 = ysv_ext1 = ysb + 1 + dy_ext0 = dy_ext1 = dy0 - 1 + + if (c & 0x04) == 0: + zsv_ext0 = zsb + zsv_ext1 = zsb - 1 + dz_ext0 = dz0 + dz_ext1 = dz0 + 1 + else: + zsv_ext0 = zsv_ext1 = zsb + 1 + dz_ext0 = dz_ext1 = dz0 - 1 + else: # (0,0,0) is not one of the closest two tetrahedral vertices. + c = (a_point | b_point) # Our two extra vertices are determined by the closest two. + + if (c & 0x01) == 0: + xsv_ext0 = xsb + xsv_ext1 = xsb - 1 + dx_ext0 = dx0 - 2 * SQUISH_CONSTANT_3D + dx_ext1 = dx0 + 1 - SQUISH_CONSTANT_3D + else: + xsv_ext0 = xsv_ext1 = xsb + 1 + dx_ext0 = dx0 - 1 - 2 * SQUISH_CONSTANT_3D + dx_ext1 = dx0 - 1 - SQUISH_CONSTANT_3D + + if (c & 0x02) == 0: + ysv_ext0 = ysb + ysv_ext1 = ysb - 1 + dy_ext0 = dy0 - 2 * SQUISH_CONSTANT_3D + dy_ext1 = dy0 + 1 - SQUISH_CONSTANT_3D + else: + ysv_ext0 = ysv_ext1 = ysb + 1 + dy_ext0 = dy0 - 1 - 2 * SQUISH_CONSTANT_3D + dy_ext1 = dy0 - 1 - SQUISH_CONSTANT_3D + + if (c & 0x04) == 0: + zsv_ext0 = zsb + zsv_ext1 = zsb - 1 + dz_ext0 = dz0 - 2 * SQUISH_CONSTANT_3D + dz_ext1 = dz0 + 1 - SQUISH_CONSTANT_3D + else: + zsv_ext0 = zsv_ext1 = zsb + 1 + dz_ext0 = dz0 - 1 - 2 * SQUISH_CONSTANT_3D + dz_ext1 = dz0 - 1 - SQUISH_CONSTANT_3D + + # Contribution (0,0,0) + attn0 = 2 - dx0 * dx0 - dy0 * dy0 - dz0 * dz0 + if attn0 > 0: + attn0 *= attn0 + value += attn0 * attn0 * extrapolate(xsb + 0, ysb + 0, zsb + 0, dx0, dy0, dz0) + + # Contribution (1,0,0) + dx1 = dx0 - 1 - SQUISH_CONSTANT_3D + dy1 = dy0 - 0 - SQUISH_CONSTANT_3D + dz1 = dz0 - 0 - SQUISH_CONSTANT_3D + attn1 = 2 - dx1 * dx1 - dy1 * dy1 - dz1 * dz1 + if attn1 > 0: + attn1 *= attn1 + value += attn1 * attn1 * extrapolate(xsb + 1, ysb + 0, zsb + 0, dx1, dy1, dz1) + + # Contribution (0,1,0) + dx2 = dx0 - 0 - SQUISH_CONSTANT_3D + dy2 = dy0 - 1 - SQUISH_CONSTANT_3D + dz2 = dz1 + attn2 = 2 - dx2 * dx2 - dy2 * dy2 - dz2 * dz2 + if attn2 > 0: + attn2 *= attn2 + value += attn2 * attn2 * extrapolate(xsb + 0, ysb + 1, zsb + 0, dx2, dy2, dz2) + + # Contribution (0,0,1) + dx3 = dx2 + dy3 = dy1 + dz3 = dz0 - 1 - SQUISH_CONSTANT_3D + attn3 = 2 - dx3 * dx3 - dy3 * dy3 - dz3 * dz3 + if attn3 > 0: + attn3 *= attn3 + value += attn3 * attn3 * extrapolate(xsb + 0, ysb + 0, zsb + 1, dx3, dy3, dz3) + elif in_sum >= 2: # We're inside the tetrahedron (3-Simplex) at (1,1,1) + + # Determine which two tetrahedral vertices are the closest, out of (1,1,0), (1,0,1), (0,1,1) but not (1,1,1). + a_point = 0x06 + a_score = xins + b_point = 0x05 + b_score = yins + if a_score <= b_score and zins < b_score: + b_score = zins + b_point = 0x03 + elif a_score > b_score and zins < a_score: + a_score = zins + a_point = 0x03 + + # Now we determine the two lattice points not part of the tetrahedron that may contribute. + # This depends on the closest two tetrahedral vertices, including (1,1,1) + wins = 3 - in_sum + if wins < a_score or wins < b_score: # (1,1,1) is one of the closest two tetrahedral vertices. + c = b_point if (b_score < a_score) else a_point # Our other closest vertex is the closest out of a and b. + + if (c & 0x01) != 0: + xsv_ext0 = xsb + 2 + xsv_ext1 = xsb + 1 + dx_ext0 = dx0 - 2 - 3 * SQUISH_CONSTANT_3D + dx_ext1 = dx0 - 1 - 3 * SQUISH_CONSTANT_3D + else: + xsv_ext0 = xsv_ext1 = xsb + dx_ext0 = dx_ext1 = dx0 - 3 * SQUISH_CONSTANT_3D + + if (c & 0x02) != 0: + ysv_ext0 = ysv_ext1 = ysb + 1 + dy_ext0 = dy_ext1 = dy0 - 1 - 3 * SQUISH_CONSTANT_3D + if (c & 0x01) != 0: + ysv_ext1 += 1 + dy_ext1 -= 1 + else: + ysv_ext0 += 1 + dy_ext0 -= 1 + else: + ysv_ext0 = ysv_ext1 = ysb + dy_ext0 = dy_ext1 = dy0 - 3 * SQUISH_CONSTANT_3D + + if (c & 0x04) != 0: + zsv_ext0 = zsb + 1 + zsv_ext1 = zsb + 2 + dz_ext0 = dz0 - 1 - 3 * SQUISH_CONSTANT_3D + dz_ext1 = dz0 - 2 - 3 * SQUISH_CONSTANT_3D + else: + zsv_ext0 = zsv_ext1 = zsb + dz_ext0 = dz_ext1 = dz0 - 3 * SQUISH_CONSTANT_3D + else: # (1,1,1) is not one of the closest two tetrahedral vertices. + c = (a_point & b_point) # Our two extra vertices are determined by the closest two. + + if (c & 0x01) != 0: + xsv_ext0 = xsb + 1 + xsv_ext1 = xsb + 2 + dx_ext0 = dx0 - 1 - SQUISH_CONSTANT_3D + dx_ext1 = dx0 - 2 - 2 * SQUISH_CONSTANT_3D + else: + xsv_ext0 = xsv_ext1 = xsb + dx_ext0 = dx0 - SQUISH_CONSTANT_3D + dx_ext1 = dx0 - 2 * SQUISH_CONSTANT_3D + + if (c & 0x02) != 0: + ysv_ext0 = ysb + 1 + ysv_ext1 = ysb + 2 + dy_ext0 = dy0 - 1 - SQUISH_CONSTANT_3D + dy_ext1 = dy0 - 2 - 2 * SQUISH_CONSTANT_3D + else: + ysv_ext0 = ysv_ext1 = ysb + dy_ext0 = dy0 - SQUISH_CONSTANT_3D + dy_ext1 = dy0 - 2 * SQUISH_CONSTANT_3D + + if (c & 0x04) != 0: + zsv_ext0 = zsb + 1 + zsv_ext1 = zsb + 2 + dz_ext0 = dz0 - 1 - SQUISH_CONSTANT_3D + dz_ext1 = dz0 - 2 - 2 * SQUISH_CONSTANT_3D + else: + zsv_ext0 = zsv_ext1 = zsb + dz_ext0 = dz0 - SQUISH_CONSTANT_3D + dz_ext1 = dz0 - 2 * SQUISH_CONSTANT_3D + + # Contribution (1,1,0) + dx3 = dx0 - 1 - 2 * SQUISH_CONSTANT_3D + dy3 = dy0 - 1 - 2 * SQUISH_CONSTANT_3D + dz3 = dz0 - 0 - 2 * SQUISH_CONSTANT_3D + attn3 = 2 - dx3 * dx3 - dy3 * dy3 - dz3 * dz3 + if attn3 > 0: + attn3 *= attn3 + value += attn3 * attn3 * extrapolate(xsb + 1, ysb + 1, zsb + 0, dx3, dy3, dz3) + + # Contribution (1,0,1) + dx2 = dx3 + dy2 = dy0 - 0 - 2 * SQUISH_CONSTANT_3D + dz2 = dz0 - 1 - 2 * SQUISH_CONSTANT_3D + attn2 = 2 - dx2 * dx2 - dy2 * dy2 - dz2 * dz2 + if attn2 > 0: + attn2 *= attn2 + value += attn2 * attn2 * extrapolate(xsb + 1, ysb + 0, zsb + 1, dx2, dy2, dz2) + + # Contribution (0,1,1) + dx1 = dx0 - 0 - 2 * SQUISH_CONSTANT_3D + dy1 = dy3 + dz1 = dz2 + attn1 = 2 - dx1 * dx1 - dy1 * dy1 - dz1 * dz1 + if attn1 > 0: + attn1 *= attn1 + value += attn1 * attn1 * extrapolate(xsb + 0, ysb + 1, zsb + 1, dx1, dy1, dz1) + + # Contribution (1,1,1) + dx0 = dx0 - 1 - 3 * SQUISH_CONSTANT_3D + dy0 = dy0 - 1 - 3 * SQUISH_CONSTANT_3D + dz0 = dz0 - 1 - 3 * SQUISH_CONSTANT_3D + attn0 = 2 - dx0 * dx0 - dy0 * dy0 - dz0 * dz0 + if attn0 > 0: + attn0 *= attn0 + value += attn0 * attn0 * extrapolate(xsb + 1, ysb + 1, zsb + 1, dx0, dy0, dz0) + else: # We're inside the octahedron (Rectified 3-Simplex) in between. + # Decide between point (0,0,1) and (1,1,0) as closest + p1 = xins + yins + if p1 > 1: + a_score = p1 - 1 + a_point = 0x03 + a_is_further_side = True + else: + a_score = 1 - p1 + a_point = 0x04 + a_is_further_side = False + + # Decide between point (0,1,0) and (1,0,1) as closest + p2 = xins + zins + if p2 > 1: + b_score = p2 - 1 + b_point = 0x05 + b_is_further_side = True + else: + b_score = 1 - p2 + b_point = 0x02 + b_is_further_side = False + + # The closest out of the two (1,0,0) and (0,1,1) will replace the furthest out of the two decided above, if closer. + p3 = yins + zins + if p3 > 1: + score = p3 - 1 + if a_score <= b_score and a_score < score: + a_point = 0x06 + a_is_further_side = True + elif a_score > b_score and b_score < score: + b_point = 0x06 + b_is_further_side = True + else: + score = 1 - p3 + if a_score <= b_score and a_score < score: + a_point = 0x01 + a_is_further_side = False + elif a_score > b_score and b_score < score: + b_point = 0x01 + b_is_further_side = False + + # Where each of the two closest points are determines how the extra two vertices are calculated. + if a_is_further_side == b_is_further_side: + if a_is_further_side: # Both closest points on (1,1,1) side + + # One of the two extra points is (1,1,1) + dx_ext0 = dx0 - 1 - 3 * SQUISH_CONSTANT_3D + dy_ext0 = dy0 - 1 - 3 * SQUISH_CONSTANT_3D + dz_ext0 = dz0 - 1 - 3 * SQUISH_CONSTANT_3D + xsv_ext0 = xsb + 1 + ysv_ext0 = ysb + 1 + zsv_ext0 = zsb + 1 + + # Other extra point is based on the shared axis. + c = (a_point & b_point) + if (c & 0x01) != 0: + dx_ext1 = dx0 - 2 - 2 * SQUISH_CONSTANT_3D + dy_ext1 = dy0 - 2 * SQUISH_CONSTANT_3D + dz_ext1 = dz0 - 2 * SQUISH_CONSTANT_3D + xsv_ext1 = xsb + 2 + ysv_ext1 = ysb + zsv_ext1 = zsb + elif (c & 0x02) != 0: + dx_ext1 = dx0 - 2 * SQUISH_CONSTANT_3D + dy_ext1 = dy0 - 2 - 2 * SQUISH_CONSTANT_3D + dz_ext1 = dz0 - 2 * SQUISH_CONSTANT_3D + xsv_ext1 = xsb + ysv_ext1 = ysb + 2 + zsv_ext1 = zsb + else: + dx_ext1 = dx0 - 2 * SQUISH_CONSTANT_3D + dy_ext1 = dy0 - 2 * SQUISH_CONSTANT_3D + dz_ext1 = dz0 - 2 - 2 * SQUISH_CONSTANT_3D + xsv_ext1 = xsb + ysv_ext1 = ysb + zsv_ext1 = zsb + 2 + else:# Both closest points on (0,0,0) side + + # One of the two extra points is (0,0,0) + dx_ext0 = dx0 + dy_ext0 = dy0 + dz_ext0 = dz0 + xsv_ext0 = xsb + ysv_ext0 = ysb + zsv_ext0 = zsb + + # Other extra point is based on the omitted axis. + c = (a_point | b_point) + if (c & 0x01) == 0: + dx_ext1 = dx0 + 1 - SQUISH_CONSTANT_3D + dy_ext1 = dy0 - 1 - SQUISH_CONSTANT_3D + dz_ext1 = dz0 - 1 - SQUISH_CONSTANT_3D + xsv_ext1 = xsb - 1 + ysv_ext1 = ysb + 1 + zsv_ext1 = zsb + 1 + elif (c & 0x02) == 0: + dx_ext1 = dx0 - 1 - SQUISH_CONSTANT_3D + dy_ext1 = dy0 + 1 - SQUISH_CONSTANT_3D + dz_ext1 = dz0 - 1 - SQUISH_CONSTANT_3D + xsv_ext1 = xsb + 1 + ysv_ext1 = ysb - 1 + zsv_ext1 = zsb + 1 + else: + dx_ext1 = dx0 - 1 - SQUISH_CONSTANT_3D + dy_ext1 = dy0 - 1 - SQUISH_CONSTANT_3D + dz_ext1 = dz0 + 1 - SQUISH_CONSTANT_3D + xsv_ext1 = xsb + 1 + ysv_ext1 = ysb + 1 + zsv_ext1 = zsb - 1 + else: # One point on (0,0,0) side, one point on (1,1,1) side + if a_is_further_side: + c1 = a_point + c2 = b_point + else: + c1 = b_point + c2 = a_point + + # One contribution is a _permutation of (1,1,-1) + if (c1 & 0x01) == 0: + dx_ext0 = dx0 + 1 - SQUISH_CONSTANT_3D + dy_ext0 = dy0 - 1 - SQUISH_CONSTANT_3D + dz_ext0 = dz0 - 1 - SQUISH_CONSTANT_3D + xsv_ext0 = xsb - 1 + ysv_ext0 = ysb + 1 + zsv_ext0 = zsb + 1 + elif (c1 & 0x02) == 0: + dx_ext0 = dx0 - 1 - SQUISH_CONSTANT_3D + dy_ext0 = dy0 + 1 - SQUISH_CONSTANT_3D + dz_ext0 = dz0 - 1 - SQUISH_CONSTANT_3D + xsv_ext0 = xsb + 1 + ysv_ext0 = ysb - 1 + zsv_ext0 = zsb + 1 + else: + dx_ext0 = dx0 - 1 - SQUISH_CONSTANT_3D + dy_ext0 = dy0 - 1 - SQUISH_CONSTANT_3D + dz_ext0 = dz0 + 1 - SQUISH_CONSTANT_3D + xsv_ext0 = xsb + 1 + ysv_ext0 = ysb + 1 + zsv_ext0 = zsb - 1 + + # One contribution is a _permutation of (0,0,2) + dx_ext1 = dx0 - 2 * SQUISH_CONSTANT_3D + dy_ext1 = dy0 - 2 * SQUISH_CONSTANT_3D + dz_ext1 = dz0 - 2 * SQUISH_CONSTANT_3D + xsv_ext1 = xsb + ysv_ext1 = ysb + zsv_ext1 = zsb + if (c2 & 0x01) != 0: + dx_ext1 -= 2 + xsv_ext1 += 2 + elif (c2 & 0x02) != 0: + dy_ext1 -= 2 + ysv_ext1 += 2 + else: + dz_ext1 -= 2 + zsv_ext1 += 2 + + # Contribution (1,0,0) + dx1 = dx0 - 1 - SQUISH_CONSTANT_3D + dy1 = dy0 - 0 - SQUISH_CONSTANT_3D + dz1 = dz0 - 0 - SQUISH_CONSTANT_3D + attn1 = 2 - dx1 * dx1 - dy1 * dy1 - dz1 * dz1 + if attn1 > 0: + attn1 *= attn1 + value += attn1 * attn1 * extrapolate(xsb + 1, ysb + 0, zsb + 0, dx1, dy1, dz1) + + # Contribution (0,1,0) + dx2 = dx0 - 0 - SQUISH_CONSTANT_3D + dy2 = dy0 - 1 - SQUISH_CONSTANT_3D + dz2 = dz1 + attn2 = 2 - dx2 * dx2 - dy2 * dy2 - dz2 * dz2 + if attn2 > 0: + attn2 *= attn2 + value += attn2 * attn2 * extrapolate(xsb + 0, ysb + 1, zsb + 0, dx2, dy2, dz2) + + # Contribution (0,0,1) + dx3 = dx2 + dy3 = dy1 + dz3 = dz0 - 1 - SQUISH_CONSTANT_3D + attn3 = 2 - dx3 * dx3 - dy3 * dy3 - dz3 * dz3 + if attn3 > 0: + attn3 *= attn3 + value += attn3 * attn3 * extrapolate(xsb + 0, ysb + 0, zsb + 1, dx3, dy3, dz3) + + # Contribution (1,1,0) + dx4 = dx0 - 1 - 2 * SQUISH_CONSTANT_3D + dy4 = dy0 - 1 - 2 * SQUISH_CONSTANT_3D + dz4 = dz0 - 0 - 2 * SQUISH_CONSTANT_3D + attn4 = 2 - dx4 * dx4 - dy4 * dy4 - dz4 * dz4 + if attn4 > 0: + attn4 *= attn4 + value += attn4 * attn4 * extrapolate(xsb + 1, ysb + 1, zsb + 0, dx4, dy4, dz4) + + # Contribution (1,0,1) + dx5 = dx4 + dy5 = dy0 - 0 - 2 * SQUISH_CONSTANT_3D + dz5 = dz0 - 1 - 2 * SQUISH_CONSTANT_3D + attn5 = 2 - dx5 * dx5 - dy5 * dy5 - dz5 * dz5 + if attn5 > 0: + attn5 *= attn5 + value += attn5 * attn5 * extrapolate(xsb + 1, ysb + 0, zsb + 1, dx5, dy5, dz5) + + # Contribution (0,1,1) + dx6 = dx0 - 0 - 2 * SQUISH_CONSTANT_3D + dy6 = dy4 + dz6 = dz5 + attn6 = 2 - dx6 * dx6 - dy6 * dy6 - dz6 * dz6 + if attn6 > 0: + attn6 *= attn6 + value += attn6 * attn6 * extrapolate(xsb + 0, ysb + 1, zsb + 1, dx6, dy6, dz6) + + # First extra vertex + attn_ext0 = 2 - dx_ext0 * dx_ext0 - dy_ext0 * dy_ext0 - dz_ext0 * dz_ext0 + if attn_ext0 > 0: + attn_ext0 *= attn_ext0 + value += attn_ext0 * attn_ext0 * extrapolate(xsv_ext0, ysv_ext0, zsv_ext0, dx_ext0, dy_ext0, dz_ext0) + + # Second extra vertex + attn_ext1 = 2 - dx_ext1 * dx_ext1 - dy_ext1 * dy_ext1 - dz_ext1 * dz_ext1 + if attn_ext1 > 0: + attn_ext1 *= attn_ext1 + value += attn_ext1 * attn_ext1 * extrapolate(xsv_ext1, ysv_ext1, zsv_ext1, dx_ext1, dy_ext1, dz_ext1) + + return value / NORM_CONSTANT_3D + + def noise4d(self, x, y, z, w): + """ + Generate 4D OpenSimplex noise from X,Y,Z,W coordinates. + """ + # Place input coordinates on simplectic honeycomb. + stretch_offset = (x + y + z + w) * STRETCH_CONSTANT_4D + xs = x + stretch_offset + ys = y + stretch_offset + zs = z + stretch_offset + ws = w + stretch_offset + + # Floor to get simplectic honeycomb coordinates of rhombo-hypercube super-cell origin. + xsb = floor(xs) + ysb = floor(ys) + zsb = floor(zs) + wsb = floor(ws) + + # Skew out to get actual coordinates of stretched rhombo-hypercube origin. We'll need these later. + squish_offset = (xsb + ysb + zsb + wsb) * SQUISH_CONSTANT_4D + xb = xsb + squish_offset + yb = ysb + squish_offset + zb = zsb + squish_offset + wb = wsb + squish_offset + + # Compute simplectic honeycomb coordinates relative to rhombo-hypercube origin. + xins = xs - xsb + yins = ys - ysb + zins = zs - zsb + wins = ws - wsb + + # Sum those together to get a value that determines which region we're in. + in_sum = xins + yins + zins + wins + + # Positions relative to origin po. + dx0 = x - xb + dy0 = y - yb + dz0 = z - zb + dw0 = w - wb + + value = 0 + extrapolate = self._extrapolate4d + if in_sum <= 1: # We're inside the pentachoron (4-Simplex) at (0,0,0,0) + + # Determine which two of (0,0,0,1), (0,0,1,0), (0,1,0,0), (1,0,0,0) are closest. + a_po = 0x01 + a_score = xins + b_po = 0x02 + b_score = yins + if a_score >= b_score and zins > b_score: + b_score = zins + b_po = 0x04 + elif a_score < b_score and zins > a_score: + a_score = zins + a_po = 0x04 + + if a_score >= b_score and wins > b_score: + b_score = wins + b_po = 0x08 + elif a_score < b_score and wins > a_score: + a_score = wins + a_po = 0x08 + + # Now we determine the three lattice pos not part of the pentachoron that may contribute. + # This depends on the closest two pentachoron vertices, including (0,0,0,0) + uins = 1 - in_sum + if uins > a_score or uins > b_score: # (0,0,0,0) is one of the closest two pentachoron vertices. + c = b_po if (b_score > a_score) else a_po # Our other closest vertex is the closest out of a and b. + if (c & 0x01) == 0: + xsv_ext0 = xsb - 1 + xsv_ext1 = xsv_ext2 = xsb + dx_ext0 = dx0 + 1 + dx_ext1 = dx_ext2 = dx0 + else: + xsv_ext0 = xsv_ext1 = xsv_ext2 = xsb + 1 + dx_ext0 = dx_ext1 = dx_ext2 = dx0 - 1 + + if (c & 0x02) == 0: + ysv_ext0 = ysv_ext1 = ysv_ext2 = ysb + dy_ext0 = dy_ext1 = dy_ext2 = dy0 + if (c & 0x01) == 0x01: + ysv_ext0 -= 1 + dy_ext0 += 1 + else: + ysv_ext1 -= 1 + dy_ext1 += 1 + + else: + ysv_ext0 = ysv_ext1 = ysv_ext2 = ysb + 1 + dy_ext0 = dy_ext1 = dy_ext2 = dy0 - 1 + + if (c & 0x04) == 0: + zsv_ext0 = zsv_ext1 = zsv_ext2 = zsb + dz_ext0 = dz_ext1 = dz_ext2 = dz0 + if (c & 0x03) != 0: + if (c & 0x03) == 0x03: + zsv_ext0 -= 1 + dz_ext0 += 1 + else: + zsv_ext1 -= 1 + dz_ext1 += 1 + + else: + zsv_ext2 -= 1 + dz_ext2 += 1 + + else: + zsv_ext0 = zsv_ext1 = zsv_ext2 = zsb + 1 + dz_ext0 = dz_ext1 = dz_ext2 = dz0 - 1 + + if (c & 0x08) == 0: + wsv_ext0 = wsv_ext1 = wsb + wsv_ext2 = wsb - 1 + dw_ext0 = dw_ext1 = dw0 + dw_ext2 = dw0 + 1 + else: + wsv_ext0 = wsv_ext1 = wsv_ext2 = wsb + 1 + dw_ext0 = dw_ext1 = dw_ext2 = dw0 - 1 + + else: # (0,0,0,0) is not one of the closest two pentachoron vertices. + c = (a_po | b_po) # Our three extra vertices are determined by the closest two. + + if (c & 0x01) == 0: + xsv_ext0 = xsv_ext2 = xsb + xsv_ext1 = xsb - 1 + dx_ext0 = dx0 - 2 * SQUISH_CONSTANT_4D + dx_ext1 = dx0 + 1 - SQUISH_CONSTANT_4D + dx_ext2 = dx0 - SQUISH_CONSTANT_4D + else: + xsv_ext0 = xsv_ext1 = xsv_ext2 = xsb + 1 + dx_ext0 = dx0 - 1 - 2 * SQUISH_CONSTANT_4D + dx_ext1 = dx_ext2 = dx0 - 1 - SQUISH_CONSTANT_4D + + if (c & 0x02) == 0: + ysv_ext0 = ysv_ext1 = ysv_ext2 = ysb + dy_ext0 = dy0 - 2 * SQUISH_CONSTANT_4D + dy_ext1 = dy_ext2 = dy0 - SQUISH_CONSTANT_4D + if (c & 0x01) == 0x01: + ysv_ext1 -= 1 + dy_ext1 += 1 + else: + ysv_ext2 -= 1 + dy_ext2 += 1 + + else: + ysv_ext0 = ysv_ext1 = ysv_ext2 = ysb + 1 + dy_ext0 = dy0 - 1 - 2 * SQUISH_CONSTANT_4D + dy_ext1 = dy_ext2 = dy0 - 1 - SQUISH_CONSTANT_4D + + if (c & 0x04) == 0: + zsv_ext0 = zsv_ext1 = zsv_ext2 = zsb + dz_ext0 = dz0 - 2 * SQUISH_CONSTANT_4D + dz_ext1 = dz_ext2 = dz0 - SQUISH_CONSTANT_4D + if (c & 0x03) == 0x03: + zsv_ext1 -= 1 + dz_ext1 += 1 + else: + zsv_ext2 -= 1 + dz_ext2 += 1 + + else: + zsv_ext0 = zsv_ext1 = zsv_ext2 = zsb + 1 + dz_ext0 = dz0 - 1 - 2 * SQUISH_CONSTANT_4D + dz_ext1 = dz_ext2 = dz0 - 1 - SQUISH_CONSTANT_4D + + if (c & 0x08) == 0: + wsv_ext0 = wsv_ext1 = wsb + wsv_ext2 = wsb - 1 + dw_ext0 = dw0 - 2 * SQUISH_CONSTANT_4D + dw_ext1 = dw0 - SQUISH_CONSTANT_4D + dw_ext2 = dw0 + 1 - SQUISH_CONSTANT_4D + else: + wsv_ext0 = wsv_ext1 = wsv_ext2 = wsb + 1 + dw_ext0 = dw0 - 1 - 2 * SQUISH_CONSTANT_4D + dw_ext1 = dw_ext2 = dw0 - 1 - SQUISH_CONSTANT_4D + + # Contribution (0,0,0,0) + attn0 = 2 - dx0 * dx0 - dy0 * dy0 - dz0 * dz0 - dw0 * dw0 + if attn0 > 0: + attn0 *= attn0 + value += attn0 * attn0 * extrapolate(xsb + 0, ysb + 0, zsb + 0, wsb + 0, dx0, dy0, dz0, dw0) + + # Contribution (1,0,0,0) + dx1 = dx0 - 1 - SQUISH_CONSTANT_4D + dy1 = dy0 - 0 - SQUISH_CONSTANT_4D + dz1 = dz0 - 0 - SQUISH_CONSTANT_4D + dw1 = dw0 - 0 - SQUISH_CONSTANT_4D + attn1 = 2 - dx1 * dx1 - dy1 * dy1 - dz1 * dz1 - dw1 * dw1 + if attn1 > 0: + attn1 *= attn1 + value += attn1 * attn1 * extrapolate(xsb + 1, ysb + 0, zsb + 0, wsb + 0, dx1, dy1, dz1, dw1) + + # Contribution (0,1,0,0) + dx2 = dx0 - 0 - SQUISH_CONSTANT_4D + dy2 = dy0 - 1 - SQUISH_CONSTANT_4D + dz2 = dz1 + dw2 = dw1 + attn2 = 2 - dx2 * dx2 - dy2 * dy2 - dz2 * dz2 - dw2 * dw2 + if attn2 > 0: + attn2 *= attn2 + value += attn2 * attn2 * extrapolate(xsb + 0, ysb + 1, zsb + 0, wsb + 0, dx2, dy2, dz2, dw2) + + # Contribution (0,0,1,0) + dx3 = dx2 + dy3 = dy1 + dz3 = dz0 - 1 - SQUISH_CONSTANT_4D + dw3 = dw1 + attn3 = 2 - dx3 * dx3 - dy3 * dy3 - dz3 * dz3 - dw3 * dw3 + if attn3 > 0: + attn3 *= attn3 + value += attn3 * attn3 * extrapolate(xsb + 0, ysb + 0, zsb + 1, wsb + 0, dx3, dy3, dz3, dw3) + + # Contribution (0,0,0,1) + dx4 = dx2 + dy4 = dy1 + dz4 = dz1 + dw4 = dw0 - 1 - SQUISH_CONSTANT_4D + attn4 = 2 - dx4 * dx4 - dy4 * dy4 - dz4 * dz4 - dw4 * dw4 + if attn4 > 0: + attn4 *= attn4 + value += attn4 * attn4 * extrapolate(xsb + 0, ysb + 0, zsb + 0, wsb + 1, dx4, dy4, dz4, dw4) + + elif in_sum >= 3: # We're inside the pentachoron (4-Simplex) at (1,1,1,1) + # Determine which two of (1,1,1,0), (1,1,0,1), (1,0,1,1), (0,1,1,1) are closest. + a_po = 0x0E + a_score = xins + b_po = 0x0D + b_score = yins + if a_score <= b_score and zins < b_score: + b_score = zins + b_po = 0x0B + elif a_score > b_score and zins < a_score: + a_score = zins + a_po = 0x0B + + if a_score <= b_score and wins < b_score: + b_score = wins + b_po = 0x07 + elif a_score > b_score and wins < a_score: + a_score = wins + a_po = 0x07 + + # Now we determine the three lattice pos not part of the pentachoron that may contribute. + # This depends on the closest two pentachoron vertices, including (0,0,0,0) + uins = 4 - in_sum + if uins < a_score or uins < b_score: # (1,1,1,1) is one of the closest two pentachoron vertices. + c = b_po if (b_score < a_score) else a_po # Our other closest vertex is the closest out of a and b. + + if (c & 0x01) != 0: + xsv_ext0 = xsb + 2 + xsv_ext1 = xsv_ext2 = xsb + 1 + dx_ext0 = dx0 - 2 - 4 * SQUISH_CONSTANT_4D + dx_ext1 = dx_ext2 = dx0 - 1 - 4 * SQUISH_CONSTANT_4D + else: + xsv_ext0 = xsv_ext1 = xsv_ext2 = xsb + dx_ext0 = dx_ext1 = dx_ext2 = dx0 - 4 * SQUISH_CONSTANT_4D + + if (c & 0x02) != 0: + ysv_ext0 = ysv_ext1 = ysv_ext2 = ysb + 1 + dy_ext0 = dy_ext1 = dy_ext2 = dy0 - 1 - 4 * SQUISH_CONSTANT_4D + if (c & 0x01) != 0: + ysv_ext1 += 1 + dy_ext1 -= 1 + else: + ysv_ext0 += 1 + dy_ext0 -= 1 + + else: + ysv_ext0 = ysv_ext1 = ysv_ext2 = ysb + dy_ext0 = dy_ext1 = dy_ext2 = dy0 - 4 * SQUISH_CONSTANT_4D + + if (c & 0x04) != 0: + zsv_ext0 = zsv_ext1 = zsv_ext2 = zsb + 1 + dz_ext0 = dz_ext1 = dz_ext2 = dz0 - 1 - 4 * SQUISH_CONSTANT_4D + if (c & 0x03) != 0x03: + if (c & 0x03) == 0: + zsv_ext0 += 1 + dz_ext0 -= 1 + else: + zsv_ext1 += 1 + dz_ext1 -= 1 + + else: + zsv_ext2 += 1 + dz_ext2 -= 1 + + else: + zsv_ext0 = zsv_ext1 = zsv_ext2 = zsb + dz_ext0 = dz_ext1 = dz_ext2 = dz0 - 4 * SQUISH_CONSTANT_4D + + if (c & 0x08) != 0: + wsv_ext0 = wsv_ext1 = wsb + 1 + wsv_ext2 = wsb + 2 + dw_ext0 = dw_ext1 = dw0 - 1 - 4 * SQUISH_CONSTANT_4D + dw_ext2 = dw0 - 2 - 4 * SQUISH_CONSTANT_4D + else: + wsv_ext0 = wsv_ext1 = wsv_ext2 = wsb + dw_ext0 = dw_ext1 = dw_ext2 = dw0 - 4 * SQUISH_CONSTANT_4D + + else: # (1,1,1,1) is not one of the closest two pentachoron vertices. + c = (a_po & b_po) # Our three extra vertices are determined by the closest two. + + if (c & 0x01) != 0: + xsv_ext0 = xsv_ext2 = xsb + 1 + xsv_ext1 = xsb + 2 + dx_ext0 = dx0 - 1 - 2 * SQUISH_CONSTANT_4D + dx_ext1 = dx0 - 2 - 3 * SQUISH_CONSTANT_4D + dx_ext2 = dx0 - 1 - 3 * SQUISH_CONSTANT_4D + else: + xsv_ext0 = xsv_ext1 = xsv_ext2 = xsb + dx_ext0 = dx0 - 2 * SQUISH_CONSTANT_4D + dx_ext1 = dx_ext2 = dx0 - 3 * SQUISH_CONSTANT_4D + + if (c & 0x02) != 0: + ysv_ext0 = ysv_ext1 = ysv_ext2 = ysb + 1 + dy_ext0 = dy0 - 1 - 2 * SQUISH_CONSTANT_4D + dy_ext1 = dy_ext2 = dy0 - 1 - 3 * SQUISH_CONSTANT_4D + if (c & 0x01) != 0: + ysv_ext2 += 1 + dy_ext2 -= 1 + else: + ysv_ext1 += 1 + dy_ext1 -= 1 + + else: + ysv_ext0 = ysv_ext1 = ysv_ext2 = ysb + dy_ext0 = dy0 - 2 * SQUISH_CONSTANT_4D + dy_ext1 = dy_ext2 = dy0 - 3 * SQUISH_CONSTANT_4D + + if (c & 0x04) != 0: + zsv_ext0 = zsv_ext1 = zsv_ext2 = zsb + 1 + dz_ext0 = dz0 - 1 - 2 * SQUISH_CONSTANT_4D + dz_ext1 = dz_ext2 = dz0 - 1 - 3 * SQUISH_CONSTANT_4D + if (c & 0x03) != 0: + zsv_ext2 += 1 + dz_ext2 -= 1 + else: + zsv_ext1 += 1 + dz_ext1 -= 1 + + else: + zsv_ext0 = zsv_ext1 = zsv_ext2 = zsb + dz_ext0 = dz0 - 2 * SQUISH_CONSTANT_4D + dz_ext1 = dz_ext2 = dz0 - 3 * SQUISH_CONSTANT_4D + + if (c & 0x08) != 0: + wsv_ext0 = wsv_ext1 = wsb + 1 + wsv_ext2 = wsb + 2 + dw_ext0 = dw0 - 1 - 2 * SQUISH_CONSTANT_4D + dw_ext1 = dw0 - 1 - 3 * SQUISH_CONSTANT_4D + dw_ext2 = dw0 - 2 - 3 * SQUISH_CONSTANT_4D + else: + wsv_ext0 = wsv_ext1 = wsv_ext2 = wsb + dw_ext0 = dw0 - 2 * SQUISH_CONSTANT_4D + dw_ext1 = dw_ext2 = dw0 - 3 * SQUISH_CONSTANT_4D + + # Contribution (1,1,1,0) + dx4 = dx0 - 1 - 3 * SQUISH_CONSTANT_4D + dy4 = dy0 - 1 - 3 * SQUISH_CONSTANT_4D + dz4 = dz0 - 1 - 3 * SQUISH_CONSTANT_4D + dw4 = dw0 - 3 * SQUISH_CONSTANT_4D + attn4 = 2 - dx4 * dx4 - dy4 * dy4 - dz4 * dz4 - dw4 * dw4 + if attn4 > 0: + attn4 *= attn4 + value += attn4 * attn4 * extrapolate(xsb + 1, ysb + 1, zsb + 1, wsb + 0, dx4, dy4, dz4, dw4) + + # Contribution (1,1,0,1) + dx3 = dx4 + dy3 = dy4 + dz3 = dz0 - 3 * SQUISH_CONSTANT_4D + dw3 = dw0 - 1 - 3 * SQUISH_CONSTANT_4D + attn3 = 2 - dx3 * dx3 - dy3 * dy3 - dz3 * dz3 - dw3 * dw3 + if attn3 > 0: + attn3 *= attn3 + value += attn3 * attn3 * extrapolate(xsb + 1, ysb + 1, zsb + 0, wsb + 1, dx3, dy3, dz3, dw3) + + # Contribution (1,0,1,1) + dx2 = dx4 + dy2 = dy0 - 3 * SQUISH_CONSTANT_4D + dz2 = dz4 + dw2 = dw3 + attn2 = 2 - dx2 * dx2 - dy2 * dy2 - dz2 * dz2 - dw2 * dw2 + if attn2 > 0: + attn2 *= attn2 + value += attn2 * attn2 * extrapolate(xsb + 1, ysb + 0, zsb + 1, wsb + 1, dx2, dy2, dz2, dw2) + + # Contribution (0,1,1,1) + dx1 = dx0 - 3 * SQUISH_CONSTANT_4D + dz1 = dz4 + dy1 = dy4 + dw1 = dw3 + attn1 = 2 - dx1 * dx1 - dy1 * dy1 - dz1 * dz1 - dw1 * dw1 + if attn1 > 0: + attn1 *= attn1 + value += attn1 * attn1 * extrapolate(xsb + 0, ysb + 1, zsb + 1, wsb + 1, dx1, dy1, dz1, dw1) + + # Contribution (1,1,1,1) + dx0 = dx0 - 1 - 4 * SQUISH_CONSTANT_4D + dy0 = dy0 - 1 - 4 * SQUISH_CONSTANT_4D + dz0 = dz0 - 1 - 4 * SQUISH_CONSTANT_4D + dw0 = dw0 - 1 - 4 * SQUISH_CONSTANT_4D + attn0 = 2 - dx0 * dx0 - dy0 * dy0 - dz0 * dz0 - dw0 * dw0 + if attn0 > 0: + attn0 *= attn0 + value += attn0 * attn0 * extrapolate(xsb + 1, ysb + 1, zsb + 1, wsb + 1, dx0, dy0, dz0, dw0) + + elif in_sum <= 2: # We're inside the first dispentachoron (Rectified 4-Simplex) + a_is_bigger_side = True + b_is_bigger_side = True + + # Decide between (1,1,0,0) and (0,0,1,1) + if xins + yins > zins + wins: + a_score = xins + yins + a_po = 0x03 + else: + a_score = zins + wins + a_po = 0x0C + + # Decide between (1,0,1,0) and (0,1,0,1) + if xins + zins > yins + wins: + b_score = xins + zins + b_po = 0x05 + else: + b_score = yins + wins + b_po = 0x0A + + # Closer between (1,0,0,1) and (0,1,1,0) will replace the further of a and b, if closer. + if xins + wins > yins + zins: + score = xins + wins + if a_score >= b_score and score > b_score: + b_score = score + b_po = 0x09 + elif a_score < b_score and score > a_score: + a_score = score + a_po = 0x09 + + else: + score = yins + zins + if a_score >= b_score and score > b_score: + b_score = score + b_po = 0x06 + elif a_score < b_score and score > a_score: + a_score = score + a_po = 0x06 + + # Decide if (1,0,0,0) is closer. + p1 = 2 - in_sum + xins + if a_score >= b_score and p1 > b_score: + b_score = p1 + b_po = 0x01 + b_is_bigger_side = False + elif a_score < b_score and p1 > a_score: + a_score = p1 + a_po = 0x01 + a_is_bigger_side = False + + # Decide if (0,1,0,0) is closer. + p2 = 2 - in_sum + yins + if a_score >= b_score and p2 > b_score: + b_score = p2 + b_po = 0x02 + b_is_bigger_side = False + elif a_score < b_score and p2 > a_score: + a_score = p2 + a_po = 0x02 + a_is_bigger_side = False + + # Decide if (0,0,1,0) is closer. + p3 = 2 - in_sum + zins + if a_score >= b_score and p3 > b_score: + b_score = p3 + b_po = 0x04 + b_is_bigger_side = False + elif a_score < b_score and p3 > a_score: + a_score = p3 + a_po = 0x04 + a_is_bigger_side = False + + # Decide if (0,0,0,1) is closer. + p4 = 2 - in_sum + wins + if a_score >= b_score and p4 > b_score: + b_po = 0x08 + b_is_bigger_side = False + elif a_score < b_score and p4 > a_score: + a_po = 0x08 + a_is_bigger_side = False + + # Where each of the two closest pos are determines how the extra three vertices are calculated. + if a_is_bigger_side == b_is_bigger_side: + if a_is_bigger_side: # Both closest pos on the bigger side + c1 = (a_po | b_po) + c2 = (a_po & b_po) + if (c1 & 0x01) == 0: + xsv_ext0 = xsb + xsv_ext1 = xsb - 1 + dx_ext0 = dx0 - 3 * SQUISH_CONSTANT_4D + dx_ext1 = dx0 + 1 - 2 * SQUISH_CONSTANT_4D + else: + xsv_ext0 = xsv_ext1 = xsb + 1 + dx_ext0 = dx0 - 1 - 3 * SQUISH_CONSTANT_4D + dx_ext1 = dx0 - 1 - 2 * SQUISH_CONSTANT_4D + + if (c1 & 0x02) == 0: + ysv_ext0 = ysb + ysv_ext1 = ysb - 1 + dy_ext0 = dy0 - 3 * SQUISH_CONSTANT_4D + dy_ext1 = dy0 + 1 - 2 * SQUISH_CONSTANT_4D + else: + ysv_ext0 = ysv_ext1 = ysb + 1 + dy_ext0 = dy0 - 1 - 3 * SQUISH_CONSTANT_4D + dy_ext1 = dy0 - 1 - 2 * SQUISH_CONSTANT_4D + + if (c1 & 0x04) == 0: + zsv_ext0 = zsb + zsv_ext1 = zsb - 1 + dz_ext0 = dz0 - 3 * SQUISH_CONSTANT_4D + dz_ext1 = dz0 + 1 - 2 * SQUISH_CONSTANT_4D + else: + zsv_ext0 = zsv_ext1 = zsb + 1 + dz_ext0 = dz0 - 1 - 3 * SQUISH_CONSTANT_4D + dz_ext1 = dz0 - 1 - 2 * SQUISH_CONSTANT_4D + + if (c1 & 0x08) == 0: + wsv_ext0 = wsb + wsv_ext1 = wsb - 1 + dw_ext0 = dw0 - 3 * SQUISH_CONSTANT_4D + dw_ext1 = dw0 + 1 - 2 * SQUISH_CONSTANT_4D + else: + wsv_ext0 = wsv_ext1 = wsb + 1 + dw_ext0 = dw0 - 1 - 3 * SQUISH_CONSTANT_4D + dw_ext1 = dw0 - 1 - 2 * SQUISH_CONSTANT_4D + + # One combination is a _permutation of (0,0,0,2) based on c2 + xsv_ext2 = xsb + ysv_ext2 = ysb + zsv_ext2 = zsb + wsv_ext2 = wsb + dx_ext2 = dx0 - 2 * SQUISH_CONSTANT_4D + dy_ext2 = dy0 - 2 * SQUISH_CONSTANT_4D + dz_ext2 = dz0 - 2 * SQUISH_CONSTANT_4D + dw_ext2 = dw0 - 2 * SQUISH_CONSTANT_4D + if (c2 & 0x01) != 0: + xsv_ext2 += 2 + dx_ext2 -= 2 + elif (c2 & 0x02) != 0: + ysv_ext2 += 2 + dy_ext2 -= 2 + elif (c2 & 0x04) != 0: + zsv_ext2 += 2 + dz_ext2 -= 2 + else: + wsv_ext2 += 2 + dw_ext2 -= 2 + + else: # Both closest pos on the smaller side + # One of the two extra pos is (0,0,0,0) + xsv_ext2 = xsb + ysv_ext2 = ysb + zsv_ext2 = zsb + wsv_ext2 = wsb + dx_ext2 = dx0 + dy_ext2 = dy0 + dz_ext2 = dz0 + dw_ext2 = dw0 + + # Other two pos are based on the omitted axes. + c = (a_po | b_po) + + if (c & 0x01) == 0: + xsv_ext0 = xsb - 1 + xsv_ext1 = xsb + dx_ext0 = dx0 + 1 - SQUISH_CONSTANT_4D + dx_ext1 = dx0 - SQUISH_CONSTANT_4D + else: + xsv_ext0 = xsv_ext1 = xsb + 1 + dx_ext0 = dx_ext1 = dx0 - 1 - SQUISH_CONSTANT_4D + + if (c & 0x02) == 0: + ysv_ext0 = ysv_ext1 = ysb + dy_ext0 = dy_ext1 = dy0 - SQUISH_CONSTANT_4D + if (c & 0x01) == 0x01: + ysv_ext0 -= 1 + dy_ext0 += 1 + else: + ysv_ext1 -= 1 + dy_ext1 += 1 + + else: + ysv_ext0 = ysv_ext1 = ysb + 1 + dy_ext0 = dy_ext1 = dy0 - 1 - SQUISH_CONSTANT_4D + + if (c & 0x04) == 0: + zsv_ext0 = zsv_ext1 = zsb + dz_ext0 = dz_ext1 = dz0 - SQUISH_CONSTANT_4D + if (c & 0x03) == 0x03: + zsv_ext0 -= 1 + dz_ext0 += 1 + else: + zsv_ext1 -= 1 + dz_ext1 += 1 + + else: + zsv_ext0 = zsv_ext1 = zsb + 1 + dz_ext0 = dz_ext1 = dz0 - 1 - SQUISH_CONSTANT_4D + + if (c & 0x08) == 0: + wsv_ext0 = wsb + wsv_ext1 = wsb - 1 + dw_ext0 = dw0 - SQUISH_CONSTANT_4D + dw_ext1 = dw0 + 1 - SQUISH_CONSTANT_4D + else: + wsv_ext0 = wsv_ext1 = wsb + 1 + dw_ext0 = dw_ext1 = dw0 - 1 - SQUISH_CONSTANT_4D + + else: # One po on each "side" + if a_is_bigger_side: + c1 = a_po + c2 = b_po + else: + c1 = b_po + c2 = a_po + + # Two contributions are the bigger-sided po with each 0 replaced with -1. + if (c1 & 0x01) == 0: + xsv_ext0 = xsb - 1 + xsv_ext1 = xsb + dx_ext0 = dx0 + 1 - SQUISH_CONSTANT_4D + dx_ext1 = dx0 - SQUISH_CONSTANT_4D + else: + xsv_ext0 = xsv_ext1 = xsb + 1 + dx_ext0 = dx_ext1 = dx0 - 1 - SQUISH_CONSTANT_4D + + if (c1 & 0x02) == 0: + ysv_ext0 = ysv_ext1 = ysb + dy_ext0 = dy_ext1 = dy0 - SQUISH_CONSTANT_4D + if (c1 & 0x01) == 0x01: + ysv_ext0 -= 1 + dy_ext0 += 1 + else: + ysv_ext1 -= 1 + dy_ext1 += 1 + + else: + ysv_ext0 = ysv_ext1 = ysb + 1 + dy_ext0 = dy_ext1 = dy0 - 1 - SQUISH_CONSTANT_4D + + if (c1 & 0x04) == 0: + zsv_ext0 = zsv_ext1 = zsb + dz_ext0 = dz_ext1 = dz0 - SQUISH_CONSTANT_4D + if (c1 & 0x03) == 0x03: + zsv_ext0 -= 1 + dz_ext0 += 1 + else: + zsv_ext1 -= 1 + dz_ext1 += 1 + + else: + zsv_ext0 = zsv_ext1 = zsb + 1 + dz_ext0 = dz_ext1 = dz0 - 1 - SQUISH_CONSTANT_4D + + if (c1 & 0x08) == 0: + wsv_ext0 = wsb + wsv_ext1 = wsb - 1 + dw_ext0 = dw0 - SQUISH_CONSTANT_4D + dw_ext1 = dw0 + 1 - SQUISH_CONSTANT_4D + else: + wsv_ext0 = wsv_ext1 = wsb + 1 + dw_ext0 = dw_ext1 = dw0 - 1 - SQUISH_CONSTANT_4D + + # One contribution is a _permutation of (0,0,0,2) based on the smaller-sided po + xsv_ext2 = xsb + ysv_ext2 = ysb + zsv_ext2 = zsb + wsv_ext2 = wsb + dx_ext2 = dx0 - 2 * SQUISH_CONSTANT_4D + dy_ext2 = dy0 - 2 * SQUISH_CONSTANT_4D + dz_ext2 = dz0 - 2 * SQUISH_CONSTANT_4D + dw_ext2 = dw0 - 2 * SQUISH_CONSTANT_4D + if (c2 & 0x01) != 0: + xsv_ext2 += 2 + dx_ext2 -= 2 + elif (c2 & 0x02) != 0: + ysv_ext2 += 2 + dy_ext2 -= 2 + elif (c2 & 0x04) != 0: + zsv_ext2 += 2 + dz_ext2 -= 2 + else: + wsv_ext2 += 2 + dw_ext2 -= 2 + + # Contribution (1,0,0,0) + dx1 = dx0 - 1 - SQUISH_CONSTANT_4D + dy1 = dy0 - 0 - SQUISH_CONSTANT_4D + dz1 = dz0 - 0 - SQUISH_CONSTANT_4D + dw1 = dw0 - 0 - SQUISH_CONSTANT_4D + attn1 = 2 - dx1 * dx1 - dy1 * dy1 - dz1 * dz1 - dw1 * dw1 + if attn1 > 0: + attn1 *= attn1 + value += attn1 * attn1 * extrapolate(xsb + 1, ysb + 0, zsb + 0, wsb + 0, dx1, dy1, dz1, dw1) + + # Contribution (0,1,0,0) + dx2 = dx0 - 0 - SQUISH_CONSTANT_4D + dy2 = dy0 - 1 - SQUISH_CONSTANT_4D + dz2 = dz1 + dw2 = dw1 + attn2 = 2 - dx2 * dx2 - dy2 * dy2 - dz2 * dz2 - dw2 * dw2 + if attn2 > 0: + attn2 *= attn2 + value += attn2 * attn2 * extrapolate(xsb + 0, ysb + 1, zsb + 0, wsb + 0, dx2, dy2, dz2, dw2) + + # Contribution (0,0,1,0) + dx3 = dx2 + dy3 = dy1 + dz3 = dz0 - 1 - SQUISH_CONSTANT_4D + dw3 = dw1 + attn3 = 2 - dx3 * dx3 - dy3 * dy3 - dz3 * dz3 - dw3 * dw3 + if attn3 > 0: + attn3 *= attn3 + value += attn3 * attn3 * extrapolate(xsb + 0, ysb + 0, zsb + 1, wsb + 0, dx3, dy3, dz3, dw3) + + # Contribution (0,0,0,1) + dx4 = dx2 + dy4 = dy1 + dz4 = dz1 + dw4 = dw0 - 1 - SQUISH_CONSTANT_4D + attn4 = 2 - dx4 * dx4 - dy4 * dy4 - dz4 * dz4 - dw4 * dw4 + if attn4 > 0: + attn4 *= attn4 + value += attn4 * attn4 * extrapolate(xsb + 0, ysb + 0, zsb + 0, wsb + 1, dx4, dy4, dz4, dw4) + + # Contribution (1,1,0,0) + dx5 = dx0 - 1 - 2 * SQUISH_CONSTANT_4D + dy5 = dy0 - 1 - 2 * SQUISH_CONSTANT_4D + dz5 = dz0 - 0 - 2 * SQUISH_CONSTANT_4D + dw5 = dw0 - 0 - 2 * SQUISH_CONSTANT_4D + attn5 = 2 - dx5 * dx5 - dy5 * dy5 - dz5 * dz5 - dw5 * dw5 + if attn5 > 0: + attn5 *= attn5 + value += attn5 * attn5 * extrapolate(xsb + 1, ysb + 1, zsb + 0, wsb + 0, dx5, dy5, dz5, dw5) + + # Contribution (1,0,1,0) + dx6 = dx0 - 1 - 2 * SQUISH_CONSTANT_4D + dy6 = dy0 - 0 - 2 * SQUISH_CONSTANT_4D + dz6 = dz0 - 1 - 2 * SQUISH_CONSTANT_4D + dw6 = dw0 - 0 - 2 * SQUISH_CONSTANT_4D + attn6 = 2 - dx6 * dx6 - dy6 * dy6 - dz6 * dz6 - dw6 * dw6 + if attn6 > 0: + attn6 *= attn6 + value += attn6 * attn6 * extrapolate(xsb + 1, ysb + 0, zsb + 1, wsb + 0, dx6, dy6, dz6, dw6) + + # Contribution (1,0,0,1) + dx7 = dx0 - 1 - 2 * SQUISH_CONSTANT_4D + dy7 = dy0 - 0 - 2 * SQUISH_CONSTANT_4D + dz7 = dz0 - 0 - 2 * SQUISH_CONSTANT_4D + dw7 = dw0 - 1 - 2 * SQUISH_CONSTANT_4D + attn7 = 2 - dx7 * dx7 - dy7 * dy7 - dz7 * dz7 - dw7 * dw7 + if attn7 > 0: + attn7 *= attn7 + value += attn7 * attn7 * extrapolate(xsb + 1, ysb + 0, zsb + 0, wsb + 1, dx7, dy7, dz7, dw7) + + # Contribution (0,1,1,0) + dx8 = dx0 - 0 - 2 * SQUISH_CONSTANT_4D + dy8 = dy0 - 1 - 2 * SQUISH_CONSTANT_4D + dz8 = dz0 - 1 - 2 * SQUISH_CONSTANT_4D + dw8 = dw0 - 0 - 2 * SQUISH_CONSTANT_4D + attn8 = 2 - dx8 * dx8 - dy8 * dy8 - dz8 * dz8 - dw8 * dw8 + if attn8 > 0: + attn8 *= attn8 + value += attn8 * attn8 * extrapolate(xsb + 0, ysb + 1, zsb + 1, wsb + 0, dx8, dy8, dz8, dw8) + + # Contribution (0,1,0,1) + dx9 = dx0 - 0 - 2 * SQUISH_CONSTANT_4D + dy9 = dy0 - 1 - 2 * SQUISH_CONSTANT_4D + dz9 = dz0 - 0 - 2 * SQUISH_CONSTANT_4D + dw9 = dw0 - 1 - 2 * SQUISH_CONSTANT_4D + attn9 = 2 - dx9 * dx9 - dy9 * dy9 - dz9 * dz9 - dw9 * dw9 + if attn9 > 0: + attn9 *= attn9 + value += attn9 * attn9 * extrapolate(xsb + 0, ysb + 1, zsb + 0, wsb + 1, dx9, dy9, dz9, dw9) + + # Contribution (0,0,1,1) + dx10 = dx0 - 0 - 2 * SQUISH_CONSTANT_4D + dy10 = dy0 - 0 - 2 * SQUISH_CONSTANT_4D + dz10 = dz0 - 1 - 2 * SQUISH_CONSTANT_4D + dw10 = dw0 - 1 - 2 * SQUISH_CONSTANT_4D + attn10 = 2 - dx10 * dx10 - dy10 * dy10 - dz10 * dz10 - dw10 * dw10 + if attn10 > 0: + attn10 *= attn10 + value += attn10 * attn10 * extrapolate(xsb + 0, ysb + 0, zsb + 1, wsb + 1, dx10, dy10, dz10, dw10) + + else: # We're inside the second dispentachoron (Rectified 4-Simplex) + a_is_bigger_side = True + b_is_bigger_side = True + + # Decide between (0,0,1,1) and (1,1,0,0) + if xins + yins < zins + wins: + a_score = xins + yins + a_po = 0x0C + else: + a_score = zins + wins + a_po = 0x03 + + # Decide between (0,1,0,1) and (1,0,1,0) + if xins + zins < yins + wins: + b_score = xins + zins + b_po = 0x0A + else: + b_score = yins + wins + b_po = 0x05 + + # Closer between (0,1,1,0) and (1,0,0,1) will replace the further of a and b, if closer. + if xins + wins < yins + zins: + score = xins + wins + if a_score <= b_score and score < b_score: + b_score = score + b_po = 0x06 + elif a_score > b_score and score < a_score: + a_score = score + a_po = 0x06 + + else: + score = yins + zins + if a_score <= b_score and score < b_score: + b_score = score + b_po = 0x09 + elif a_score > b_score and score < a_score: + a_score = score + a_po = 0x09 + + # Decide if (0,1,1,1) is closer. + p1 = 3 - in_sum + xins + if a_score <= b_score and p1 < b_score: + b_score = p1 + b_po = 0x0E + b_is_bigger_side = False + elif a_score > b_score and p1 < a_score: + a_score = p1 + a_po = 0x0E + a_is_bigger_side = False + + # Decide if (1,0,1,1) is closer. + p2 = 3 - in_sum + yins + if a_score <= b_score and p2 < b_score: + b_score = p2 + b_po = 0x0D + b_is_bigger_side = False + elif a_score > b_score and p2 < a_score: + a_score = p2 + a_po = 0x0D + a_is_bigger_side = False + + # Decide if (1,1,0,1) is closer. + p3 = 3 - in_sum + zins + if a_score <= b_score and p3 < b_score: + b_score = p3 + b_po = 0x0B + b_is_bigger_side = False + elif a_score > b_score and p3 < a_score: + a_score = p3 + a_po = 0x0B + a_is_bigger_side = False + + # Decide if (1,1,1,0) is closer. + p4 = 3 - in_sum + wins + if a_score <= b_score and p4 < b_score: + b_po = 0x07 + b_is_bigger_side = False + elif a_score > b_score and p4 < a_score: + a_po = 0x07 + a_is_bigger_side = False + + # Where each of the two closest pos are determines how the extra three vertices are calculated. + if a_is_bigger_side == b_is_bigger_side: + if a_is_bigger_side: # Both closest pos on the bigger side + c1 = (a_po & b_po) + c2 = (a_po | b_po) + + # Two contributions are _permutations of (0,0,0,1) and (0,0,0,2) based on c1 + xsv_ext0 = xsv_ext1 = xsb + ysv_ext0 = ysv_ext1 = ysb + zsv_ext0 = zsv_ext1 = zsb + wsv_ext0 = wsv_ext1 = wsb + dx_ext0 = dx0 - SQUISH_CONSTANT_4D + dy_ext0 = dy0 - SQUISH_CONSTANT_4D + dz_ext0 = dz0 - SQUISH_CONSTANT_4D + dw_ext0 = dw0 - SQUISH_CONSTANT_4D + dx_ext1 = dx0 - 2 * SQUISH_CONSTANT_4D + dy_ext1 = dy0 - 2 * SQUISH_CONSTANT_4D + dz_ext1 = dz0 - 2 * SQUISH_CONSTANT_4D + dw_ext1 = dw0 - 2 * SQUISH_CONSTANT_4D + if (c1 & 0x01) != 0: + xsv_ext0 += 1 + dx_ext0 -= 1 + xsv_ext1 += 2 + dx_ext1 -= 2 + elif (c1 & 0x02) != 0: + ysv_ext0 += 1 + dy_ext0 -= 1 + ysv_ext1 += 2 + dy_ext1 -= 2 + elif (c1 & 0x04) != 0: + zsv_ext0 += 1 + dz_ext0 -= 1 + zsv_ext1 += 2 + dz_ext1 -= 2 + else: + wsv_ext0 += 1 + dw_ext0 -= 1 + wsv_ext1 += 2 + dw_ext1 -= 2 + + # One contribution is a _permutation of (1,1,1,-1) based on c2 + xsv_ext2 = xsb + 1 + ysv_ext2 = ysb + 1 + zsv_ext2 = zsb + 1 + wsv_ext2 = wsb + 1 + dx_ext2 = dx0 - 1 - 2 * SQUISH_CONSTANT_4D + dy_ext2 = dy0 - 1 - 2 * SQUISH_CONSTANT_4D + dz_ext2 = dz0 - 1 - 2 * SQUISH_CONSTANT_4D + dw_ext2 = dw0 - 1 - 2 * SQUISH_CONSTANT_4D + if (c2 & 0x01) == 0: + xsv_ext2 -= 2 + dx_ext2 += 2 + elif (c2 & 0x02) == 0: + ysv_ext2 -= 2 + dy_ext2 += 2 + elif (c2 & 0x04) == 0: + zsv_ext2 -= 2 + dz_ext2 += 2 + else: + wsv_ext2 -= 2 + dw_ext2 += 2 + + else: # Both closest pos on the smaller side + # One of the two extra pos is (1,1,1,1) + xsv_ext2 = xsb + 1 + ysv_ext2 = ysb + 1 + zsv_ext2 = zsb + 1 + wsv_ext2 = wsb + 1 + dx_ext2 = dx0 - 1 - 4 * SQUISH_CONSTANT_4D + dy_ext2 = dy0 - 1 - 4 * SQUISH_CONSTANT_4D + dz_ext2 = dz0 - 1 - 4 * SQUISH_CONSTANT_4D + dw_ext2 = dw0 - 1 - 4 * SQUISH_CONSTANT_4D + + # Other two pos are based on the shared axes. + c = (a_po & b_po) + if (c & 0x01) != 0: + xsv_ext0 = xsb + 2 + xsv_ext1 = xsb + 1 + dx_ext0 = dx0 - 2 - 3 * SQUISH_CONSTANT_4D + dx_ext1 = dx0 - 1 - 3 * SQUISH_CONSTANT_4D + else: + xsv_ext0 = xsv_ext1 = xsb + dx_ext0 = dx_ext1 = dx0 - 3 * SQUISH_CONSTANT_4D + + if (c & 0x02) != 0: + ysv_ext0 = ysv_ext1 = ysb + 1 + dy_ext0 = dy_ext1 = dy0 - 1 - 3 * SQUISH_CONSTANT_4D + if (c & 0x01) == 0: + ysv_ext0 += 1 + dy_ext0 -= 1 + else: + ysv_ext1 += 1 + dy_ext1 -= 1 + + else: + ysv_ext0 = ysv_ext1 = ysb + dy_ext0 = dy_ext1 = dy0 - 3 * SQUISH_CONSTANT_4D + + if (c & 0x04) != 0: + zsv_ext0 = zsv_ext1 = zsb + 1 + dz_ext0 = dz_ext1 = dz0 - 1 - 3 * SQUISH_CONSTANT_4D + if (c & 0x03) == 0: + zsv_ext0 += 1 + dz_ext0 -= 1 + else: + zsv_ext1 += 1 + dz_ext1 -= 1 + + else: + zsv_ext0 = zsv_ext1 = zsb + dz_ext0 = dz_ext1 = dz0 - 3 * SQUISH_CONSTANT_4D + + if (c & 0x08) != 0: + wsv_ext0 = wsb + 1 + wsv_ext1 = wsb + 2 + dw_ext0 = dw0 - 1 - 3 * SQUISH_CONSTANT_4D + dw_ext1 = dw0 - 2 - 3 * SQUISH_CONSTANT_4D + else: + wsv_ext0 = wsv_ext1 = wsb + dw_ext0 = dw_ext1 = dw0 - 3 * SQUISH_CONSTANT_4D + + else: # One po on each "side" + if a_is_bigger_side: + c1 = a_po + c2 = b_po + else: + c1 = b_po + c2 = a_po + + # Two contributions are the bigger-sided po with each 1 replaced with 2. + if (c1 & 0x01) != 0: + xsv_ext0 = xsb + 2 + xsv_ext1 = xsb + 1 + dx_ext0 = dx0 - 2 - 3 * SQUISH_CONSTANT_4D + dx_ext1 = dx0 - 1 - 3 * SQUISH_CONSTANT_4D + else: + xsv_ext0 = xsv_ext1 = xsb + dx_ext0 = dx_ext1 = dx0 - 3 * SQUISH_CONSTANT_4D + + if (c1 & 0x02) != 0: + ysv_ext0 = ysv_ext1 = ysb + 1 + dy_ext0 = dy_ext1 = dy0 - 1 - 3 * SQUISH_CONSTANT_4D + if (c1 & 0x01) == 0: + ysv_ext0 += 1 + dy_ext0 -= 1 + else: + ysv_ext1 += 1 + dy_ext1 -= 1 + + else: + ysv_ext0 = ysv_ext1 = ysb + dy_ext0 = dy_ext1 = dy0 - 3 * SQUISH_CONSTANT_4D + + if (c1 & 0x04) != 0: + zsv_ext0 = zsv_ext1 = zsb + 1 + dz_ext0 = dz_ext1 = dz0 - 1 - 3 * SQUISH_CONSTANT_4D + if (c1 & 0x03) == 0: + zsv_ext0 += 1 + dz_ext0 -= 1 + else: + zsv_ext1 += 1 + dz_ext1 -= 1 + + else: + zsv_ext0 = zsv_ext1 = zsb + dz_ext0 = dz_ext1 = dz0 - 3 * SQUISH_CONSTANT_4D + + if (c1 & 0x08) != 0: + wsv_ext0 = wsb + 1 + wsv_ext1 = wsb + 2 + dw_ext0 = dw0 - 1 - 3 * SQUISH_CONSTANT_4D + dw_ext1 = dw0 - 2 - 3 * SQUISH_CONSTANT_4D + else: + wsv_ext0 = wsv_ext1 = wsb + dw_ext0 = dw_ext1 = dw0 - 3 * SQUISH_CONSTANT_4D + + # One contribution is a _permutation of (1,1,1,-1) based on the smaller-sided po + xsv_ext2 = xsb + 1 + ysv_ext2 = ysb + 1 + zsv_ext2 = zsb + 1 + wsv_ext2 = wsb + 1 + dx_ext2 = dx0 - 1 - 2 * SQUISH_CONSTANT_4D + dy_ext2 = dy0 - 1 - 2 * SQUISH_CONSTANT_4D + dz_ext2 = dz0 - 1 - 2 * SQUISH_CONSTANT_4D + dw_ext2 = dw0 - 1 - 2 * SQUISH_CONSTANT_4D + if (c2 & 0x01) == 0: + xsv_ext2 -= 2 + dx_ext2 += 2 + elif (c2 & 0x02) == 0: + ysv_ext2 -= 2 + dy_ext2 += 2 + elif (c2 & 0x04) == 0: + zsv_ext2 -= 2 + dz_ext2 += 2 + else: + wsv_ext2 -= 2 + dw_ext2 += 2 + + # Contribution (1,1,1,0) + dx4 = dx0 - 1 - 3 * SQUISH_CONSTANT_4D + dy4 = dy0 - 1 - 3 * SQUISH_CONSTANT_4D + dz4 = dz0 - 1 - 3 * SQUISH_CONSTANT_4D + dw4 = dw0 - 3 * SQUISH_CONSTANT_4D + attn4 = 2 - dx4 * dx4 - dy4 * dy4 - dz4 * dz4 - dw4 * dw4 + if attn4 > 0: + attn4 *= attn4 + value += attn4 * attn4 * extrapolate(xsb + 1, ysb + 1, zsb + 1, wsb + 0, dx4, dy4, dz4, dw4) + + # Contribution (1,1,0,1) + dx3 = dx4 + dy3 = dy4 + dz3 = dz0 - 3 * SQUISH_CONSTANT_4D + dw3 = dw0 - 1 - 3 * SQUISH_CONSTANT_4D + attn3 = 2 - dx3 * dx3 - dy3 * dy3 - dz3 * dz3 - dw3 * dw3 + if attn3 > 0: + attn3 *= attn3 + value += attn3 * attn3 * extrapolate(xsb + 1, ysb + 1, zsb + 0, wsb + 1, dx3, dy3, dz3, dw3) + + # Contribution (1,0,1,1) + dx2 = dx4 + dy2 = dy0 - 3 * SQUISH_CONSTANT_4D + dz2 = dz4 + dw2 = dw3 + attn2 = 2 - dx2 * dx2 - dy2 * dy2 - dz2 * dz2 - dw2 * dw2 + if attn2 > 0: + attn2 *= attn2 + value += attn2 * attn2 * extrapolate(xsb + 1, ysb + 0, zsb + 1, wsb + 1, dx2, dy2, dz2, dw2) + + # Contribution (0,1,1,1) + dx1 = dx0 - 3 * SQUISH_CONSTANT_4D + dz1 = dz4 + dy1 = dy4 + dw1 = dw3 + attn1 = 2 - dx1 * dx1 - dy1 * dy1 - dz1 * dz1 - dw1 * dw1 + if attn1 > 0: + attn1 *= attn1 + value += attn1 * attn1 * extrapolate(xsb + 0, ysb + 1, zsb + 1, wsb + 1, dx1, dy1, dz1, dw1) + + # Contribution (1,1,0,0) + dx5 = dx0 - 1 - 2 * SQUISH_CONSTANT_4D + dy5 = dy0 - 1 - 2 * SQUISH_CONSTANT_4D + dz5 = dz0 - 0 - 2 * SQUISH_CONSTANT_4D + dw5 = dw0 - 0 - 2 * SQUISH_CONSTANT_4D + attn5 = 2 - dx5 * dx5 - dy5 * dy5 - dz5 * dz5 - dw5 * dw5 + if attn5 > 0: + attn5 *= attn5 + value += attn5 * attn5 * extrapolate(xsb + 1, ysb + 1, zsb + 0, wsb + 0, dx5, dy5, dz5, dw5) + + # Contribution (1,0,1,0) + dx6 = dx0 - 1 - 2 * SQUISH_CONSTANT_4D + dy6 = dy0 - 0 - 2 * SQUISH_CONSTANT_4D + dz6 = dz0 - 1 - 2 * SQUISH_CONSTANT_4D + dw6 = dw0 - 0 - 2 * SQUISH_CONSTANT_4D + attn6 = 2 - dx6 * dx6 - dy6 * dy6 - dz6 * dz6 - dw6 * dw6 + if attn6 > 0: + attn6 *= attn6 + value += attn6 * attn6 * extrapolate(xsb + 1, ysb + 0, zsb + 1, wsb + 0, dx6, dy6, dz6, dw6) + + # Contribution (1,0,0,1) + dx7 = dx0 - 1 - 2 * SQUISH_CONSTANT_4D + dy7 = dy0 - 0 - 2 * SQUISH_CONSTANT_4D + dz7 = dz0 - 0 - 2 * SQUISH_CONSTANT_4D + dw7 = dw0 - 1 - 2 * SQUISH_CONSTANT_4D + attn7 = 2 - dx7 * dx7 - dy7 * dy7 - dz7 * dz7 - dw7 * dw7 + if attn7 > 0: + attn7 *= attn7 + value += attn7 * attn7 * extrapolate(xsb + 1, ysb + 0, zsb + 0, wsb + 1, dx7, dy7, dz7, dw7) + + # Contribution (0,1,1,0) + dx8 = dx0 - 0 - 2 * SQUISH_CONSTANT_4D + dy8 = dy0 - 1 - 2 * SQUISH_CONSTANT_4D + dz8 = dz0 - 1 - 2 * SQUISH_CONSTANT_4D + dw8 = dw0 - 0 - 2 * SQUISH_CONSTANT_4D + attn8 = 2 - dx8 * dx8 - dy8 * dy8 - dz8 * dz8 - dw8 * dw8 + if attn8 > 0: + attn8 *= attn8 + value += attn8 * attn8 * extrapolate(xsb + 0, ysb + 1, zsb + 1, wsb + 0, dx8, dy8, dz8, dw8) + + # Contribution (0,1,0,1) + dx9 = dx0 - 0 - 2 * SQUISH_CONSTANT_4D + dy9 = dy0 - 1 - 2 * SQUISH_CONSTANT_4D + dz9 = dz0 - 0 - 2 * SQUISH_CONSTANT_4D + dw9 = dw0 - 1 - 2 * SQUISH_CONSTANT_4D + attn9 = 2 - dx9 * dx9 - dy9 * dy9 - dz9 * dz9 - dw9 * dw9 + if attn9 > 0: + attn9 *= attn9 + value += attn9 * attn9 * extrapolate(xsb + 0, ysb + 1, zsb + 0, wsb + 1, dx9, dy9, dz9, dw9) + + # Contribution (0,0,1,1) + dx10 = dx0 - 0 - 2 * SQUISH_CONSTANT_4D + dy10 = dy0 - 0 - 2 * SQUISH_CONSTANT_4D + dz10 = dz0 - 1 - 2 * SQUISH_CONSTANT_4D + dw10 = dw0 - 1 - 2 * SQUISH_CONSTANT_4D + attn10 = 2 - dx10 * dx10 - dy10 * dy10 - dz10 * dz10 - dw10 * dw10 + if attn10 > 0: + attn10 *= attn10 + value += attn10 * attn10 * extrapolate(xsb + 0, ysb + 0, zsb + 1, wsb + 1, dx10, dy10, dz10, dw10) + + # First extra vertex + attn_ext0 = 2 - dx_ext0 * dx_ext0 - dy_ext0 * dy_ext0 - dz_ext0 * dz_ext0 - dw_ext0 * dw_ext0 + if attn_ext0 > 0: + attn_ext0 *= attn_ext0 + value += attn_ext0 * attn_ext0 * extrapolate(xsv_ext0, ysv_ext0, zsv_ext0, wsv_ext0, dx_ext0, dy_ext0, dz_ext0, dw_ext0) + + # Second extra vertex + attn_ext1 = 2 - dx_ext1 * dx_ext1 - dy_ext1 * dy_ext1 - dz_ext1 * dz_ext1 - dw_ext1 * dw_ext1 + if attn_ext1 > 0: + attn_ext1 *= attn_ext1 + value += attn_ext1 * attn_ext1 * extrapolate(xsv_ext1, ysv_ext1, zsv_ext1, wsv_ext1, dx_ext1, dy_ext1, dz_ext1, dw_ext1) + + # Third extra vertex + attn_ext2 = 2 - dx_ext2 * dx_ext2 - dy_ext2 * dy_ext2 - dz_ext2 * dz_ext2 - dw_ext2 * dw_ext2 + if attn_ext2 > 0: + attn_ext2 *= attn_ext2 + value += attn_ext2 * attn_ext2 * extrapolate(xsv_ext2, ysv_ext2, zsv_ext2, wsv_ext2, dx_ext2, dy_ext2, dz_ext2, dw_ext2) + + return value / NORM_CONSTANT_4D diff --git a/python/quiver_plot.py b/python/quiver_plot.py new file mode 100644 index 0000000..c10f73d --- /dev/null +++ b/python/quiver_plot.py @@ -0,0 +1,121 @@ +# -*- coding: utf-8 -*- + +# Scratch code for visualizing a 3D vector field with arrows + +import functools + +import numpy as np +import vispy +import vispy.scene +from vispy.scene import visuals + +import opensimplex + +def gradient(fn, eps=1e-3): + eps_inv = 1.0 / eps + # Numerical gradient by finite differences + def _grad(x, y, z): + p = fn(x, y, z) + p_dx = fn(x+eps, y, z) + p_dy = fn(x, y+eps, z) + p_dz = fn(x, y, z+eps) + return (p_dx - p)*eps_inv, (p_dy - p)*eps_inv, (p_dz - p)*eps_inv + return _grad + +def curl_3d(grads): + # 'grads' should be of shape (..., 3, 3); + # 2nd-to-last dimension is gradient of potential x, y, z. + # Last dimension is gradient of that w.r.t. x, y, and z. + cx = grads[..., 2, 1] - grads[..., 1, 2] + cy = grads[..., 0, 2] - grads[..., 2, 0] + cz = grads[..., 1, 0] - grads[..., 0, 1] + return np.stack((cx, cy, cz), axis=-1) + +def generate(grid, t=0): + x_spx = opensimplex.OpenSimplex(seed=0) + y_spx = opensimplex.OpenSimplex(seed=12345) + z_spx = opensimplex.OpenSimplex(seed=45678) + # grad_x = gradient(x_spx.noise3d) + # grad_y = gradient(y_spx.noise3d) + # grad_z = gradient(z_spx.noise3d) + # grad_x = gradient(functools.partial(x_spx.noise4d, t)) + # grad_y = gradient(functools.partial(y_spx.noise4d, t)) + # grad_z = gradient(functools.partial(z_spx.noise4d, t)) + grad_x = gradient(lambda x,y,z: x_spx.noise3d(x + t, y, z)) + grad_y = gradient(lambda x,y,z: y_spx.noise3d(x + t, y, z)) + grad_z = gradient(lambda x,y,z: z_spx.noise3d(x + t, y, z)) + grads = np.array([ + (grad_x(x,y,z), grad_y(x,y,z), grad_z(x,y,z)) + for (x,y,z) in grid + ]) + curl = curl_3d(grads) + return grid, curl + #pos = np.random.normal(size=(count, 3), scale=0.2) + # one could stop here for the data generation, the rest is just to make the + # data look more interesting. Copied over from magnify.py + #return pos + # centers = np.random.normal(size=(50, 3)) + # indexes = np.random.normal(size=count, loc=centers.shape[0]/2., + # scale=centers.shape[0]/3.) + # indexes = np.clip(indexes, 0, centers.shape[0]-1).astype(int) + # scales = 2**(np.linspace(-2, 0.5, centers.shape[0]))[indexes][:, np.newaxis] + # pos *= scales + # pos += centers[indexes] + # return pos + +class Data(object): + def __init__(self, view): + self.view = view + self.visual = vispy.scene.visuals.Arrow( + color=(1,1,1,0.75), + connect='segments', + arrow_size=8, + arrow_type="triangle_30", + ) + self.view.add(self.visual) + self.view.camera = 'turntable' # or try 'arcball' + #view.camera = 'arcball' + s = 1 + count = 8 + xs = zs = np.linspace(-s, s, count) + ys = np.array([0]) + self.grid = np.array([i.flatten() for i in np.meshgrid(xs,ys,zs)]).T + self.update() + def update(self, ev=None): + t = 0 if ev is None else ev.elapsed + grid, curl = generate(self.grid, t*0.4) + a1 = grid + a2 = grid + curl*0.1 + # So I guess if I give argument 'pos', it needs to be of shape + # (2*N, 3) and consist of alternating starting & ending vertices? + # The docs don't say a thing about this. + self.visual.set_data( + pos=np.hstack((a1, a2)).reshape(grid.shape[0]*2, -1), + arrows=np.hstack((a1, a2)), + ) + # self.scatter = visuals.Markers() + # self.scatter.set_data(self.d, edge_color=None, face_color=(1, 0.5, 1, .5), size=4) + # view.add(self.scatter) + #m = np.array([[np.cos(t), np.sin(t*1.01), np.cos(t*1.02)]]) + #d2 = self.d*m + #self.scatter.set_data(d2, edge_color=None, face_color=(1, 0.5, 1, .5), size=4) + +def main(): + # + # Make a canvas and add simple view + # + canvas = vispy.scene.SceneCanvas(keys='interactive', show=True) + view = canvas.central_widget.add_view() + import sys + # add a colored 3D axis for orientation + axis = visuals.XYZAxis(parent=view.scene) + timer = vispy.app.Timer() + da = Data(view) + # Problem is how slow update() is: + timer.connect(da.update) + timer.start(0.05) + if sys.flags.interactive != 1: + vispy.app.run() + +if __name__ == '__main__': + main() diff --git a/python/shell.nix b/python/shell.nix new file mode 100644 index 0000000..811810a --- /dev/null +++ b/python/shell.nix @@ -0,0 +1,15 @@ +{ pkgs ? import {} }: +pkgs.stdenv.mkDerivation rec { + name = "vispy-test"; + + # PyQt5 should be fine but seems to have other issues + buildInputs = [ + (pkgs.python3.withPackages + (ps: [ + ps.vispy + ps.wxPython_4_0 + ps.sympy + ps.jupyterlab + ])) + ]; +}