(************** Content-type: application/mathematica **************
                     CreatedBy='Mathematica 5.1'

                    Mathematica-Compatible Notebook

This notebook can be used with any Mathematica-compatible
application, such as Mathematica, MathReader or Publicon. The data
for the notebook starts with the line containing stars above.

To get the notebook into a Mathematica-compatible application, do
one of the following:

* Save the data starting with the line of stars above into a file
  with a name ending in .nb, then open the file inside the
  application;

* Copy the data starting with the line of stars above to the
  clipboard, then use the Paste menu command inside the application.

Data for notebooks contains only printable 7-bit ASCII and can be
sent directly in email or through ftp in text mode.  Newlines can be
CR, LF or CRLF (Unix, Macintosh or MS-DOS style).

NOTE: If you modify the data for this notebook not in a Mathematica-
compatible application, you must delete the line below containing
the word CacheID, otherwise Mathematica-compatible applications may
try to use invalid cache data.

For more information on notebooks and Mathematica-compatible 
applications, contact Wolfram Research:
  web: http://www.wolfram.com
  email: info@wolfram.com
  phone: +1-217-398-0700 (U.S.)

Notebook reader applications are available free of charge from 
Wolfram Research.
*******************************************************************)

(*CacheID: 232*)


(*NotebookFileLineBreakTest
NotebookFileLineBreakTest*)
(*NotebookOptionsPosition[     42647,       1193]*)
(*NotebookOutlinePosition[     44123,       1236]*)
(*  CellTagsIndexPosition[     44079,       1232]*)
(*WindowFrame->Normal*)



Notebook[{
Cell["Transformada coseno discreta. Ampliaci\[OAcute]n de im\[AAcute]genes", \
"Title"],

Cell["\<\
Francisco Javier P\[EAcute]rez Gonz\[AAcute]lez
Departamento de An\[AAcute]lisis Matem\[AAcute]tico
Universidad de Granada\
\>", "Author"],

Cell[CellGroupData[{

Cell["Inicializaci\[OAcute]n", "Section"],

Cell[BoxData[{
    \(<< Graphics`\), "\n", 
    \(<< LinearAlgebra`MatrixManipulation`\), "\n", 
    \(\(SetOptions[ListDensityPlot, Mesh \[Rule] False, 
        Frame \[Rule] False];\)\)}], "Input",
  InitializationCell->True]
}, Closed]],

Cell[CellGroupData[{

Cell["Transformada coseno discreta", "Section",
  CounterAssignments->{{"Section", 3}, {"Subsection", 0}, {
      "Subsubsection", 0}, {"NumberedEquation", 0}}],

Cell["\<\
La transformada coseno discreta (DCT) es una herramienta muy \
\[UAcute]til en el procesamiento de se\[NTilde]ales bidimensionales y es \
usada en los est\[AAcute]ndares de compresi\[OAcute]n de imagen JPEG (Joint \
Photographic Experts Group)  y de v\[IAcute]deo MPEG (Motion Picture Experts \
Group).  Como veremos, la DCT puede calcularse mediante la transformada de \
Fourier discreta (DFT) y, por tanto, la DCT puede calcularse con el algoritmo \
FFT (Fast Fourier Transform) con todas las ventajas que esto supone. \
\>", \
"Text"],

Cell[CellGroupData[{

Cell["1D - DCT", "Subsubsection"],

Cell[TextData[{
  "La transformada coseno discreta en una variable (1D - DCT)",
  StyleBox[" ",
    FontSlant->"Italic"],
  "de un vector ",
  Cell[BoxData[
      \(TraditionalForm\`x = \((x[0], x[1],  ... , 
            x[N - 1])\) \[Element] \[DoubleStruckCapitalC]\^N\)]],
  "es el vector ",
  Cell[BoxData[
      \(TraditionalForm\`X = \((X[0], X[1],  ... , 
            X[N - 1])\) \[Element] \[DoubleStruckCapitalC]\^N\)]],
  " definido por: "
}], "Text"],

Cell[BoxData[
    FormBox[
      FrameBox[
        RowBox[{
          
          FormBox[\(X[
                n]\  = \ \ \[Sum]\+\(q = 0\)\%\(N - 1\)\(\(\[Alpha][n, 
                    N]\)\(\ \)\(x[
                    q]\)\(\ \ \)\(Cos[\(\((2  q + 1)\)\ n\ \[Pi]\)\/\(2\ N\)]\
\)\(\ \)\)\),
            "TraditionalForm"], ",", "\t", \(n\  = \ 0\), ",", "1", ",", "...",
           ",", \(N - 1\)}]], TraditionalForm]], "NumberedEquation",
  TextAlignment->Center],

Cell["\<\
Seguimos la costumbre establecida de usar \[IAcute]ndices que \
empiezan en 0 lo que obliga a modificar convenientemente las expresiones para \
definir comandos. Por ejemplo, si queremos usar la definici\[OAcute]n \
anterior para definir un comando, deberemos escribirla en la forma \
\>", \
"Text",
  TextAlignment->Left],

Cell[BoxData[
    FormBox[
      FrameBox[
        RowBox[{
          
          FormBox[\(X[\([n]\)]\  = \ \ \[Sum]\+\(q = 1\)\%N\(\( \[Alpha][
                    n - 1, N]\)\(\ \)\(x[\([q]\)]\)\(\ \ \)\(Cos[\(\((2  q - \
1)\) \((\ n - 1)\)\ \[Pi]\)\/\(2\ N\)]\)\(\ \)\)\),
            "TraditionalForm"], ",", "\t", \(n\  = \ 1\), ",", "2", ",", "...",
           ",", "N"}]], TraditionalForm]], "NumberedEquation",
  TextAlignment->Center],

Cell["La transformada inversa 1D - DCT viene dada por:", "Text"],

Cell[BoxData[
    FormBox[
      FrameBox[
        RowBox[{
          
          FormBox[\(x[
                n]\  = \ \[Sum]\+\(q = 0\)\%\(N - 1\)\[Alpha][q, N]\ X[
                    q]\ \ Cos[\(\(\ \)\(\((2  n + 1)\)\ q\ \[Pi]\)\)\/\(2\ \
N\)]\),
            "TraditionalForm"], ",", "\t", \(n\  = \ 0\), ",", "1", ",", "...",
           ",", \(N - 1\)}]], TraditionalForm]], "NumberedEquation",
  TextAlignment->Center],

Cell[BoxData[
    FormBox[
      RowBox[{"donde", 
        FormBox[\(\[Alpha][\ , ]\),
          "TraditionalForm"], "es", " ", "una", " ", "funci\[OAcute]n", "  ", 
        "definida", " ", "para", " ", "todo", " ", Cell[TextData[Cell[
        BoxData[
            \(TraditionalForm\`P \[Element] \[DoubleStruckCapitalN]\)]]]], 
        StyleBox[
          RowBox[{" ", 
            StyleBox[" ",
              FontSlant->"Plain"]}]], 
        StyleBox["y",
          FontSlant->"Plain"], 
        StyleBox[" ",
          FontSlant->"Plain"], 
        StyleBox["para",
          FontSlant->"Plain"], 
        StyleBox[" ",
          FontSlant->"Plain"], 
        StyleBox["todo",
          FontSlant->"Plain"], 
        StyleBox[" ",
          FontSlant->"Plain"], Cell[TextData[Cell[BoxData[
            \(TraditionalForm\`j \[Element] \[DoubleStruckCapitalN]\)]]]], 
        StyleBox["  ",
          FontSlant->"Plain"], 
        StyleBox["con",
          FontSlant->"Plain"], 
        StyleBox["  ",
          FontSlant->"Plain"], Cell[TextData[Cell[BoxData[
            \(TraditionalForm\`0 \[LessEqual] j \[LessEqual] P - 1\)]]]], 
        "por"}], TraditionalForm]], "Text"],

Cell[BoxData[
    FormBox[
      FrameBox[
        RowBox[{
          
          FormBox[\(\[Alpha][0, P] = 
              1\/\@P, \ \ \[Alpha][j, P] = \@\(2\/P\)\),
            "TraditionalForm"], ",", 
          RowBox[{"para", "   ", 
            RowBox[{
              FormBox[\(j = 1, 2,  ... , P - 1\),
                "TraditionalForm"], "."}]}]}]], 
      TraditionalForm]], "NumberedEquation",
  TextAlignment->Center],

Cell[TextData[{
  "Es f\[AAcute]cil definir un comando que calcule simb\[OAcute]licamente \
dicha transformada. Observa que ",
  Cell[BoxData[
      \(TraditionalForm\`X[n]\)]],
  " se obtiene haciendo el producto escalar en ",
  Cell[BoxData[
      \(TraditionalForm\`\[DoubleStruckCapitalC]\^N\)]],
  "del vector",
  Cell[BoxData[
      \(TraditionalForm\`\(\(\ \)\(x = \((x[0], x[1],  ... , 
          x[N - 1])\)\)\)\)]],
  " por el vector "
}], "Text"],

Cell[BoxData[
    \(TraditionalForm\`u[n, 
        N] = \[Alpha][n, N] \((Cos[\(\(\ \)\(n\ \[Pi]\)\)\/\(2\ N\)], 
          Cos[\(3\ n\ \[Pi]\)\/\(2\ N\)],  ... , 
          Cos[\(\(\ \)\(\((2  N - 1)\)\ n\ \[Pi]\)\)\/\(2\ N\)]\ )\)\)], \
"NumberedEquation",
  TextAlignment->Center],

Cell["Definamos estos vectores.", "Text"],

Cell[BoxData[{
    \(\(\[Alpha][0, 
          M_] := \(\[Alpha][0, M] = \@\(1/M\)\);\)\), "\[IndentingNewLine]", 
    \(\(\[Alpha][k_, 
          M_] := \(\[Alpha][k, M] = \@\(2/M\)\);\)\), "\[IndentingNewLine]", 
    \(u[k_, 
        M_] := \(u[k, M] = \[Alpha][k, M] 
          Cos[\(\((2\ Range[0, M - 1] + 1)\) k\ \[Pi]\)\/\(2\ M\)]\)\)}], \
"Input"],

Cell[TextData[{
  "Podemos comprobar que dichos vectores forman una base ortonormal de ",
  Cell[BoxData[
      \(TraditionalForm\`\[DoubleStruckCapitalC]\^N\)]],
  "."
}], "Text"],

Cell[BoxData[
    \(\((M = Table[u[k, 6], {k, 0, 5}])\) // MatrixForm\)], "Input"],

Cell[BoxData[{
    \(\(M . Transpose[M] // MatrixForm\) // 
      Simplify\), "\[IndentingNewLine]", 
    \(Clear[M]\)}], "Input"],

Cell[TextData[{
  "Algebraicamente, la DCT no es m\[AAcute]s que la representaci\[OAcute]n \
del vector ",
  Cell[BoxData[
      \(TraditionalForm\`x\)]],
  " en la base de ",
  Cell[BoxData[
      \(TraditionalForm\`\[DoubleStruckCapitalC]\^N\)]],
  " formada por los vectores ",
  Cell[BoxData[
      \(TraditionalForm\`u[k, N]\)]],
  ", ",
  Cell[BoxData[
      \(TraditionalForm\`k = 0, 1, 2,  ... , N - 1\)]],
  ". Como dicha base es ortonormal se verifica que la DCT conserva la \
distancia eucl\[IAcute]dea o, dicho de otra forma, la DCT conserva la energ\
\[IAcute]a de la se\[NTilde]al. \nLos siguientes comandos calculan la DCT de \
forma simb\[OAcute]lica y de forma num\[EAcute]rica. En el comando num\
\[EAcute]rico, en vez de usar Table[ ] usamos Range[ ] que es m\[AAcute]s r\
\[AAcute]pido; para lo que asignamos a la funci\[OAcute]n u[ ]  el atributo \
\"Listable\". Adem\[AAcute]s, para evitar que ",
  StyleBox["Mathematica",
    FontSlant->"Italic"],
  " almacene largos vectores simb\[OAcute]licos, hacemos todos los \
c\[AAcute]lculos num\[EAcute]ricamente."
}], "Text"],

Cell[BoxData[{
    \(SetAttributes[u, Listable]\), "\[IndentingNewLine]", 
    \(u[Range[0, 5], 6] // MatrixForm\)}], "Input"],

Cell[BoxData[{
    \(sbDCT[
        x_List] := \(\(u[Range[0, Length[x] - 1], Length[x]] . 
          x\)\(//\)\(Simplify\)\( (*\ 
          calcula\ la\ DCT\ de\ forma\ simb\[OAcute]lica\ *) \)\)\), "\
\[IndentingNewLine]", 
    \(numDCT[x_List] := 
      u[Range[0, Length[x] - 1], N[Length[x]]] . x // 
        Chop\), "\[IndentingNewLine]", 
    \( (*\ calcula\ la\ DCT\ de\ forma\ num\[EAcute]rica\ *) \)}], "Input"],

Cell[BoxData[{
    \(Clear[a, b, c, d, e, f]\), "\[IndentingNewLine]", 
    \(sbDCT[{a, b, c, d, e, f}]\)}], "Input"],

Cell[BoxData[
    \(numDCT[{1, 2, 3, 4, 5, 6}]\)], "Input"],

Cell["\<\
Los siguientes comandos calculan la inversa de la DCT de forma simb\
\[OAcute]lica y num\[EAcute]rica. Para hacer m\[AAcute]s r\[AAcute]pido este \
comando damos el atributo \"Listable\" a la funci\[OAcute]n \[Alpha] y usamos \
Range[ ]. Observa que para que la funci\[OAcute]n \[Alpha] trabaje num\
\[EAcute]ricamente su segundo argumento debe ser un valor aproximado.\
\>", \
"Text"],

Cell[BoxData[{
    \(SetAttributes[\[Alpha], Listable]\), "\[IndentingNewLine]", 
    \(\[Alpha][Range[0, 3], 4. ]\)}], "Input"],

Cell[BoxData[{
    \(\(\(inversesbDCT[x_] := 
        Table[\[Alpha][Range[0, Length[x] - 1], 
                Length[x]] . \((Cos[\(\((2\ j + 1)\) Range[0, Length[x] - 1]\ \
\[Pi]\)\/\(2\ Length[x]\)]*x)\), {j, 0, Length[x] - 1}] // 
          Simplify;\)\(\[IndentingNewLine]\)
    \) (*\ 
      calcula\ la\ inversa\ DCT\ de\ forma\ simb\[OAcute]lica\ *) \), "\
\[IndentingNewLine]", 
    \(\(inversenumDCT[
          x_] := \(Table[\[Alpha][Range[0, Length[x] - 1], 
                  N[Length[
                      x]]] . \((Cos[\(\((2\ j + 1)\) Range[0, Length[x] - 1]\ \
N[\[Pi]]\)\/\(2\ Length[x]\)]*x)\), {j, 0, Length[x] - 1}] // N\) // 
          Chop;\)\), "\[IndentingNewLine]", 
    \( (*\ 
      calcula\ la\ inversa\ DCT\ de\ forma\ num\[EAcute]rica\ *) \)}], "Input"],

Cell[BoxData[{
    \(sbDCT[{a, b, c, d, e, f}]\), "\[IndentingNewLine]", 
    \(inversesbDCT[%]\)}], "Input"],

Cell[BoxData[{
    \(numDCT[Range[6]]\), "\[IndentingNewLine]", 
    \(inversenumDCT[%]\)}], "Input"],

Cell["Para calcular la DCT usando la DFT se procede como sigue.", "Text"],

Cell[BoxData[
    \(\(\(\(x = {1, 2, 3, 4}\  (*\ 
        este\ es\ el\ vector\ cuya\ DCT\ se\ quiere\ calcular\ *) ;\)\
\[IndentingNewLine]
    y = Join[x, Reverse[x]]\ \  (*\ 
      se\ a\[NTilde]ade\ al\ final\ de\ x\ su\ revertido\ \
*) \[IndentingNewLine]
    z = Fourier[y, FourierParameters \[Rule] {1, \(-1\)}]\  (*\ 
      calculamos\ la\ DFT\ del\ vector\ y\ *) \[IndentingNewLine]
    Take[z, Length[x]]*
        Table[\[Alpha][k, Length[x]]\/2*
            Exp[\(-\(\(I\ Pi\ k\)\/\(2  Length[x]\)\)\)], {k, 0, 
            Length[x] - 1}]\  // Chop (*\ 
      nos\ quedamos\ con\ la\ primera\ mitad\ de\ z\ y\ multiplicamos\ por\ \
coeficientes\ apropiados\ *) \[IndentingNewLine]
    % \[Equal] numDCT[x]\)\(\ \)\( (*\ 
      comprobamos\ el\ resultado\ obtenido\ *) \)\)\)], "Input"],

Cell["\<\
Podemos definir un comando que haga estas operaciones el cual \
solamente trabajar\[AAcute] num\[EAcute]ricamente. \
\>", "Text"],

Cell[BoxData[
    \(fastDCT[list_List] := 
      Module[{temp}, temp = Join[list, Reverse[list]]; 
        Take[Fourier[temp, FourierParameters \[Rule] {1, \(-1\)}], 
              Length[list]]*\((\((\[Alpha][#, N[Length[list]]]\/2*
                      Exp[\(-\(\(\(I\)\(\ \)\(N[Pi]\)\(\ \)\)\/\(2 
                                   Length[
                                    list]\)\)\) #] &)\)\[IndentingNewLine][
                Range[0, Length[list] - 1]])\) // Chop]\)], "Input"],

Cell["\<\
Comprobamos la coincidencia de este comando con el antes \
definido.\
\>", "Text"],

Cell[BoxData[{
    \(Clear[a]\), "\[IndentingNewLine]", 
    \(\(a = Table[Random[], {200}];\)\), "\[IndentingNewLine]", 
    \(Simplify[\((\((numDCT[a] - fastDCT[a])\) // Chop)\) \[Equal] 
        0]\)}], "Input"],

Cell["\<\
Naturalmente, el comando fstDCT[ ] es much\[IAcute]simo m\[AAcute]s \
r\[AAcute]pido que numDCT[ ].\
\>", "Text"],

Cell[BoxData[
    \(\(\(Timing[\(numDCT[
          Table[Random[], {4096}]];\)] (*\ \[DownExclamation]tiene\ que\ \
calcular\ una\ matriz\ con\ 4096\[Times]4096 = 
        16777216\ de\ \(entradas!\)\ *) \[IndentingNewLine]
    Timing[\(fastDCT[Table[Random[], {4096}]];\)]\)\( (*\ 
      usa\ FFT\ *) \)\( (*\ 
      la\ fastDCT\ es\ mucho\ m\[AAcute]s\ r\[AAcute]pida\ que\ la\ numDCT\ *) \
\)\)\)], "Input"],

Cell["El siguiente comando calcula la inversa de la fastDCT[ ].", "Text"],

Cell[BoxData[
    \(inversefastDCT[list_] := 
      Module[{temp}, 
        temp = \((\((2\/\[Alpha][#, N[Length[list]]]*
                      Exp[\(\(I\ N[Pi]\)\/\(2 
                                 Length[list]\)\) #] &)\)\[IndentingNewLine][
                Range[0, Length[list] - 1]])\)*list; 
        Take[InverseFourier[
              Join[Append[temp, 0], Reverse[Conjugate[Rest[temp]]]], 
              FourierParameters \[Rule] {1, \(-1\)}], Length[list]] // 
          Chop]\)], "Input"],

Cell[BoxData[
    \(inversefastDCT[fastDCT[Range[100]]]\)], "Input"],

Cell[BoxData[
    \(inversefastDCT[
      numDCT[{1, 2, 3, 4, 5, 6, \(-1\),  .25, \(-4.5\), 8, \(- .0003\), 
          4.657, 3.1416, N[Pi]}]]\)], "Input"],

Cell[TextData[{
  "Aqu\[IAcute] est\[AAcute]n representados los vectores b\[AAcute]sicos para \
",
  Cell[BoxData[
      \(TraditionalForm\`N = 16\)]],
  "."
}], "Text",
  GridBoxOptions->{RowLines->False}],

Cell[BoxData[{
    \(makeplotlist[x_List] := 
      Line /@ Table[{{k, 0}, {k, x[\([\)\(k\)\(]\)]}}, {k, 1, 
            Length[x]}]\), "\[IndentingNewLine]", 
    \(\(Show[
        GraphicsArray[
          Partition[\((ListPlot[#1\ , PlotStyle \[Rule] PointSize[0.04], 
                    Ticks \[Rule] None, Epilog \[Rule] makeplotlist[#1], 
                    PlotRange \[Rule] {\(- .5\),  .5}, 
                    AxesOrigin \[Rule] {1, 0}, 
                    DisplayFunction \[Rule] Identity] &)\) /@ \ 
              Table[u[k, 16. ], {k, 0, 15}], 4]], 
        ImageSize \[Rule] {426.625, 263.5}];\)\)}], "Input"],

Cell[TextData[{
  "Observa que los vectores b\[AAcute]sicos son muestras de funciones \
cosenoidales del tipo ",
  Cell[BoxData[
      \(TraditionalForm\`Cos[2  k\ \[Pi]\ t + k\ \[Pi]/2  N]\)]],
  " tomadas en los momentos ",
  Cell[BoxData[
      \(TraditionalForm\`t\_n = n/2  N\)]],
  " por lo que la frecuencia de muestreo es ",
  Cell[BoxData[
      \(TraditionalForm\`2  N\)]],
  ". Podemos comprobarlo gr\[AAcute]ficamente como sigue."
}], "Text"],

Cell[BoxData[
    RowBox[{
      RowBox[{"Show", "[", 
        RowBox[{
          RowBox[{"GraphicsArray", "[", 
            RowBox[{"Partition", "[", 
              RowBox[{
                
                RowBox[{\((Plot[#, {t, 0,  .5 - 1/32}, 
                        PlotRange \[Rule] {\(- .5\),  .5}, 
                        Ticks \[Rule] None, 
                        Epilog \[Rule] {PointSize[0.04], 
                            Point /@ 
                              Table[{j, # /. t \[Rule] j}, {j, 0, 15/32, 
                                  1/32}]}, 
                        DisplayFunction \[Rule] Identity] &)\), "/@", 
                  RowBox[{"Table", "[", 
                    RowBox[{
                      RowBox[{\(\[Alpha][k, 16. ]\), "*", " ", 
                        
                        FormBox[\(Cos[2  k\ \[Pi]\ t + k\ \[Pi]/\((2*16)\)]\),
                          
                          "TraditionalForm"]}], ",", \({k, 0, 15}\)}], 
                    "]"}]}], ",", "4"}], "]"}], "]"}], 
          ",", \(ImageSize \[Rule] {426.625, 263.5}\)}], "]"}], 
      ";"}]], "Input"]
}, Open  ]],

Cell[CellGroupData[{

Cell["2D - DCT", "Subsubsection"],

Cell[TextData[{
  "La transformada coseno discreta en dos variables (2D - DCT)",
  StyleBox[" ",
    FontSlant->"Italic"],
  "de una matriz ",
  Cell[BoxData[
      \(TraditionalForm\`x = \(\(x[p, q]\)\_\(0 \[LessEqual] p \[LessEqual] M \
- 1, 0 \[LessEqual] q \[LessEqual] N - 1\) \[Element] \[ScriptCapitalM]\_\(M\
\[Times]N\) = \[DoubleStruckCapitalC]\^\(M\[Times]N\)\)\)]],
  " es la matriz ",
  Cell[BoxData[
      \(TraditionalForm\`X = \(X[m, n]\)\_\(0 \[LessEqual] m \[LessEqual] M - \
1, 0 \[LessEqual] n \[LessEqual] N - 1\) \[Element] \[ScriptCapitalM]\_\(M\
\[Times]N\)\)]],
  " definida por:"
}], "Text"],

Cell[BoxData[
    FormBox[
      FrameBox[\(X[m, \ 
            n]\  = \ \[Alpha][m, M]\ \ \[Alpha][n, 
              N]\ \ \(\[Sum]\+\(p = 0\)\%\(M - 1\)\(\[Sum]\+\(q = 0\)\%\(N - \
1\)x[p, \ 
                    q]\ \ Cos[\(\(\ \)\(\((2  p + 1)\)\ \ m\ \[Pi]\)\)\/\(2\ \
M\)]\ \ Cos[\(\(\ \)\(\((2  q + 1)\)\ n\ \[Pi]\)\)\/\(2\ N\)]\)\)\)], 
      TraditionalForm]], "NumberedEquation",
  TextAlignment->Left],

Cell[TextData[{
  "para  ",
  Cell[BoxData[
      \(TraditionalForm\`m\  = \ 0, 1,  ... , M - 1, \ n = 0, 1,  ... , 
      N - 1\)]],
  ". Donde los coeficientes ",
  Cell[BoxData[
      \(TraditionalForm\`\[Alpha][\ , ]\)]],
  " est\[AAcute]n definidos como en (3)."
}], "Text"],

Cell[TextData[{
  "Representando por ",
  Cell[BoxData[
      \(TraditionalForm\`x\_p\)]],
  " la fila p-\[EAcute]sima de la matriz ",
  Cell[BoxData[
      \(TraditionalForm\`x[p, q]\)]],
  " y por ",
  Cell[BoxData[
      \(TraditionalForm\`X\_p[n]\)]],
  " la componente n-\[EAcute]sima de su 1D",
  StyleBox[" - ",
    FontSize->10],
  "DCT, podemos escribir"
}], "Text"],

Cell[BoxData[
    \(TraditionalForm\`X[m, \ 
        n]\  = \ \ \(\[Sum]\+\(p = 0\)\%\(M - 1\)\ \[Alpha][m, M] 
            Cos[\(\(\ \)\(\((2  p + 1)\)\ m\ \[Pi]\)\)\/\(2\ M\)]\ \
\((\[Sum]\+\(q = 0\)\%\(N - 1\)\[Alpha][n, N]\ \ x[p, \ 
                    q]\ \ \ Cos[\(\(\ \)\(\((2  q + 1)\)\ n\ \[Pi]\)\)\/\(2\ \
N\)])\) = \[Sum]\+\(p = 0\)\%\(M - 1\)\ \[Alpha][m, M] 
            Cos[\(\(\ \)\(\((2  p + 1)\)\ m\ \[Pi]\)\)\/\(2\ M\)] 
            X\_p[n]\)\)], "Text"],

Cell[TextData[{
  "Definiendo ahora ",
  Cell[BoxData[
      \(TraditionalForm\`z[n] = \((X\_0[n], \ X\_1[n],  ... , 
          X\_\(M - 1\)[n])\)\)]],
  ", con lo que ",
  Cell[BoxData[
      \(TraditionalForm\`\(z[n]\)[p] = X\_p[n]\)]],
  ", y notando por ",
  Cell[BoxData[
      \(TraditionalForm\`Z[n]\)]],
  " la 1D-DCT de ",
  Cell[BoxData[
      \(TraditionalForm\`z[n]\)]],
  ", tenemos que"
}], "Text"],

Cell[BoxData[
    FormBox[
      RowBox[{\(X[m, \ n]\), " ", "=", 
        RowBox[{
          
          FormBox[\(\[Sum]\+\(p = 0\)\%\(M - 1\)\ \[Alpha][m, M] 
                Cos[\(\(\ \)\(\((2  p + 1)\)\ m\ \[Pi]\)\)\/\(2\ M\)] 
                X\_p[n]\),
            "TraditionalForm"], 
          "=", \(\[Sum]\+\(p = 0\)\%\(M - 1\)\ \[Alpha][m, M] 
                Cos[\(\(\ \)\(\((2  p + 1)\)\ m\ \[Pi]\)\)\/\(2\ M\)] \(z[
                    n]\)[p] = \(Z[n]\)[m]\)}]}], TraditionalForm]], "Text"],

Cell[TextData[{
  StyleBox["Lo que nos indica que la 2D - DCT de una matriz puede calcularse \
con sucesivas 1D - DCT. Primero se calculan las 1D - DCT de las ",
    FontWeight->"Bold"],
  StyleBox["filas",
    FontWeight->"Bold",
    FontSlant->"Italic"],
  StyleBox[" de la matriz y se obtienen ",
    FontWeight->"Bold"],
  Cell[BoxData[
      \(TraditionalForm\`M\)],
    FontWeight->"Bold"],
  StyleBox[" vectores ",
    FontWeight->"Bold"],
  Cell[BoxData[
      \(TraditionalForm\`X\_p \[Element] \[DoubleStruckCapitalC]\^N\)],
    FontWeight->"Bold"],
  StyleBox["; se forma ahora la matriz cuyas ",
    FontWeight->"Bold"],
  StyleBox["filas",
    FontWeight->"Bold",
    FontSlant->"Italic"],
  StyleBox[" son estos vectores y se calcula la 1D - DCT de las ",
    FontWeight->"Bold"],
  StyleBox["columnas",
    FontWeight->"Bold",
    FontSlant->"Italic"],
  StyleBox[" de la misma y as\[IAcute] se obtienen las ",
    FontWeight->"Bold"],
  StyleBox["columnas",
    FontWeight->"Bold",
    FontSlant->"Italic"],
  StyleBox[" de la matriz transformada de la matriz dada.",
    FontWeight->"Bold"],
  " A continuaci\[OAcute]n definimos dos comandos que calculan la 2D - DCT \
transformadad de una matriz. Uno de ellos lo hace de forma directa usando la \
definici\[OAcute]n dada arriba y el otro usa la la transformada \
r\[AAcute]pida 1D - DCT, fastDCT[ ], antes definida."
}], "Text"],

Cell[BoxData[{
    \(Off[General::spell1]\), "\[IndentingNewLine]", 
    \(\(twoDCT[X_?MatrixQ] := \ 
        Module[{M, K}, M = N[\(Dimensions[X]\)[\([1]\)]]; 
          K = N[\(Dimensions[X]\)[\([2]\)]]; 
          Table[\[Alpha][m - 1, M]\ \ \[Alpha][n - 1, 
                  K]\ \ \(\[Sum]\+\(p = 1\)\%M\(\[Sum]\+\(q = 1\)\%K 
                        X[\([p, \ 
                          q]\)]\ \ Cos[\(\(\ \)\(\((2  p - 1)\)\ \ \((m - \
1)\)\ N[\[Pi]]\)\)\/\(2\ M\)]\ \ Cos[\(\(\ \)\(\((2  q - 1)\)\ \((n - 1)\)\ \
N[\[Pi]]\)\)\/\(2\ K\)]\)\) // Chop, {m, 1, M}, {n, 1, 
              K}]];\)\), "\[IndentingNewLine]", 
    \(fast2DCT[X_?MatrixQ] := 
      Map[fastDCT, Transpose[Map[fastDCT, X]]] // 
        Transpose\), "\[IndentingNewLine]", 
    \(On[General::spell1]\)}], "Input"],

Cell["Comprobamos que ambos comandos producen igual resultado.", "Text"],

Cell[BoxData[{
    \(\(a = 
        Table[Random[], {j, 1, 10}, {k, 1, 10}];\)\), "\[IndentingNewLine]", 
    \(Simplify[Chop[twoDCT[a] - fast2DCT[a]] \[Equal] 0]\)}], "Input"],

Cell["\<\
Naturalmente, el comando fast2DCT[ ] es much\[IAcute]simo \
m\[AAcute]s r\[AAcute]pido.\
\>", "Text"],

Cell[BoxData[{
    \(Clear[a, b]\), "\[IndentingNewLine]", 
    \(a = Table[Random[], {j, 1, 16}, {k, 1, 16}];  (*\ 
      una\ matriz\ 16\[Times]16\ *) \[IndentingNewLine]b = 
      Table[Random[], {j, 1, 128}, {k, 1, 128}];  (*\ 
      una\ matriz\ 128\[Times]128\ *) \[IndentingNewLine]Timing[\(twoDCT[
          a];\)]\), "\[IndentingNewLine]", 
    \(Timing[\(fast2DCT[b];\)]\), "\[IndentingNewLine]", 
    \(Clear[a, b]\)}], "Input"],

Cell["Podemos escribir la igualdad (6) en la forma", "Text"],

Cell[BoxData[
    \(TraditionalForm\`X[m, \ 
        n]\  = \ \[Alpha][m, M]\ \ \[Alpha][n, 
          N]\ \ \(\[Sum]\+\(p = 0\)\%\(M - 1\)Cos[\(\(\ \)\(\((2  p + \
1)\)\ \ m\ \[Pi]\)\)\/\(2\ M\)]\ \(\[Sum]\+\(q = 0\)\%\(N - 1\)x[p, \ 
                  q]\ \ \ Cos[\(\(\ \)\(\((2  q + 1)\)\ n\ \[Pi]\)\)\/\(2\ \
N\)]\)\)\)], "NumberedEquation",
  TextAlignment->Center],

Cell["lo que nos lleva a definir las matrices siguientes.", "Text"],

Cell[BoxData[{
    \(R[m_, n_, M_, K_] := Table[\(u[m, M]\)[\([j]\)]*u[n, K], {j, 1, M}] (*\ 
      M\ filas, \ K\ columnas\ *) \), "\[IndentingNewLine]", 
    \(R[m_, n_, K_] := R[m, n, K, K]\)}], "Input"],

Cell[TextData[{
  "El siguiente resultado prueba que las matrices R[i,j,M,N] son linealmente \
independientes (lo hacemos en el caso M=N=4) y, por tanto, constituyen una \
base de ",
  Cell[BoxData[
      \(TraditionalForm\`\[ScriptCapitalM]\_\(M\[Times]N\)\)]],
  "."
}], "Text"],

Cell[BoxData[
    \(Solve[
        Sum[Flatten[t[i, j]*N[R[i, j, 4]]], {i, 0, 3}, {j, 0, 3}] \[Equal] 
          Flatten[ZeroMatrix[4]]] /. \((z_ \[Rule] u_)\) \[Rule] u\)], "Input"],

Cell[TextData[{
  "Estas matrices tambi\[EAcute]n tienen algunas propiedades de \
\"ortogonalidad matricial\". El siguiente comando nos dice que el producto \
matricial R[j,k,M].R[p,q,M] es la matriz nula si p\[NotEqual]k y, cuando p=k \
es una matriz cuya norma eucl\[IAcute]dea (vista como vector de ",
  Cell[BoxData[
      \(TraditionalForm\`\[DoubleStruckCapitalC]\^\(\[ScriptCapitalM]\^2\)\)]]\
,
  ") es igual a 1. Lo hacemos para M=4 y la raz\[OAcute]n de calcular las \
normas es para que no salgan en pantalla 256 matrices 4\[Times]4."
}], "Text"],

Cell[BoxData[
    \(\(\(Table[
      Norm[Flatten[N[R[j, k, 4] . R[p, q, 4]]]] // Chop, {j, 0, 3}, {k, 0, 
        3}, {p, 0, 3}, {q, 0, 3}]\)\( (*\ 
      es\ igual\ a\ 1\ cuando\ p = 
        k\ y\ es\ 0\ en\ otro\ caso\ *) \)\)\)], "Input"],

Cell[BoxData[
    \(Table[
      Norm[Flatten[N[R[j, k, 4] . R[k, q, 4]]]] // Chop, {j, 0, 3}, {k, 0, 
        3}, {q, 0, 3}]\)], "Input"],

Cell["\<\
La matriz R[m,n,K] proporciona el elemento X[m,n] de la matriz \
transformada. \
\>", "Text"],

Cell[BoxData[
    \(new2DCT[X_?MatrixQ] := 
      Module[{M, K}, M = \(Dimensions[X]\)[\([1]\)]; 
        K = \(Dimensions[X]\)[\([2]\)]; 
        Table[Sum[\(R[m, n, M, K]\)[\([j]\)] . X[\([j]\)], {j, 1, M}], {m, 0, 
            M - 1}, {n, 0, K - 1}]]\)], "Input"],

Cell[TextData[{
  "Podemos comprobar f\[AAcute]cilmente que este comando coincide con el \
antes definido \"twoDCT[ ]\". Definamos antes las matrices de la base can\
\[OAcute]nica de ",
  Cell[BoxData[
      \(TraditionalForm\`\[ScriptCapitalM]\_\(M\[Times]N\)\)]],
  ". "
}], "Text"],

Cell[BoxData[{
    \(H[i_, j_, M_, K_] := 
      Module[{a}, a = Table[0, {m, 1, M}, {n, 1, K}]; 
        For[m = 1, m \[LessEqual] M, \(m++\), 
          For[n = 1, n \[LessEqual] K, \(n++\), 
            If[{m, n} \[Equal] {i, j}, a[\([m, n]\)] = 1, 
              a[\([m, n]\)] = 0]]]; a]\), "\[IndentingNewLine]", 
    \(H[i_, j_, M_] := H[i, j, M, M]\)}], "Input"],

Cell[TextData[{
  "Comprobemos que new2DCT[ ] y twoDCT[ ] coinciden en todas las matrices de \
la base can\[OAcute]nica de ",
  Cell[BoxData[
      \(TraditionalForm\`M\_\(4\[Times]4\)\)]],
  " con lo cual, por linealidad, coincidir\[AAcute]n en toda matriz de ",
  Cell[BoxData[
      \(TraditionalForm\`M\_\(4\[Times]4\)\)]],
  "."
}], "Text"],

Cell[BoxData[
    \(Table[
      Simplify[Chop[twoDCT[H[i, j, 4]] - new2DCT[H[i, j, 4]]] \[Equal] 
          0], {i, 1, 4}, {j, 1, 4}]\)], "Input"],

Cell["La transformada inversa 2D-DCT viene dada por: ", "Text",
  TextAlignment->Left],

Cell[BoxData[
    FormBox[
      FrameBox[\(x[p, \ 
            q]\  = \ \ \[Sum]\+\(m = 0\)\%\(M - 1\)\ \(\[Sum]\+\(n = 0\)\%\(N \
- 1\)\[Alpha][m, M]\ \ \[Alpha][n, N] 
                X[m, \ n]\ \ Cos[\(\(\ \)\(\((2  p + 1)\)\ m\ \[Pi]\)\)\/\(2\ \
M\)]\ \ Cos[\(\((2  q + 1)\)\ n\ \[Pi]\)\/\(2\ N\)]\)\)], 
      TraditionalForm]], "NumberedEquation",
  TextAlignment->Left],

Cell[TextData[{
  "para  ",
  Cell[BoxData[
      \(TraditionalForm\`p = \ 0, 1,  ... , M - 1, \ q = 0, 1,  ... , 
      N - 1\)]],
  ". "
}], "Text"],

Cell[TextData[{
  "Realmente, lo que estamos haciendo es usar la base \
{R[m,n,K],{m,0,M-1},{n,0,K-1}} en el espacio de las matrices ",
  Cell[BoxData[
      \(TraditionalForm\`\[ScriptCapitalM]\_\(M\[Times]N\) = \
\[DoubleStruckCapitalC]\^\(M\[Times]N\)\)]],
  " para expresar en ella la matriz dada ",
  Cell[BoxData[
      \(TraditionalForm\`x = \(x[p, q]\)\_\(0 \[LessEqual] p \[LessEqual] M - \
1, 0 \[LessEqual] q \[LessEqual] N - 1\)\)]],
  ". "
}], "Text"],

Cell["\<\
El siguiente comando, que se deduce de la igualdad (8),  define la \
transformaci\[OAcute]n inversa 2D - DCT por medio de una combinaci\[OAcute]n \
lineal de las matrices b\[AAcute]sicas R[m,n,K].\
\>", "Text"],

Cell[BoxData[
    \(inverse2DCT[X_?MatrixQ] := 
      Module[{M, K, temp}, M = \(Dimensions[X]\)[\([1]\)]; 
        K = \(Dimensions[X]\)[\([2]\)]; 
        Sum[X[\([m + 1, n + 1]\)] R[m, n, M, K], {m, 0, M - 1}, {n, 0, 
              K - 1}] // Chop]\)], "Input"],

Cell["\<\
Cuando aplicamos este comando a la transformada 2D - DCT de una \
matriz obtenemos una identidad (justamente la expresada en la igualdad (8)). \
Es decir, expresamos la matriz original en la base \
{R[m,n,K],{m,0,M-1},{n,0,K-1}}. Observa que las coordenadas de la matriz dada \
en dicha base son precisamente los elementos de la matriz transformada.\
\>", \
"Text"],

Cell[BoxData[{
    \(twoDCT[H[1, 1, 4]]\), "\[IndentingNewLine]", 
    \(inverse2DCT[%]\)}], "Input"],

Cell["\<\
Podemos definir c\[OAcute]modamente un comando r\[AAcute]pido para \
calcular la inversa de la 2D - DCT.\
\>", "Text"],

Cell[BoxData[
    \(fastinverse2DCT[X_?MatrixQ] := 
      Map[inversefastDCT, 
        Transpose[Map[inversefastDCT, Transpose[X]]]]\)], "Input"],

Cell[BoxData[
    \(fastinverse2DCT[twoDCT[H[3, 3, 5]]]\)], "Input"],

Cell[TextData[{
  "Aqu\[IAcute] puedes ver una representaci\[OAcute]n gr\[AAcute]fica de las \
64 matrices 2D-DCT b\[AAcute]sicas  correspondientes a ",
  Cell[BoxData[
      \(TraditionalForm\`M = \(N = 8\)\)]],
  ". "
}], "Text"],

Cell[BoxData[
    \(\(Show[
        GraphicsArray[
          Map[\((ListDensityPlot[Reverse[#], Mesh \[Rule] False, 
                  FrameTicks \[Rule] None, 
                  DisplayFunction \[Rule] Identity]\  &)\), 
            Table[R[p, q, 8], {p, 0, 7}, {q, 0, 7}], {2}]]];\)\)], "Input"],

Cell["\<\
Veamos qu\[EAcute] pasa si consideramos las matrices \
b\[AAcute]sicas R[p,q,M,K] como vectores.\
\>", "Text"],

Cell[BoxData[{
    \(V[p_, q_, M_, K_] := Flatten[R[p, q, M, K]] (*\ 
      la\ matriz\ como\ vector\ de\ \[DoubleStruckCapitalC]\^\(M\[Times]K\)\ *) \
\), "\n", 
    \(V[p_, q_, M_] := Flatten[R[p, q, M]]\)}], "Input"],

Cell[TextData[{
  "A continuaci\[OAcute]n comprobamos que dichos vectores forman una base \
ortonormal de ",
  Cell[BoxData[
      \(TraditionalForm\`\[DoubleStruckCapitalC]\^\(M\[Times]N\)\)]],
  "."
}], "Text"],

Cell[BoxData[
    \(\(\(Table[
      Norm[N[V[p, q, 4] . V[j, k, 4]]] // Chop, {p, 0, 3}, {q, 0, 3}, {j, 0, 
        3}, {k, 0, 3}]\)\( (*\ todos\ los\ productos\ escalares, \ 
      en\ total\ 16\[Times]16 = 256, \ 
      son\ nulos\ excepto\ cuanado\ p = \(j\ y\ q = k\)\ *) \)\)\)], "Input"]
}, Open  ]]
}, Open  ]],

Cell[CellGroupData[{

Cell["Ampliaci\[OAcute]n de im\[AAcute]genes", "Section"],

Cell["\<\
Como sabemos, una imagen no es m\[AAcute]s que una matriz num\
\[EAcute]rica. A continuaci\[OAcute]n vamos a importar una imagen muy famosa.\
\
\>", "Text"],

Cell[BoxData[
    \(\(SetDirectory["\<D:/Trabajo/02_Docencia \
Informatica/Ingenieria_Informatica/Curso2006_2007/Calculo \
Avanzado/practicas/DCT&imagen/\>"];\)\)], "Input"],

Cell[BoxData[
    \(\(malena = Import["\<malena.dat\>"];\)\)], "Input"],

Cell[BoxData[{
    \(InputForm[malena] // Shallow\  (*\ 
      es\ una\ matriz\ *) \), "\[IndentingNewLine]", 
    \(Dimensions[malena]\), "\[IndentingNewLine]", 
    \(Take[malena[\([289]\)], 17] (*\ 
      sus\ elementos\ son\ n\[UAcute]meros\ enteros\ comprendidos\ entre\ 0\ \
y\ 255\ *) \), "\[IndentingNewLine]", 
    \(\(\(malena = malena/255 // N;\)\(\[IndentingNewLine]\)
    \) (*\ 
      dividimos\ por\ 255\ y\ convertimos\ los\ datos\ en\ valores\ reales\ \
comprendidos\ entre\ 0\ y\ 1\ *) \), "\[IndentingNewLine]", 
    \(Take[malena[\([289]\)], 17]\)}], "Input"],

Cell[BoxData[{
    \(\(foto = 
        Show[Graphics[Raster[malena]], AspectRatio \[Rule] Automatic, 
          ImageSize \[Rule] 0.5*{512, 512}];\)\), "\[IndentingNewLine]", 
    \( (*\ mitad\ de\ su\ tama\[NTilde]o\ natural\ *) \)}], "Input"],

Cell["\<\
Vamos a seleccionar un rect\[AAcute]ngulo en esta imagen que \
comprenda la zona de los ojos. Para ello, se selecciona la imagen y se pulsa \
\"Control\" y, entonces, al mover el puntero sobre la imagen nos va mostrando \
las coordenadas de cada punto de la misma. Las coordenadas que se ven son de \
la forma {n\.bacolumna, n\.bafila} y, una vez elegidas las coordenadas \
{n\.bacolumna1, n\.bafila1}, {n\.bacolumna2, n\.bafila2} de dos \
v\[EAcute]rtices opuestos del rect\[AAcute]ngulo que se quiere selelccionar, \
hay que trasponerlas para elegir las filas y columnas correspondientes de la \
matriz {n\.bafila1,n\.bafila2},{n\.bacolumna1,n\.bacolumna2}.\
\>", "Text"],

Cell[BoxData[
    \(\(\(\(ojos = Take[malena, \ {229, 281}, {237, \ 373}];\) (*\ 
      que\ corresponde\ a\ las\ coordendas\ del\ v\[EAcute]rtice\ inferior \
{237, 229}\ y\ a\ las\ del\ v\[EAcute]rtice\ superior\ {373, 
          281}\ del\ rect\[AAcute]ngulo\ seleccionado\ *) \[IndentingNewLine]
    Dimensions[ojos]\)\( (*\ 53\  = \ 281 - 228, \ 
      137\  = \ 373 - 236\ *) \)\)\)], "Input"],

Cell["\<\
Representamos los ojos ampliando cuatro veces su tama\[NTilde]o \
natural.\
\>", "Text"],

Cell[BoxData[
    \(\(\(Show[Graphics[Raster[ojos]], AspectRatio \[Rule] Automatic, 
        ImageSize \[Rule] 4*{137, 53}];\)\( (*\ 
      ampliamos\ 4\ veces\ *) \)\)\)], "Input"],

Cell["\<\
Al ampliar la imagen se pierde calidad. Podemos obtener un \
resultado mejor usando la 2D - DCT para interpolar la imagen que queremos \
ampliar por medio de un polinomio trigonom\[EAcute]trico en dos variables. \
Seg\[UAcute]n sabemos, por la igualdad (8), se verifica que\
\>", "Text"],

Cell[BoxData[
    FormBox[
      FrameBox[\(ojos\[LeftDoubleBracket]p, \ 
            q\[RightDoubleBracket] = \[Sum]\+\(m = 1\)\%53\ \(\[Sum]\+\(n = 1\
\)\%137 \[Alpha][m - 1, 53]\ \ \[Alpha][n - 1, 137] 
                OJOS\[LeftDoubleBracket]m, \ 
                  n\[RightDoubleBracket]\ Cos[\(\((2  p - 1)\) \((\ m - 1)\)\ \
\[Pi]\)\/\(2\ *53\)]\ Cos[\(\((2  q - 1)\) \((\ n - 1)\)\ \[Pi]\)\/\(2\ \
*137\)]\)\)], TraditionalForm]], "NumberedEquation",
  TextAlignment->Left],

Cell["\<\
Donde \"OJOS\" es la transformada 2D - DCT de \"ojos\". Esto nos \
lleva a definir el siguiente polinomio trigonom\[EAcute]trico en dos \
variables.\
\>", "Text"],

Cell[BoxData[
    \(\(OJOS = fast2DCT[ojos];\)\)], "Input"],

Cell[BoxData[
    \(\(poli[x_, 
          y_] = \[Sum]\+\(m = 1\)\%53\(\[Sum]\+\(n = 1\)\%137 \[Alpha][m - 1, 
                53. ]*\[Alpha][n - 1, 137. ] 
              OJOS\[LeftDoubleBracket]m, n\[RightDoubleBracket]*
              Cos[\(\(\((\ m - 1)\) N[\[Pi]]\)\/53\) x]\ *
              Cos[\(\(\((\ n - 1)\)\ N[\[Pi]]\)\/137\) y]\);\)\)], "Input"],

Cell[TextData[{
  "Este polinomio coincide en los puntos de la forma {x,y}={p-1/2,q-1/2} \
donde ",
  Cell[BoxData[
      \(TraditionalForm\`1 \[LessEqual] p \[LessEqual] 53\)]],
  ", ",
  Cell[BoxData[
      \(TraditionalForm\`1 \[LessEqual] q \[LessEqual] 137\)]],
  " con ojos[[p,q]] y, en consecuencia, es un polinomio que interpola a la \
matriz \"ojos\" en dichos puntos. Comprobemos algunos valores."
}], "Text"],

Cell[BoxData[
    \(Table[
      poli[p - 1/2, q - 1/2] \[Equal] ojos[\([p, q]\)], {p, 1, 53, 10}, {q, 
        1, 137, 25}]\)], "Input"],

Cell["Podemos representar gr\[AAcute]ficamente la funci\[OAcute]n \
poli[x,y].", "Text"],

Cell[BoxData[
    \(\(Plot3D[Evaluate[poli[x, y]], {x, 0, 53}, {y, 0, 137}, 
        PlotPoints \[Rule] 60];\)\)], "Input"],

Cell["Tambi\[EAcute]n podemos representra en 3D la matriz \"ojos\".", "Text"],

Cell[BoxData[
    \(\(ListPlot3D[Transpose[ojos]];\)\)], "Input"],

Cell["\<\
Aunuqe las dos gr\[AAcute]ficas parecen iguales, la primera es una \
funci\[OAcute]n continua y la segunda es una matriz. Podemos ahora ampliar la \
gr\[AAcute]fica de la funci\[OAcute]n poli[x,y] cuatro veces sin perder \
calidad. Para ello, \"dilatamos\" la funci\[OAcute]n poli[x,y].\
\>", "Text"],

Cell[BoxData[
    \(newpoli[x_, y_] := poli[0.25*x, 0.25*y]\)], "Input"],

Cell["\<\
La funci\[OAcute]n newpoli[x,y] toma en el intervalo \
[0,4\[Times]53]\[Times][0,4\[Times]137] los mismos valores que toma poli[x,y] \
en el intervalo [0,53]\[Times][0,137]. Aqu\[IAcute] podemos ver la misma gr\
\[AAcute]fica anterior, sin p\[EAcute]rdida de calidad, cuatro veces m\
\[AAcute]s grande (observa los valores en los ejes).\
\>", "Text"],

Cell[BoxData[
    \(\(Plot3D[Evaluate[newpoli[x, y]], {x, 0, 4*53}, {y, 0, 4*137}, 
        PlotPoints \[Rule] 60];\)\)], "Input"],

Cell[TextData[{
  "Podemos utilizar ahora \[NTilde]a funci\[OAcute]n newpoli[x,y] para \
ampliar nuestra matriz \"ojos\" a una matriz 4 veces mayor. Para ello \
discretizamos la funci\[OAcute]n newpoli[x,y] calculando sus valores en los \
puntos de la forma {p-1/2,q-1/2} para ",
  Cell[BoxData[
      \(TraditionalForm\`1 \[LessEqual] p \[LessEqual] 4*53, \ 
      1 \[LessEqual] q \[LessEqual] 4*137\)]],
  ". Es decir, tendr\[IAcute]amos que calcular la siguiente matriz."
}], "Text"],

Cell["Table[newpoli[m-1/2,n-1/2],{m,1,4*53},{n,1,4*137}]", "Text",
  TextAlignment->Center],

Cell[TextData[{
  StyleBox["Mathematica",
    FontSlant->"Italic"],
  " (o tu CPU) no puede hacer este c\[AAcute]lculo. Observa que queremos \
hacer 4*53*4*137=116176 evaluaciones de la funci\[OAcute]n newpoli[x,y] y \
cada una de ellas supone hacer 53*137*2=14522 multiplicaciones y 53*137=7261 \
sumas. En total 116176*14522=1687107872 multiplicaciones y \
116176*7261=843553936 sumas. Demasiado. Pero, observando lo que queremos \
hacer, caemos en la cuenta de que podemos usar la transformaci\[OAcute]n \
inversa de la 2D - DCT para lograr nuestro objetivo y esto tiene la ventaja \
de que podemos usar el comando fastinverse2DCT[ ] antes definido. Tenemos \
que"
}], "Text"],

Cell[BoxData[
    \(TraditionalForm\`\(newpoli[p - 1/2, 
          q - 1/2] = \[Sum]\+\(m = 1\)\%53\ \(\[Sum]\+\(n = 1\)\%137 \[Alpha][
                m - 1, 53]\ \ \[Alpha][n - 1, 137] 
              OJOS\[LeftDoubleBracket]m, \ 
                n\[RightDoubleBracket]\ Cos[\(\((2  p - 1)\) \((\ m - 1)\) N[\
\[Pi]]\)\/\(2*53*4\)]\ Cos[\(\((2  q - 1)\) \((\ n - 1)\)\ \
N[\[Pi]]\)\/\(2*137*4\)]\);\)\)], "Text"],

Cell["\<\
Esto puede scribirse como una transformada inversa 2D - DCT. Para \
ello ampliamos la matriz OJOS a\[NTilde]adi\[EAcute]ndole ceros hasta \
convertirla en una matriz cuatro veces mayor.\
\>", "Text"],

Cell[BoxData[
    \(\(newOJOS = 
        Module[{temp}, temp = Table[0. , {m, 1, 4*53}, {n, 1, 4*137}]; 
          For[m = 1, m \[LessEqual] 4*53, \(m++\), 
            For[n = 1, n \[LessEqual] 4*137, \(n++\), 
              If[m \[LessEqual] 53 && n \[LessEqual] 137, 
                temp[\([m, n]\)] = OJOS[\([m, n]\)], 
                temp[\([m, n]\)] = 0. ]]]; temp];\)\)], "Input"],

Cell["\<\
Ahora ya es f\[AAcute]cil darse cuenta de que se verifica la \
siguiente igualdad.\
\>", "Text"],

Cell[BoxData[
    \(TraditionalForm\`\(newpoli[p - 1/2, q - 1/2] = 
        4 \(\[Sum]\+\(m = 1\)\%\(4*53\)\ \(\[Sum]\+\(n = 1\)\%\(4*137\)\
\[Alpha][m - 1, 4*53]\ \ \[Alpha][n - 1, 4*137] 
                newOJOS\[LeftDoubleBracket]m, \ 
                  n\[RightDoubleBracket]\ Cos[\(\((2  p - 1)\) \((\ m - 1)\) \
N[\[Pi]]\)\/\(2*53*4\)]\ Cos[\(\((2  q - 1)\) \((\ n - 1)\)\ \
N[\[Pi]]\)\/\(2*137*4\)]\)\);\)\)], "Text"],

Cell["\<\
Y la expresi\[OAcute]n a la derecha del igual se calcula por medio \
de la transformada inversa 2D - DCT de la matriz \"newOJOS\".\
\>", "Text"],

Cell[BoxData[
    \(\(zoomojos = 4*fastinverse2DCT[newOJOS];\)\)], "Input"],

Cell["\<\
Por si tienes alguna duda, podemos comprobar algunos valores.\
\>", \
"Text"],

Cell[BoxData[
    \(Table[
      zoomojos[\([p, q]\)] == newpoli[p - 1/2, q - 1/2], {p, 100, 200, 
        25}, {q, 300, 500, 50}]\)], "Input"],

Cell["Finalmente, s\[OAcute]lo queda mostrar la imagen ampliada.", "Text"],

Cell[BoxData[
    \(\(Show[Graphics[Raster[zoomojos]], AspectRatio \[Rule] Automatic, 
        ImageSize \[Rule] 4*{137, 53}];\)\)], "Input"],

Cell["\<\
Al compara esta imagen con la anterior se aprecia enseguida que la \
hemos mejorado mucho (ya no se observa \"pixelaci\[OAcute]n\").\
\>", "Text"],

Cell[BoxData[
    \(\(Show[Graphics[Raster[ojos]], AspectRatio \[Rule] Automatic, 
        ImageSize \[Rule] 4*{137, 53}];\)\)], "Input"],

Cell["\<\
Podemos ahora convertir el proceso anterior en un comando que \
acepte como entrada una imagen y un factor de ampliaci\[OAcute]n y \
proporcione como salida la imagen ampliada.\
\>", "Text"],

Cell[BoxData[
    \(zoomimagen[X_?MatrixQ, factor_?IntegerQ] := 
      Module[{M, K, P, Q, temp, tempampliada}, 
        M = \(Dimensions[X]\)[\([1]\)]; K = \(Dimensions[X]\)[\([2]\)]; 
        P = factor*M; Q = factor*K; \ 
        temp = fast2DCT[X]; \[IndentingNewLine]tempampliada = 
          Module[{tempceros}, tempceros = Table[0. , {m, 1, P}, {n, 1, Q}]; 
            For[m = 1, m \[LessEqual] P, \(m++\), 
              For[n = 1, n \[LessEqual] Q, \(n++\), 
                If[m \[LessEqual] M && n \[LessEqual] K, 
                  tempceros[\([m, n]\)] = temp[\([m, n]\)], 
                  tempceros[\([m, n]\)] = 0. ]]]; 
            tempceros]; \[IndentingNewLine]Show[
          Graphics[Raster[factor*fastinverse2DCT[tempampliada]]], 
          AspectRatio \[Rule] Automatic, 
          ImageSize \[Rule] {Q, P}];]\)], "Input"],

Cell[BoxData[
    \(zoomimagen[ojos, 4]\)], "Input"],

Cell["Hagamos otra prueba. Seleecionamos los labios.", "Text"],

Cell[BoxData[{
    \(\(labios = 
        Take[malena, {144, 175}, {258, 327}];\)\), "\[IndentingNewLine]", 
    \(\(Show[Graphics[Raster[labios]], AspectRatio \[Rule] Automatic, 
        ImageSize \[Rule] 4*Reverse[Dimensions[labios]]];\)\)}], "Input"],

Cell[BoxData[
    \(zoomimagen[labios, 4]\)], "Input"],

Cell["\<\
Lo llamativo es que, a simple vista, parece que la calidad se \
conserva al aumentar el tama\[NTilde]o de la ampliaci\[OAcute]n.\
\>", "Text"],

Cell[BoxData[
    \(\(Do[zoomimagen[labios, k], {k, 1, 6}];\)\)], "Input"],

Cell["\<\
Lo hecho no ser\[IAcute]a posible sin el algoritmo FFT... \
\[DownExclamation]Las matem\[AAcute]ticas son eficaces!\
\>", "Text"]
}, Open  ]]
},
FrontEndVersion->"5.1 for Microsoft Windows",
ScreenRectangle->{{0, 1024}, {0, 685}},
AutoGeneratedPackage->None,
WindowSize->{1016, 651},
WindowMargins->{{0, Automatic}, {Automatic, 0}},
PageHeaders->{{Cell[ 
        TextData[ 
        "F. Javier P\[EAcute]rez Gonz\[AAcute]lez\nDpto. An\[AAcute]lisis \
Matem\[AAcute]tico - Universidad de Granada"], "Header"], Cell[ 
        TextData[ {
          CounterBox[ "Page"]}], "PageNumber"], Cell[ 
        TextData[ 
        "Transformada Coseno Discreta. Ampliaci\[OAcute]n de \
im\[AAcute]genes"], "Header"]}, {Cell[ 
        TextData[ 
        "F. Javier P\[EAcute]rez Gonz\[AAcute]lez\nDpto. An\[AAcute]lisis \
Matem\[AAcute]tico - Universidad de Granada"], "Header"], Cell[ 
        TextData[ {
          CounterBox[ "Page"]}], "PageNumber"], Cell[ 
        TextData[ 
        "Transformada Coseno Discreta. Ampliaci\[OAcute]n de \
im\[AAcute]genes"], "Header"]}},
ShowSelection->True,
Magnification->1,
StyleDefinitions -> "Classroom.nb"
]

(*******************************************************************
Cached data follows.  If you edit this Notebook file directly, not
using Mathematica, you must remove the line containing CacheID at
the top of  the file.  The cache data will then be recreated when
you save this file from within Mathematica.
*******************************************************************)

(*CellTagsOutline
CellTagsIndex->{}
*)

(*CellTagsIndex
CellTagsIndex->{}
*)

(*NotebookFileOutline
Notebook[{
Cell[1754, 51, 87, 1, 67, "Title"],
Cell[1844, 54, 148, 4, 87, "Author"],

Cell[CellGroupData[{
Cell[2017, 62, 41, 0, 62, "Section"],
Cell[2061, 64, 227, 5, 82, "Input",
  InitializationCell->True]
}, Closed]],

Cell[CellGroupData[{
Cell[2325, 74, 160, 2, 42, "Section",
  CounterAssignments->{{"Section", 3}, {"Subsection", 0}, {
      "Subsubsection", 0}, {"NumberedEquation", 0}}],
Cell[2488, 78, 547, 9, 67, "Text"],

Cell[CellGroupData[{
Cell[3060, 91, 33, 0, 45, "Subsubsection"],
Cell[3096, 93, 461, 13, 29, "Text"],
Cell[3560, 108, 467, 12, 70, "NumberedEquation"],
Cell[4030, 122, 332, 7, 48, "Text"],
Cell[4365, 131, 443, 10, 70, "NumberedEquation"],
Cell[4811, 143, 64, 0, 29, "Text"],
Cell[4878, 145, 424, 11, 70, "NumberedEquation"],
Cell[5305, 158, 1181, 32, 27, "Text"],
Cell[6489, 192, 427, 13, 58, "NumberedEquation"],
Cell[6919, 207, 457, 13, 29, "Text"],
Cell[7379, 222, 283, 6, 42, "NumberedEquation"],
Cell[7665, 230, 41, 0, 29, "Text"],
Cell[7709, 232, 354, 8, 108, "Input"],
Cell[8066, 242, 180, 5, 29, "Text"],
Cell[8249, 249, 82, 1, 50, "Input"],
Cell[8334, 252, 130, 3, 70, "Input"],
Cell[8467, 257, 1092, 25, 98, "Text"],
Cell[9562, 284, 126, 2, 70, "Input"],
Cell[9691, 288, 421, 9, 90, "Input"],
Cell[10115, 299, 117, 2, 70, "Input"],
Cell[10235, 303, 59, 1, 50, "Input"],
Cell[10297, 306, 395, 7, 48, "Text"],
Cell[10695, 315, 128, 2, 70, "Input"],
Cell[10826, 319, 785, 16, 181, "Input"],
Cell[11614, 337, 109, 2, 70, "Input"],
Cell[11726, 341, 101, 2, 70, "Input"],
Cell[11830, 345, 73, 0, 29, "Text"],
Cell[11906, 347, 797, 16, 166, "Input"],
Cell[12706, 365, 139, 3, 29, "Text"],
Cell[12848, 370, 489, 8, 140, "Input"],
Cell[13340, 380, 92, 3, 29, "Text"],
Cell[13435, 385, 214, 4, 90, "Input"],
Cell[13652, 391, 123, 3, 29, "Text"],
Cell[13778, 396, 410, 8, 70, "Input"],
Cell[14191, 406, 73, 0, 29, "Text"],
Cell[14267, 408, 500, 10, 121, "Input"],
Cell[14770, 420, 68, 1, 50, "Input"],
Cell[14841, 423, 155, 3, 50, "Input"],
Cell[14999, 428, 206, 7, 29, "Text"],
Cell[15208, 437, 625, 12, 150, "Input"],
Cell[15836, 451, 454, 12, 48, "Text"],
Cell[16293, 465, 1118, 25, 150, "Input"]
}, Open  ]],

Cell[CellGroupData[{
Cell[17448, 495, 33, 0, 45, "Subsubsection"],
Cell[17484, 497, 617, 15, 30, "Text"],
Cell[18104, 514, 411, 9, 70, "NumberedEquation"],
Cell[18518, 525, 279, 9, 29, "Text"],
Cell[18800, 536, 375, 14, 29, "Text"],
Cell[19178, 552, 473, 8, 64, "Text"],
Cell[19654, 562, 412, 15, 29, "Text"],
Cell[20069, 579, 505, 11, 64, "Text"],
Cell[20577, 592, 1396, 38, 86, "Text"],
Cell[21976, 632, 792, 15, 159, "Input"],
Cell[22771, 649, 72, 0, 29, "Text"],
Cell[22846, 651, 176, 3, 70, "Input"],
Cell[23025, 656, 111, 3, 29, "Text"],
Cell[23139, 661, 439, 8, 150, "Input"],
Cell[23581, 671, 60, 0, 29, "Text"],
Cell[23644, 673, 370, 7, 62, "NumberedEquation"],
Cell[24017, 682, 67, 0, 29, "Text"],
Cell[24087, 684, 206, 3, 70, "Input"],
Cell[24296, 689, 280, 7, 29, "Text"],
Cell[24579, 698, 182, 3, 50, "Input"],
Cell[24764, 703, 557, 10, 50, "Text"],
Cell[25324, 715, 243, 5, 70, "Input"],
Cell[25570, 722, 138, 3, 50, "Input"],
Cell[25711, 727, 103, 3, 29, "Text"],
Cell[25817, 732, 266, 5, 70, "Input"],
Cell[26086, 739, 284, 7, 29, "Text"],
Cell[26373, 748, 369, 7, 90, "Input"],
Cell[26745, 757, 345, 9, 29, "Text"],
Cell[27093, 768, 147, 3, 50, "Input"],
Cell[27243, 773, 86, 1, 29, "Text"],
Cell[27332, 776, 377, 8, 68, "NumberedEquation"],
Cell[27712, 786, 150, 6, 29, "Text"],
Cell[27865, 794, 464, 11, 50, "Text"],
Cell[28332, 807, 220, 4, 29, "Text"],
Cell[28555, 813, 264, 5, 70, "Input"],
Cell[28822, 820, 375, 7, 48, "Text"],
Cell[29200, 829, 101, 2, 70, "Input"],
Cell[29304, 833, 128, 3, 29, "Text"],
Cell[29435, 838, 145, 3, 50, "Input"],
Cell[29583, 843, 68, 1, 50, "Input"],
Cell[29654, 846, 231, 6, 29, "Text"],
Cell[29888, 854, 297, 6, 70, "Input"],
Cell[30188, 862, 120, 3, 29, "Text"],
Cell[30311, 867, 219, 4, 70, "Input"],
Cell[30533, 873, 212, 6, 29, "Text"],
Cell[30748, 881, 294, 5, 70, "Input"]
}, Open  ]]
}, Open  ]],

Cell[CellGroupData[{
Cell[31091, 892, 57, 0, 62, "Section"],
Cell[31151, 894, 166, 4, 29, "Text"],
Cell[31320, 900, 173, 3, 70, "Input"],
Cell[31496, 905, 71, 1, 50, "Input"],
Cell[31570, 908, 579, 11, 150, "Input"],
Cell[32152, 921, 244, 4, 70, "Input"],
Cell[32399, 927, 683, 10, 67, "Text"],
Cell[33085, 939, 398, 6, 110, "Input"],
Cell[33486, 947, 98, 3, 29, "Text"],
Cell[33587, 952, 181, 3, 50, "Input"],
Cell[33771, 957, 297, 5, 48, "Text"],
Cell[34071, 964, 481, 9, 68, "NumberedEquation"],
Cell[34555, 975, 172, 4, 29, "Text"],
Cell[34730, 981, 59, 1, 50, "Input"],
Cell[34792, 984, 356, 6, 71, "Input"],
Cell[35151, 992, 419, 10, 48, "Text"],
Cell[35573, 1004, 137, 3, 50, "Input"],
Cell[35713, 1009, 88, 1, 29, "Text"],
Cell[35804, 1012, 123, 2, 50, "Input"],
Cell[35930, 1016, 77, 0, 29, "Text"],
Cell[36010, 1018, 65, 1, 50, "Input"],
Cell[36078, 1021, 310, 5, 48, "Text"],
Cell[36391, 1028, 72, 1, 62, "Input"],
Cell[36466, 1031, 360, 6, 62, "Text"],
Cell[36829, 1039, 130, 2, 48, "Input"],
Cell[36962, 1043, 487, 9, 62, "Text"],
Cell[37452, 1054, 91, 1, 27, "Text"],
Cell[37546, 1057, 679, 12, 96, "Text"],
Cell[38228, 1071, 413, 7, 71, "Text"],
Cell[38644, 1080, 209, 4, 44, "Text"],
Cell[38856, 1086, 389, 7, 97, "Input"],
Cell[39248, 1095, 106, 3, 27, "Text"],
Cell[39357, 1100, 424, 7, 71, "Text"],
Cell[39784, 1109, 154, 3, 27, "Text"],
Cell[39941, 1114, 75, 1, 62, "Input"],
Cell[40019, 1117, 87, 3, 27, "Text"],
Cell[40109, 1122, 143, 3, 80, "Input"],
Cell[40255, 1127, 74, 0, 27, "Text"],
Cell[40332, 1129, 141, 2, 80, "Input"],
Cell[40476, 1133, 156, 3, 27, "Text"],
Cell[40635, 1138, 137, 2, 62, "Input"],
Cell[40775, 1142, 200, 4, 44, "Text"],
Cell[40978, 1148, 847, 15, 201, "Input"],
Cell[41828, 1165, 52, 1, 62, "Input"],
Cell[41883, 1168, 62, 0, 27, "Text"],
Cell[41948, 1170, 252, 4, 97, "Input"],
Cell[42203, 1176, 54, 1, 62, "Input"],
Cell[42260, 1179, 152, 3, 27, "Text"],
Cell[42415, 1184, 74, 1, 62, "Input"],
Cell[42492, 1187, 139, 3, 27, "Text"]
}, Open  ]]
}
]
*)



(*******************************************************************
End of Mathematica Notebook file.
*******************************************************************)

