Architektura SIMD i udostępniane przez nią techniki programowania niskopoziomowego używane są w celu efektywniejszego wykonywania obliczeń poprzez równoległe przetwarzanie wielu bajtów danych. Technologia ta bazuje na spakowanych typach danych i specjalnych instrukcjach procesora.
SSE to skrót od Streaming SIMD Extensions. SSE to technologia będąca odpowiednikiem MMX, ale umozliwiająca operacje na liczbach zmiennoprzecinkowych. Rejestry SSE są 128-bitowe, mogą być użyte do operacji na danych różnego typu i w przeciwieństwie do rejestrów MMX są niezależnymi rejestrami procesora. Dobry opis (Google -> SSE instruction set) arytmetycznych instrukcji SSE można znaleźć tutaj.
Pełne zestawienie instrukcji SEE z wygodną wyszukiwarką.
Iloczyn dwóch tablic wyraz po wyrazie
Proszę napisać funkcję asemblerową diff, która korzysta z operacji SSE w celu wyliczenia numerycznej pochodnej z tablicy jednowymiarowej T zmiennych 8-bitowych (np.piksele obrazu).
W praktyce jest to różnica sąsiadujących pikseli.
DIFF[i] = W[i+1] - W[i]
Sposób użycia funkcji diff przedstawiono na Listingu 3.
Listing 3 Opakowanie dla funkcji diff
#include <stdio.h>
#define N 100
// na wyjściu out[i] = wiersz[i+1] - wiersz[i]
extern "C" void diff(char *out,char *wiersz, int n);
int main(void)
{
char Tablica[N+1], DIFF[N];
int i;
Tablica[0]=1;
for(i=1;i<=N;i++)
Tablica[i]=Tablica[i-1]+i;
diff(DIFF, Tablica, N);
for(i=0;i<N;i++)
printf("%d ",DIFF[i]);
printf("\n");
}
// OUT: 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 -128 17 18 19 20
Aby uniknąć przepełnień dobrze jest stosować odejmowanie z nasyceniem PSUBSB.
Korzystając z instrukcji SSE proszę napisać procedurę wyliczającą amplitudę gradientu obrazu - 8-bitowej bitmapy. Na Listingu 4 przedstawiony jest kod z instrukcjami pozwalającymi czytać i pisać do bitmapy 256-kolorowej. Dla polepszenia dokładności obliczeń gradientu, zarówno dane wejściowe jak i wyjściowe są typu float. Na Listingu 4 przedstawiono również prototyp funkcji liczącej gradient. Gradient GRAD(i,j) w pikselu (i,j) liczony jest według wzoru:
GRAD(i,j)= scale * sqrt((DATA(i+1,j)-DATA(i-1,j))2+(DATA(i,j+1)-DATA(i,j-1))2)
Listing 4 Opakowanie asemblerowej funkcji liczącej gradient obrazu
#include <stdio.h> #include <math.h> //Na wyjściu: grad[i] = scale * sqrt( (data[i+1] - data[i-1])^2 + (data[i+N] - data[i-N])^2) // dla i=0,1,...,N-2 extern "C" void gradientSSE(float *grad, float * data, int N, float scale); int main(void) { float data[400][400],header[1078]; float grad[400][400]; int N=400,HL=1078; int i,j; FILE *strm; strm=fopen("circle3.bmp","rb"); for(i=0;i<HL;i++) header[i]=fgetc(strm); for(i=0;i<N;i++) for(j=0;j<N;j++) data[i][j]=(float)fgetc(strm); fclose(strm); for(i=1;i<N-1;i++) gradientSSE(grad[i]+1,data[i]+1,N,7); for(i=0;i<N;i++){ grad[0][i]=255; grad[N-1][i]=255; grad[i][0]=255; grad[i][N-1]=255; } strm=fopen("dum.bmp","wb"); for(i=0;i<HL;i++) fputc(header[i],strm); for(i=0;i<N;i++) for(j=0;j<N;j++) fputc((unsigned char)grad[i][j],strm); fclose(strm); }
Przykładowe obrazy wejściowy (Rys. 1) i wyjściowy (Rys. 2) przedstawiono poniżej.
Uwaga: Obrazek należy otworzyć w nowym oknie a następnie zapisać! Bezpośrednie zapisanie może powodować błędny format pliku.
Rys. 1 Obraz wejściowy Rys. 2 Obraz wyjściowyUWAGA: Zaproponowana w zadaniu procedura liczenia gradientu z obrazu nie ma szczególnie dobrych własności i raczej nie jest wykorzystywana w praktyce. Filtry gradientowe można jednak zaimplementować, korzystając z analogicznych metod, jak rozważane w zadaniu.
Korzystając z instrukcji SSE proszę napisać procedurę zmniejszającą dwukrotnie rozmiar bitmapy przedstawionej na Rys. 1. Każdy stopień szarości piksela OUT(i,j) wynikowej bitmapy ma być średnią z czterech odpowiednich stopni szarości pikseli bitmapy wejściowej IN według wzoru:
OUT(i,j) = (IN(2*i,2*j) + IN(2*i +1,2*j) + IN(2*i,2*j + 1) + IN(2*i + 1,2*j + 1))/4
Na Listingu 5 przedstawiony jest kod opakowania dla procedury skalującej bitmapę wraz z prototypem i sposobem użycia tej funkcji.
#include <stdio.h>
void scaleSSE(float *,float *,int);
int main(void)
{
float data[400][400],dum[200][200];
unsigned char header[1078];
unsigned char ch;
int N=400,HL=1078;
int i,j;
FILE *strm;
strm=fopen("circle.bmp","rb");
for(i=0;i<HL;i++) header[i]=fgetc(strm);
for(i=0;i<N;i++)
for(j=0;j<N;j++)
data[i][j]=(float)fgetc(strm);
fclose(strm);
for(i=0;i<N/2;i++) scaleSSE(dum[i],data[2*i],N);
header[4]=0;
header[3]=160;
header[2]=118;
header[18]=200;
header[19]=0;
header[22]=200;
header[23]=0;
strm=fopen("dum.bmp","wb");
for(i=0;i<HL;i++) fputc(header[i],strm);
for(i=0;i<N/2;i++)
for(j=0;j<N/2;j++)
fputc((unsigned char)dum[i][j],strm);
fclose(strm);
}
Informacje o wielkości bitmapy są zapisane w modyfikowanych bajtach nagłówka. W bajtach od 2 poczynając zakodowany jest rozmiar bitmapy w bajtach (równy wysokość*szerokość zaokrąglona do najbliższej liczby podzielnej przez 4 + ilość bajtów nagłówka). W bajtach od 18 poczynając zakodowane są wysokość i szerokość bitmapy.