ПОМОГИТЕ!!!! Поворот изображения!!! , Необходим алгоритм, а лучше исходник!!!

Пожалуйста!!! Кто-нибудь помогите!!! У меня есть 2-мерная матрица 1500х1500. В ней хранятся значения пикселей изображения (RGB). Эту матрицу надо повернуть на любой угол. Проблема в том, что повернуть нажо качественно, как в Photoshop!!! Обычный поворот не подходит!!! Нужна интерполяция или что-то ещё. Программ в c++. Нельзя использовать ни DirecX, ничего другого!!! Заранее огромное спасибо!!!

см. alglib.sources.ru раздел "интерполяция", там ищи бикубическую интерполяцию. это то, что нужно.

можно ещё подумать над тем, как представлять данные (в RGB или в системе цвет-яркость-оттенок - разные варианты по разному могут вести себя при интерполяции, цвета могут искажаться в большей или меньшей степени, в зависимости от того, что именно интерполируется).

Я там уже был :) Все исходники непонятно, что делают!!! Я конечно не гений, но и не чайник. Однако разобраться не могу...
В матиматике я дуб дубом! Мне нужен либо простая процедура поворота матрицы с интерполяцией, либо конкретная формула, которая использует:
xc,yc - центр матрицы (Точка поворота)
x0,y0 - значения повёрнутых координат
x1.y1 - текущие координаты цикла (for i=0; i<=width...)
phi - угол
Надеюсь на помощь!!!

Поворот в плоскости на угол A описывается матрицей

cos A  -sin A
sin A   cos A

т.е., чтобы повернуть точку (x,y) в плоскости, умножаем ветор (x,y)^T на матрицу:

(x') = ( cos A  -sin A  ) * (x)
(y')   ( sin A    cos A )   (y)

т.е.

x' = x * cos A - y * sin A
y' = x * sin A + y * cos A

В твоём случае, поворачиваемык координаты, это (x₁-x_c, y₁-y_c).
В виде простой интерполяции можно при повороте накапливать значения цветовых компонент R,G,B во всех пикселях, куда попадает поворачиваемый + счётчик, затем разделить их все на счётчик. На Pascal примерно вот так:

{ Rotate }
  phi := Random*Pi/2-Pi/4;
  SinCos(phi, s,c);
  xc := w div 2;
  yc := h div 2;
  for i := 0 to w-1 do
    for j := 0 to h-1 do begin
      x := i-xc;
      y := j-yc;
      rx := x*c - y*s + xc;
      ry := x*s + y*c + yc;
      l := Min(Max(floor(rx),0),w-1);
      r := Min(Max(ceil(rx),0),w-1);
      t := Min(Max(floor(ry),0),h-1);
      b := Min(Max(ceil(ry),0),h-1);
      with rm[l,t] do begin
        Inc(r, m[i,j].r);
        Inc(g, m[i,j].g);
        Inc(b, m[i,j].b);
        Inc(c);
      end;
      with rm[l,b] do begin
        Inc(r, m[i,j].r);
        Inc(g, m[i,j].g);
        Inc(b, m[i,j].b);
        Inc(c);
      end;
      with rm[r,t] do begin
        Inc(r, m[i,j].r);
        Inc(g, m[i,j].g);
        Inc(b, m[i,j].b);
        Inc(c);
      end;
      with rm[r,b] do begin
        Inc(r, m[i,j].r);
        Inc(g, m[i,j].g);
        Inc(b, m[i,j].b);
        Inc(c);
      end;
    end;
{ Calc Average }
  for i := 0 to w-1 do
    for j := 0 to h-1 do
      with rm[i,j] do if c > 0 then begin
        r := r div c;
        g := g div c;
        b := b div c;
      end;

А что за rm? Ты не мог бы написать на c++? Очень прошу!!!
И где вторая матрица, в которую заносится повёрнутое изображение? И сам мехонизм?

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

Аппроксимацию??? Неужели ни у кого нет кого-нибудь исходника на эту тему???

Цитата lupash @ 22.09.07, 23:10

Аппроксимацию??? Неужели ни у кого нет кого-нибудь исходника на эту тему???

Тут собственно и писать нечего, есть 2 формулы преобразования системы координат (приведены выше), далее организуешь циклы по X,Y, расчитываешь для каждой точки их новые координаты.
Единственный момент - задачу нужно решать с конца, если ты будешь двигаться по исходному изображению, вычисляя место где должна находится данная точка, то на новом изображении из-за округления появятся дыры - в некоторые точки мы вообще не попадем, а в некоторые попадем несколько раз. Чтоб этого не случилось - ты должен пройтись по всем точкам нового изображения и для каждой из них взять цвет из исходного, тогда изображение получается целостным.
Грубо это можно записать так:
1. Вычисляем по формуле поворота координаты 4-х углов после поворота, необходимо для определения новых размеров области, так после поворота вашего массива (1500х1500)на 45 градусов, его высота и ширина станут соответственно 2121х2121.
Формула поворота для этого случая запишется как:

XDest = x*cos(alpha)-y*sin(alpha)
YDest = x*sin(alpha)+y*cos(alpha)

2. Пройдясь по новой области (от 0 до 2120) вычисляем где была наша точка на исходном изображении, тут нужно уже применить обратную формулу:

Xsource = x*cos(alpha)+y*sin(alpha)
Ysource = -x*sin(alpha)+y*cos(alpha)

Далее просто посмотреть попадает ли точка с координатами (Xsource, Ysource) в вашу начальную область, если не попадает - значит рисуете это место на новой картинке черным цветом, там изображения нет, если попадает - рисуете соответствующим цветом из массива.
Пример кода, массив Source содержит исходное изображение, Dest - повернутое:

  Var Source,Dest:array of array of TRGB;
      x,y,i,j:integer;
      x1,y1,x2,y2,x3,y3,x4,y4:integer; //координаты углов, нужны для определения нового размера
      Xa,Ya,Xb,Yb:integer; //окно задающее исходное изображение, для вашего случая они равны 0,0,1500,1500
      Xleft,Ytop:integer;  //верхний левый угол нового изображения
      Xc,Yc:integer; //координаты оси поворота
      Xsrc,Ysrc, Width,Height:integer;
      s,c:double;
  Begin
      s:=sin(pi/180*alpha); c:=cos(pi/180*alpha); //для простоты записи, ну и скорость увеличится
      //Вычисляем центр (может быть и относительно любой другой точки)
      Xc:=(Xb-Xa)/2; Yc:=(Xb-Xa)/2;
      //Поочередно вычисляем новые координаты для всех четырех углов, не забываем сместить ось вращения в центр
      x:=Xa - Xc; y:= Ya - Yc;
      X1 =  x*cos(alpha)-y*sin(alpha);
      Y1 =  x*sin(alpha)+y*cos(alpha);
      //аналогично и для всех остальных углов    
      x:=Xa - Xc; y:= Yb - Yc;
      X2 =  x*cos(alpha)-y*sin(alpha);
      Y2 =  x*sin(alpha)+y*cos(alpha);
      //ну и так далее
........................      
      //имея координаты всех 4-х углов можем расчитать новый размер области,
      //для этого из максимальной координаты вычтем минимальную, для сокращения кода
      //предположим что у нас уже есть 2 функции находящии минимальное и максимальное
      //значение из 4-х, MinOf4 и MaxOf4
      width:=maxof4(x1,x2,x3,x4)-min(x1,x2,x3,x4);
      height:=maxof4(y1,y2,y3,y4)-min(y1,y2,y3,y4);
      //Нам также понадобятся координаты верхнего левого угла:
      Xleft:=min(x1,x2,x3,x4); Ytop:=min(y1,y2,y3,y4);
      //Ну и последний этап - пройтись по всем точкам новой области и нарисовать там повернутое изображение:
      for i:=0 to height do
      for j:=0 to width do begin
        x:=j+Xleft; y:=i+Yleft;
        Xsrc = x*c+y*s+Xc;  Ysrc = -x*s+y*c+Yc;
        //Смотрим попали ли мы на исходное изображение
        IF (Xsrc < Xa) or (Xsrc > Xb) or (Ysrc < Ya) or (Ysrc > Yb)
        Then Dest[i,j]:=0
        Else Dest[i,j]:=Source[Xsrc,Ysrc];
      end;
  end;

Тут я умышленно убрал округление, создание массивов, изменение их размеров, создание процедур (минимума и поворота), дабы не усложнять код лишними командами (если человек не знает Делфи).

Сообщение отредактировано: a_blade - 23.09.07, 06:14

Не работает!!! И первый и второй алгоритмы работаю как мой: поворачивают, но с искажениями! Вот скрины:
Мой поворот:
[img]http://up.li.ru/?id=314264;Comp+1+%280-00-00-01%29.jpg[/img]

Поворот Adobe:
[img]http://up.li.ru/?id=314265;Comp+1+%280-00-00-01%29_1.jpg[/img]


СПАСАЙТЕ МЕНЯ!!!!

lupash
Зная координаты, на старом изоброжении. Новое значение получается интерполированием.
http://en.wikipedia.org/wiki/Interpolation
Твой алгоритм, выбирает наиближайшее значение. И соответствует Nearest interpolation.
Тебе нужно билинейная интерполяция или бикубическая. Алгоритмы, найдешь ввики.

Сообщение отредактировано: Pavia - 23.09.07, 20:41

Цитата lupash @ 23.09.07, 15:07

И первый и второй алгоритмы работаю как мой: поворачивают, но с искажениями

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

  alpha := ???;
  xa := 0;
  ya := 0;
  xb := w-1;
  yb := h-1;
  s := sin(alpha);
  c := cos(alpha);
  xc := (xb-xa) div 2;
  yc := (xb-xa) div 2;
  x := xa-xc;
  y := ya-yc;
  x1 := round(x*cos(alpha) - y*sin(alpha));
  y1 := round(x*sin(alpha) + y*cos(alpha));
  x := xa-xc;
  y := yb-yc;
  x2 := round(x*cos(alpha) - y*sin(alpha));
  y2 := round(x*sin(alpha) + y*cos(alpha));
  x := xb-xc;
  y := ya-yc;
  x3 := round(x*cos(alpha) - y*sin(alpha));
  y3 := round(x*sin(alpha) + y*cos(alpha));
  x := xb-xc;
  y := yb-yc;
  x4 := round(x*cos(alpha) - y*sin(alpha));
  y4 := round(x*sin(alpha) + y*cos(alpha));
  dl := Min4(x1,x2,x3,x4);
  dt := Min4(y1,y2,y3,y4);
  dw := Max4(x1,x2,x3,x4) - dl;
  dh := Max4(y1,y2,y3,y4) - dt;
  for i := 0 to dh-1 do
    for j := 0 to dw-1 do begin
      x := j + dl;
      y := i + dt;
      sx :=  x*c + y*s + xc;
      sy := -x*s + y*c + yc;
      if (sx<xa)or(sx>xb)or(sy<ya)or(sy>yb) then begin
        dest[j,i].r := 0;
        dest[j,i].g := 0;
        dest[j,i].b := 0;
      end else begin
        sl := Min(Max(floor(sx),0),w-1);
        sr := Min(Max(sl+1,0),w-1);
        st := Min(Max(floor(sy),0),h-1);
        sb := Min(Max(st+1,0),h-1);
        klt := 1/(sqrt(sqr(sl-sx)+sqr(st-sy))+epsilon);
        klb := 1/(sqrt(sqr(sl-sx)+sqr(sb-sy))+epsilon);
        krt := 1/(sqrt(sqr(sr-sx)+sqr(st-sy))+epsilon);
        krb := 1/(sqrt(sqr(sr-sx)+sqr(sb-sy))+epsilon);
        ksum := 1/(klt + klb + krt + krb);
        klt := klt*ksum;
        klb := klb*ksum;
        krt := krt*ksum;
        krb := krb*ksum;
        with dest[j,i] do begin
          r := r + klt*source[sl,st].r
                 + klb*source[sl,sb].r
                 + krt*source[sr,st].r
                 + krb*source[sr,sb].r;
          g := g + klt*source[sl,st].g
                 + klb*source[sl,sb].g
                 + krt*source[sr,st].g
                 + krb*source[sr,sb].g;
          b := b + klt*source[sl,st].b
                 + klb*source[sl,sb].b
                 + krt*source[sr,st].b
                 + krb*source[sr,sb].b;
        end;
      end;
    end;

epsilon = 0.001 используется, чтобы не делить на 0.

Сейчас попробую... И на будущее, если можно пишите пожалуйста на c++, а то я запарился всё переписывать... Конечно если можно

Добавлено 23.09.07, 22:35
ОГРОМНОЕ СПАСИБО!!! ors_archangel - МЕГА ЧУВАК!!! Поклон тебе и уважуха!!!