Перевод поста Джона Маклуна "Droste Effect with Mathematica". Код, приведенный в статье, можно скачать в конце поста.
Выражаю огромную благодарность Кириллу Гузенко за помощь в переводе.

Эффект Дросте (wiki) представляет собой рекурсивное включение изображением самого в себя. Название происходит от какао-порошка Droste, который в 1904 году продавался в упаковке, на которой была изображена медсестра, которая держала коробку, на которой была медсестра, ну и так далее. Самая простая реализация — отмасштабировать и трансформировать изображение, а затем поместить его на свою немодифицированную точную копию, затем начать процесс снова. Взгляните на демонстрацию, в которой используется оригинальные иллюстрации упаковки Droste. Однако значительно более интересных результатов можно достичь, если использовать теорию функций ко?мплексного переменного (ТФКП). Эшер М. К. был первым, кто популяризировал идею конформных отображений применительно к изображениям, однако с помощью компьютеров мы легко можем реализовать эту идею на фотографиях для получения чего-то подобного:

A photograph conformally mapped in Mathematica

Да, идея не нова. Однако, когда я решил реализовать подобный эффект, те методы, что я находил в сети, казались мне неудовлетворительными. Одни предполагали много ручной работы типа «скопируй-вставь» над изображениями, в других были проблемы с несоответствием разрешений на стыках частей изображений. И, как обычно, я пришёл к тому, что тут надо использовать Mathematica, или в данном случае gridMathematica.

По сути идея проста. Мы создаем изображение, где пиксель в позиции {x,y} получаемого изображения получается из пикселя в позиции f[{x,y}] нашего исходного изображения. Магия заключается в выборе самой функции f. Но прежде чем мы добираемся до этого, нам нужно кое-что подготовить.

В настоящее время Mathematica не содержит операции по произвольной трансформации изображений, так что для начала я должен их реализовать (Пост написан Джоном в 7-й версии Wolfram Mathematica, в 8-й версии для этих целей появилась специальная функция ImageTransformation, которая заметно уменьшает объем кода в статье Джона — прим. ред.). Важный момент заключается в постоянном вычислении цвета между пикселями для предотвращения проблем соответствия разрешения и ненужной пикселизации в увеличенных областях изображения. Мой метод заключается в линейной интерполяции всех цветов пикселей в каждом канале RGB для изображения. Это наиболее вычислительно сложный подход, особенно с исходными изображениями в 10 мегапикселей, однако он дает необходимое мне качество.

Creating a linear interpolation of all pixel colors

Чтобы компенсировать подобный подход, я настроил всё для параллельных вычислений. Для начала, я использовал встроенные в Mathematica инструменты распараллеливания. Таким образом, вместо того, чтобы просто использовать Table для генерации сетки пикселей, я создал функцию, вызывающую Table для генерации слоя окончательного изображения, которая вызывается параллельно посредством ParallelTable для перераспределения работы между несколькими ядрами процессора. Так на 3-4 строчки больше, однако это лишь половина работы по распараллеливанию.

Generating a slice of the final image

Distributing the task across multiple CPUs

Затем, я передаю исходное изображение каждому процессору и задаю функцию интерполяции. Эта часть весьма требовательна к вычислениям, так что я хочу проделать это лишь раз и после этого иметь возможность делать любые виды изображений из исходного без проделывая повторно этой работы. Распараллеливание осуществляется здесь довольно просто: запускаются доступные ядра, осуществляется распределение определения программы на все ядра, а затем используется функция ParallelEvaluate для отправки изображения по ссылке с запросом на создание функции интерполяции.

Parallelizing and interpolating

Существует изящный приём в передаче изображения в виде строки с содержанием файла JPG вместо того, чтобы передавать реальные несжатые данные. Получается гораздо меньший объект, который, соответственно, быстрее передаётся.

При подобной настройке я легко могу привлекать дополнительные ядра компьютеров из своего офиса для ускорения вычислений посредством gridMathematica.

В код нужно ещё кое-что добавить. Так я избавляюсь от необходимости в том, чтобы исходное изображение было обрезано в правильном соотношении и с правильным центрированием.

Cropping the image

Вот то, на что похоже мое исходное изображение после обрезки:

The cropped image

Что ж, всё скучное позади, теперь самое время поразвлечься.

Скручивание, которое я использую, основывается на специфических свойствах степенной функции Power на комплексной плоскости. Мы можем представить наши пары координат {x,y} как части комплексного числа p=x+Iy, использовать отображение f[p]:=pc и вернуться обратно к декартовым координатам. Довольно трудно разобраться, каким брать значение c, особенно если использовать более общую модель a+(b+c p)d, где {a,b,c,d} являются комплексными числами. Так что я решил обратиться к Wolfram Demonstrations Project, нашёл там демонстрацию с конформными отображениями и немного изменил код, чтобы получить интересующее меня отображение.

Поэкспериментировав, я выяснил влияния параметров, настроил калибровку, получил хорошие начальные значения и задал именованные опции в необходимых местах той самой магической формулы. В итоге, у меня получилось такая конструкция для скручивания.

Placing parameters and default values into the magic formula

Я обеспечиваю подобное самовоспроизведение перескакиванием внутрь при попытке получить доступ к пикселю за пределами изображения, а затем наружу во время получения доступа к пикселю внутри определённой области.

Из соображений эстетики я бы хотел закрыть соединения между изображениями для придания им большей естественности. В своём примере я использовал картинную рамку, однако я встречал примеры, где использовались дверные проёмы, окна, компьютерные мониторы и тому подобное.

Автоматизирование этого процесса позволяет избежать ручного копирования и последующей вставки, а также позволяет избежать проблем с несоответствием разрешений различных частей изображения.

Programatically covering the image join

Всё, вся трудоёмкая работа уже позади, теперь самое время побаловаться. Инициализируйте ядра с оригинальным изображением, задав координаты внутренней области (вы легко можете получить их, используя coordinate picker в палитре 2D Drawing):

Initialize the kernels with the source image

Теперь сгенерируем изображение высотой в 400 пикселей в высоту:

Generating a 400-pixel image

Существует множество сочетаний параметров, например, как этот с двойной спиралью:

Image created from double-helix parameters

Одна спираль в противоположном направлении:

Spiral in the opposite direction

Без спиралей, одна лишь репликация:

No spiral, just replication

Этот код создаёт восьмиугольники путем создания двух копий на спираль:

Two copies per spiral, resulting in octagonal forms

А этот без рекурсии и спиралей — всего две копии, свёрнутые вместе:

Just two copies wrapped together

И вот, кульминация — целый фильм в DVD качестве, показывающий эффект Дросте закручивания изображения в рекурсивную спираль.



Создание такого видео требует большого количества вычислений — вызов трёх интерполяционных функций для каждого из 10 млн. пикселей примерно по 400 тыс. раз для каждого из 60 кадров. Это тот момент, когда мои усилия по распараллеливанию действительно окупаются. Конечно, будучи в Wolfram Research, у меня есть под под рукой несколько лицензий gridMathematica, сообщающихся посредством Wolfram Lightweight Grid Manager (Начиная с 7-й версии Wolfram Mathematica, функционал по распараллеливанию включен в ядро Wolfram Language и не требует дополнительных программ при работе, он способен использовать все ядра, доступные на вашем компьютере. Однако, если есть необходимость подключиться к кластеру, для этого потребуется gridMathematica — прим. ред.). Я открываю Parallel Configuration preferences, и все они волшебным образом появляются. Несколько кликов, и вот у меня появляются ещё 16 ядер, вдобавок к двум на моём маленьком ноутбуке, и без каких бы то ни было изменений в коде фильм генерируется в 8 раз быстрее (и в 16 раз быстрее чем без распараллеливания). Я точно не знаю, где именно работает этот код, но так как программа и исходные изображения передаются автоматически, то мне это не особо и интересно. Можете посмотреть экранную демонстрацию того, как всё это настраивается.

Exporting the movie

Этот код можно по разному модифицировать; например, изменение ReplicateRegion позволит использовать непрямоугольные кадры, в виде кругов, например, или можно даже изменять цвета и прозрачность репликаций. Попробуйте самостоятельно поэкспериментировать и посмотреть, что может из этого получиться.

Код, использованный в статье
InitializeSources[source_,p1_,p2_]:=Quiet[Block[{imgbytes=Import[source,"Byte"],sourcebytestream},If[imgbytes===$Failed,$Failed,sourcebytestream=FromCharacterCode[imgbytes];
LaunchKernels[];
DistributeDefinitions[OriginalValueFns,DrostifyRegion,ReplicateRegion,TransformCoordinates,CropData];
ParallelEvaluate[$ImageInterpolationFn=OriginalValueFns[CropData[Reverse@Developer`ToPackedArray[N[ImportString[#,{"JPG","Data"}]]],p1,p2]/255.];]&[sourcebytestream]]]]

OriginalValueFns[data_]:=($AspectRatio=1/Apply[Divide,Most[Dimensions[data]]];
Apply[Function,{{x,y},If[Abs[x]>1||Abs[y]>$AspectRatio,{1.,1.,1.},#]&[Table[ListInterpolation[data[[All,All,channel]],{{-1,1},$AspectRatio {-1,1}},InterpolationOrder->1][x,y],{channel,1,3}]]}]);

DrostifyRegion[start_,end_,res_,opts:OptionsPattern[]]:=Image[Table[Apply[$ImageInterpolationFn,TransformCoordinates[{x,y},opts]],{x,end,start,-2/(res-1)},{y,-$AspectRatio/$AspectRatio,2/res}]];

DrosteImage[resolution_,opts:OptionsPattern[]]:=ImageAssemble[ParallelTable[{DrostifyRegion[-1+2 (slice-1)/#,-1.+2 slice/#,resolution,opts]},{slice,#,1,-1}]]&[Max[Length[Kernels[]],1]];

CropData[data_,r1_,r2_]:=Block[{center,xlo,xhi,ylo,yhi,innerdims},center=Mean[{r1,r2}];
(*Find the center of the selected rectangle*)
$DrosteScale=Max[Flatten[{Abs[r2-center]/(center-{1,1}),Abs[r1-center]/Abs[Reverse@Take[Dimensions[data],2]-center]}]];

(*Find the scaling of the cropped image to the rectangle*)innerdims=Abs[r1-r2]/$DrosteScale;
{{ylo,xlo},{yhi,xhi}}=Round[{center-innerdims/2,center+innerdims/2}];
Return[data[[xlo;;xhi,ylo;;yhi,All]]]];

TransformCoordinates[{x_,y_},opts:OptionsPattern[]]:=FixedPoint[ReplicateRegion,{Re[#],Im[#]}&@((OptionValue[Zoom] E^(I OptionValue[Rotation])) (OptionValue[XShift]+I OptionValue[YShift]+x+I y)^(OptionValue[CopiesPerRotation]+OptionValue[Spirals] I Log[$DrosteScale]/(2 \[DoubledPi]))),OptionValue[MaxRecursion]];
Options[TransformCoordinates]={Zoom->1,XShift->0,YShift->0,Rotation->0,CopiesPerRotation->1,MaxRecursion->10,Spirals->1};

ReplicateRegion[{x_,y_}]:=Which[(*If outside the image area,move closer*)Abs[x]>1||Abs[y]>$AspectRatio,{x,y} $DrosteScale,(*If inside the frame move out to the main image*)Abs[x]<$DrosteScale&&Abs[y]<$DrosteScale $AspectRatio,{x,y}/$DrosteScale,(*otherwise use the calculated coordinates*)True,{x,y}];

Комментарии (6)


  1. opckSheff
    01.07.2015 16:37
    +3

    Мозг рвет по полной программе. Математика и Mathematica — это нечто, конечно. Вроде и проходил ТФКП в универе, но поверхностно, и в такие моменты понимаю, что нихрена не понимаю.


  1. jamepock
    02.07.2015 11:45

    теорию функций ко?мплексного переменного — теорию функций комплексной переменной поправьте, так будет вернее)


    1. OsipovRoman Автор
      02.07.2015 12:25

      Нет, «теория функций комплексного переменного», устоявшееся название.

      image


  1. brainick
    02.07.2015 12:54
    -3

    Очередная пальба из пушки по воробьям, бессмысленная и беспощадная.

    Первый сертифицированный инструктор технологий и учебных курсов компании Wolfram Research на территории Восточной Европы уже опубликовал с десяток таких переводов — и большая часть всякая чушь — то обработка комиксов, то гора математического мусора в статье об арбелосе.

    Зато когда нашего местного гуру Mathematiсi когда просят показать реальный проект, в котором можно увидеть реальную силу Mathematiсi он сразу же сливается. (как я и предсказывал)

    Цитирую отсюда habrahabr.ru/company/wolfram/blog/261413

    >>Все кто могут их писать, и я в том числе, в настоящий момент очень сильно заняты в своих проектах. Написание одного поста отнимает очень много времени, а писать отписки на пол странички не хочется.<<


  1. zverok
    02.07.2015 13:55

    Мы можем представить наши пары координат {x,y} как части комплексного числа p=x+I

    В последней формуле y протерялся. У автора
    If instead of thinking about our coordinates as an {x,y} pair, we think of them as a complex p=x+I y


    1. OsipovRoman Автор
      02.07.2015 14:08

      Спасибо, исправил.