SIMD का उपयोग करके C ++ में इमेज प्रोसेसिंग का अनुकूलन। मेडियन फिल्टर

परिचय

इससे पहले परिचयात्मक लेख में, मैंने उन समस्याओं की एक सूची तैयार की, जिन्हें एक डेवलपर को सामना करना पड़ेगा यदि वह SIMI निर्देशों का उपयोग करके छवि प्रसंस्करण अनुकूलन का अनुकूलन करना चाहता है। अब समय आ गया है कि एक ठोस उदाहरण दिया जाए कि उपरोक्त समस्याओं को कैसे हल किया जा सकता है। मैंने एक लंबे समय के लिए सोचा था कि कौन सा एल्गोरिथम पहले उदाहरण के लिए चुनना है, और मंझला फ़िल्टरिंग पर ध्यान केंद्रित करने का निर्णय लिया। मेडियन फ़िल्टरिंग शोर को दबाने का एक प्रभावी तरीका है जो अनिवार्य रूप से कम रोशनी वाले दृश्यों में डिजिटल कैमरों पर दिखाई देता है। यह एल्गोरिथ्म काफी संसाधन-गहन है - उदाहरण के लिए, जब एक 3x3 माध्य फ़िल्टर के साथ एक ग्रे छवि को संसाधित किया जाता है, तो प्रति छवि बिंदु पर लगभग 50 संचालन की आवश्यकता होती है। लेकिन एक ही समय में, वह केवल 8-बिट संख्या के साथ काम करता है और उसे काम करने के लिए अपेक्षाकृत कम इनपुट की आवश्यकता होती है। ये परिस्थितियाँ SIMD अनुकूलन के लिए एल्गोरिथ्म को काफी सरल बनाती हैं और साथ ही साथ इससे बहुत महत्वपूर्ण त्वरण प्राप्त करना संभव बनाती हैं।




संदर्भ के लिए, मुझे एल्गोरिथ्म का सार याद है:

• मूल छवि के प्रत्येक बिंदु के लिए, एक निश्चित पड़ोस लिया जाता है (हमारे मामले में, 3x3)।
• इस पड़ोस में अंक चमक को बढ़ाकर हल किए जाते हैं।
• सॉर्ट किए गए पड़ोस का मिडपॉइंट (एक 3x3 फ़िल्टर के लिए 5 वां) अंतिम छवि में दर्ज किया गया है।

पहले उदाहरण के लिए, मैंने कार्य को आसान बनाने का फैसला किया और एक ग्रे छवि (8 बिट प्रति बिंदु) के मामले पर विचार करने के लिए सीमा के प्रभावों को ध्यान में रखे बिना - इसके लिए, थोड़ी बड़ी छवि के मध्य भाग को छानने के लिए लिया गया था। बेशक, वास्तविक समस्याओं में, छवि के सीमा बिंदुओं की गणना एक विशेष एल्गोरिथ्म का उपयोग करके की जानी चाहिए, और कभी-कभी इन विशेष एल्गोरिदम के कार्यान्वयन में अधिकांश कोड हो सकते हैं, लेकिन हम इन मुद्दों पर एक और समय का ध्यान रखेंगे। इसके अलावा, इसी तरह के कारणों के लिए, यह माना जाता है कि छवि की चौड़ाई 32 से अधिक है।

एल्गोरिथ्म के स्केलर संस्करण

1 स्केलर संस्करण

इस प्रारंभिक कार्यान्वयन में, समस्या को सिर पर हल किया गया था: मूल छवि के प्रत्येक बिंदु के लिए पड़ोस को एक सहायक सरणी में कॉपी किया गया था, जिसे तब std :: सॉर्ट () फ़ंक्शन का उपयोग करके सॉर्ट किया गया था।

inline void Load(const uint8_t * src, int a[3]) { a[0] = src[-1]; a[1] = src[0]; a[2] = src[1]; } inline void Load(const uint8_t * src, size_t stride, int a[9]) { Load(src - stride, a + 0); Load(src, a + 3); Load(src + stride, a + 6); } inline void Sort(int a[9]) { std::sort(a, a + 9); } void MedianFilter(const uint8_t * src, size_t srcStride, size_t width, size_t height, uint8_t * dst, size_t dstStride) { int a[9]; for(size_t y = 0; y < height; ++y) { for(size_t x = 0; x < width; ++x) { Load(src + x, srcStride, a); Sort(a); dst[x] = (uint8_t)a[4]; } src += srcStride; dst += dstStride; } } 

दूसरा स्केलर संस्करण

पिछली विधि, हालांकि, तुरंत एक अड़चन दिखाती है - यह मानक सॉर्टिंग फ़ंक्शन है। यह मनमाने आकार की एक सरणी को छांटने के लिए है, यह कई आंतरिक जांच करता है, इसलिए यह सबसे अधिक संभावना है कि इस समस्या का इष्टतम समाधान नहीं होगा, जहां सरणी छोटा है और एक निश्चित आकार है। इसके अलावा, सरणी को पूरी तरह से सॉर्ट करना हमारे लिए बिल्कुल भी आवश्यक नहीं है, यह पर्याप्त है कि हम जो मूल्य चाहते हैं वह 4 वें तत्व में है। इसलिए, हम "नेटवर्क" सॉर्टिंग विधि (जिसे बिटोनिक सॉर्टिंग विधि के रूप में भी जाना जाता है) को लागू कर सकते हैं:

 inline void Sort(int & a, int & b) //    { if(a > b) { int t = a; a = b; b = t; } } inline void Sort(int a[9]) //    { Sort(a[1], a[2]); Sort(a[4], a[5]); Sort(a[7], a[8]); Sort(a[0], a[1]); Sort(a[3], a[4]); Sort(a[6], a[7]); Sort(a[1], a[2]); Sort(a[4], a[5]); Sort(a[7], a[8]); Sort(a[0], a[3]); Sort(a[5], a[8]); Sort(a[4], a[7]); Sort(a[3], a[6]); Sort(a[1], a[4]); Sort(a[2], a[5]); Sort(a[4], a[7]); Sort(a[4], a[2]); Sort(a[6], a[4]); Sort(a[4], a[2]); } 

3 स्केलर संस्करण

यद्यपि पहला तरीका शून्य विधि की तुलना में बहुत तेज निकला (नीचे परीक्षण परिणाम देखें), इसमें अभी भी एक अड़चन है - इस पद्धति में बहुत अधिक सशर्त संक्रमण हैं। जैसा कि आप जानते हैं, आधुनिक प्रोसेसर उन्हें बहुत पसंद नहीं करते हैं, क्योंकि उनके पास एक पाइपलाइन वास्तुकला है और एक घड़ी चक्र में कई निर्देशों के असाधारण निष्पादन की संभावना है। पहले और दूसरे दोनों के लिए, सशर्त संक्रमण बेहद अवांछनीय हैं, क्योंकि प्रोसेसर को रोकना चाहिए और गणना के परिणामों की प्रतीक्षा करनी चाहिए ताकि यह पता लगाया जा सके कि कार्यक्रम की किस शाखा में आगे जाना चाहिए। सौभाग्य से, सशर्त शाखाओं में बंटी बिना दो 8-बिट मानों को क्रमित किया जा सकता है:

 inline void Sort(int & a, int & b) { int d = a - b; int m = ~(d >> 8); b += d&m; a -= d&m; } 

इन प्रारंभिक स्केलर कोड अनुकूलन का परिणाम मंझला फिल्टर के निम्न संस्करण होगा:

 inline void Load(const uint8_t * src, int a[3]) { a[0] = src[-1]; a[1] = src[0]; a[2] = src[1]; } inline void Load(const uint8_t * src, size_t stride, int a[9]) { Load(src - stride, a + 0); Load(src, a + 3); Load(src + stride, a + 6); } inline void Sort(int & a, int & b) { int d = a - b; int m = ~(d >> 8); b += d&m; a -= d&m; } inline void Sort(int a[9]) //    { Sort(a[1], a[2]); Sort(a[4], a[5]); Sort(a[7], a[8]); Sort(a[0], a[1]); Sort(a[3], a[4]); Sort(a[6], a[7]); Sort(a[1], a[2]); Sort(a[4], a[5]); Sort(a[7], a[8]); Sort(a[0], a[3]); Sort(a[5], a[8]); Sort(a[4], a[7]); Sort(a[3], a[6]); Sort(a[1], a[4]); Sort(a[2], a[5]); Sort(a[4], a[7]); Sort(a[4], a[2]); Sort(a[6], a[4]); Sort(a[4], a[2]); } void MedianFilter(const uint8_t * src, size_t srcStride, size_t width, size_t height, uint8_t * dst, size_t dstStride) { int a[9]; for(size_t y = 0; y < height; ++y) { for(size_t x = 0; x < width; ++x) { Load(src + x, srcStride, a); Sort(a); dst[x] = (uint8_t)a[4]; } src += srcStride; dst += dstStride; } } 

यह संस्करण पहले से ही काफी (कहीं-कहीं 5-6 गुना) गति में एल्गोरिथ्म के हमारे शुरुआती संस्करण की तुलना में तेज है। और यह इस पर आधारित है कि हम एल्गोरिदम के SIMD अनुकूलन को पूरा करेंगे, साथ ही साथ एल्गोरिथ्म के स्केलर और वेक्टर संस्करणों की गति की तुलना करेंगे।

एल्गोरिथ्म का SIMD संस्करण

SIMD का उपयोग करके माध्य फ़िल्टरिंग एल्गोरिथ्म के अनुकूलन का वर्णन करने के लिए, मैं निर्देशों के दो सेट SSE2 और AVX2 का उपयोग करूँगा, जो MMX एक्सटेंशन को याद कर रहा है, जो वर्तमान में पुराना है और इसमें अधिक ऐतिहासिक रुचि है। सौभाग्य से, SIMD निर्देशों का उपयोग करने के लिए, कोडांतरक का उपयोग करना आवश्यक नहीं है। । अधिकांश आधुनिक सी ++ कंपाइलरों में आंतरिक (अंतर्निहित कार्य जिसके साथ आप सभी प्रकार के प्रोसेसर एक्सटेंशन का उपयोग कर सकते हैं) के लिए समर्थन है। आंतरिक का उपयोग करके प्रोग्रामिंग व्यावहारिक रूप से शुद्ध सी पर प्रोग्रामिंग से अलग नहीं है। आंतरिक कार्यों को संकलक द्वारा सीधे प्रोसेसर निर्देशों में परिवर्तित किया जाता है, हालांकि प्रोसेसर रजिस्टरों के साथ सीधे काम करना प्रोग्रामर से छिपा रहता है। ज्यादातर मामलों में, इंट्रिंसिक्स का उपयोग करने वाला एक प्रोग्राम असेंबलर में लिखे गए प्रोग्राम की गति से हीन नहीं है।

SSE2 संस्करण

SSE2 पूर्णांक कमांड को हेडर फ़ाइल <emmintrin.h> में परिभाषित किया गया है। मूल प्रकार __m128i है - एक 128 बिट वेक्टर, जो कि संदर्भ के आधार पर, 2 64 बिट, 4 32 बिट, 8 x 16 बिट, 16 x 8 बिट पर हस्ताक्षर किए या अहस्ताक्षरित संख्याओं के एक सेट के रूप में व्याख्या की जा सकती है। जैसा कि आप देख सकते हैं, वे न केवल वेक्टर अंकगणितीय संचालन का समर्थन करते हैं, बल्कि वेक्टर तार्किक संचालन भी करते हैं, साथ ही वेक्टर डेटा लोडिंग और अनलोडिंग ऑपरेशन भी करते हैं। SSE2 निर्देशों का उपयोग करके माध्यिका फ़िल्टर के लिए निम्नलिखित अनुकूलन हैं। कोड, यह मुझे लगता है, सुंदर दृश्य है।

 #include < emmintrin.h > inline void Load(const uint8_t * src, __m128i a[3]) { a[0] = _mm_loadu_si128((__m128i*)(src - 1)); // 128      16    a[1] = _mm_loadu_si128((__m128i*)(src)); a[2] = _mm_loadu_si128((__m128i*)(src + 1)); } inline void Load(const uint8_t * src, size_t stride, __m128i a[9]) { Load(src - stride, a + 0); Load(src, a + 3); Load(src + stride, a + 6); } inline void Sort(__m128i& a, __m128i& b) { __m128i t = a; a = _mm_min_epu8(t, b); //  2- 8       16   b = _mm_max_epu8(t, b); //  2- 8       16   } inline void Sort(__m128i a[9]) //    { Sort(a[1], a[2]); Sort(a[4], a[5]); Sort(a[7], a[8]); Sort(a[0], a[1]); Sort(a[3], a[4]); Sort(a[6], a[7]); Sort(a[1], a[2]); Sort(a[4], a[5]); Sort(a[7], a[8]); Sort(a[0], a[3]); Sort(a[5], a[8]); Sort(a[4], a[7]); Sort(a[3], a[6]); Sort(a[1], a[4]); Sort(a[2], a[5]); Sort(a[4], a[7]); Sort(a[4], a[2]); Sort(a[6], a[4]); Sort(a[4], a[2]); } void MedianFilter(const uint8_t * src, size_t srcStride, size_t width, size_t height, uint8_t * dst, size_t dstStride) { __m128i a[9]; for(size_t y = 0; y < height; ++y) { for(size_t x = 0; x < width; x += sizeof(__m128i)) { Load(src + x, srcStride, a); Sort(a); _mm_storeu_si128((__m128i*)(dst + x), a[4]); // 128      16    } src += srcStride; dst += dstStride; } } 

AVX2 संस्करण

पूर्णांक AVX2 कमांड को हेडर फ़ाइल <immintrin.h> में परिभाषित किया गया है। मूल प्रकार __m256i है - एक 256-बिट वेक्टर, जो कि संदर्भ के आधार पर, 4 64-बिट, 8-बिट 32, 16-बिट 16, 32-बिट 8-बिट हस्ताक्षरित या अनिश्चित संख्याओं के एक सेट के रूप में व्याख्या की जा सकती है। हालांकि AVX2 इंस्ट्रक्शन सेट काफी हद तक SSE2 इंस्ट्रक्शन सेट (वेक्टर की दोगुनी चौड़ाई को देखते हुए) को दोहराता है, लेकिन इसमें कुछ निर्देश भी होते हैं जिनका SSE2 में कोई एनालॉग नहीं है। नीचे AVX2 निर्देशों का उपयोग करके मंझला फिल्टर का अनुकूलन है।

 #include < immintrin.h > inline void Load(const uint8_t * src, __m256i a[3]) { a[0] = _mm256_loadu_si256((__m128i*)(src - 1)); // 256      32    a[1] = _mm256_loadu_si256((__m128i*)(src)); a[2] = _mm256_loadu_si256((__m128i*)(src + 1)); } inline void Load(const uint8_t * src, size_t stride, __m256i a[9]) { Load(src - stride, a + 0); Load(src, a + 3); Load(src + stride, a + 6); } inline void Sort(__m256i& a, __m256i& b) { __m256i t = a; a = _mm256_min_epu8(t, b); //  2- 8       32   b = _mm256_max_epu8(t, b); //  2- 8       32   } inline void Sort(__m256i a[9]) //    { Sort(a[1], a[2]); Sort(a[4], a[5]); Sort(a[7], a[8]); Sort(a[0], a[1]); Sort(a[3], a[4]); Sort(a[6], a[7]); Sort(a[1], a[2]); Sort(a[4], a[5]); Sort(a[7], a[8]); Sort(a[0], a[3]); Sort(a[5], a[8]); Sort(a[4], a[7]); Sort(a[3], a[6]); Sort(a[1], a[4]); Sort(a[2], a[5]); Sort(a[4], a[7]); Sort(a[4], a[2]); Sort(a[6], a[4]); Sort(a[4], a[2]); } void MedianFilter(const uint8_t * src, size_t srcStride, size_t width, size_t height, uint8_t * dst, size_t dstStride) { __m256i a[9]; for(size_t y = 0; y < height; ++y) { for(size_t x = 0; x < width; x += sizeof(__m256i)) { Load(src + x, srcStride, a); Sort(a); _mm256_storeu_si256((__m256i*)(dst + x), a[4]); // 256      32    } src += srcStride; dst += dstStride; } } 

नोट: अनुकूलन से अधिकतम त्वरण प्राप्त करने के लिए, AVX2 निर्देशों वाले कोड को निम्नलिखित कंपाइलर विकल्पों के साथ संकलित किया जाना चाहिए: (/ arch: AVX for Visual Studio 2012, -march = core-avx2 for GCC)। इसके अलावा, यह बहुत वांछनीय है कि कोड में वैकल्पिक AVX2 और SSE2 निर्देश शामिल नहीं हैं, क्योंकि इस मामले में प्रोसेसर ऑपरेटिंग मोड को स्विच करने से समग्र प्रदर्शन पर प्रतिकूल प्रभाव पड़ेगा। पूर्वगामी के आधार पर, एल्गोरिदम के AVX2 और SSE2 संस्करणों को अलग-अलग ".cpp" फाइलों में रखने की सलाह दी जाती है।

परीक्षण

परीक्षण 2 एमबी ग्रे छवियों (1920 x 1080) पर किया गया था। समय मापन की सटीकता में सुधार करने के लिए, कई बार परीक्षण चलाए गए। कुल परीक्षण निष्पादन समय को रनों की संख्या से विभाजित करके निष्पादन समय प्राप्त किया गया था। रनों की संख्या को इस तरह से चुना गया था कि कुल निष्पादन का समय कम से कम 1 सेकंड था, जो एल्गोरिदम के निष्पादन समय को मापते समय सटीकता के दो संकेत प्रदान करना चाहिए। एल्गोरिदम को Microsoft Visual Studio 2012 पर 64 बिट विंडोज 7 के लिए अधिकतम अनुकूलन के साथ संकलित किया गया था और iCore-7 4770 प्रोसेसर (3.4 गीगाहर्ट्ज) पर चलाया गया था।

निष्पादन समय, एमएस			सापेक्ष त्वरण
सबसे अच्छा स्केलर	SSE2	AVX2	सबसे अच्छा अदिश / SSE2	सर्वश्रेष्ठ स्केलर / AVX2	SSE2 / AVX2
24,814	0.565	0.424	43,920	58,566	1.333

निष्कर्ष

जैसा कि परीक्षण के परिणामों से देखा जा सकता है, आधुनिक प्रोसेसर पर SIMD निर्देशों के उपयोग से माध्य फ़िल्टरिंग एल्गोरिदम का त्वरण 40 से 60 बार तक पहुंच सकता है। यह मुझे लगता है कि यह अनुकूलन के साथ परेशान करने के लिए एक पर्याप्त राशि है (जब तक कि आपके कार्यक्रम में निष्पादन की गति महत्वपूर्ण नहीं है)। इस प्रकाशन में मैंने यथासंभव उपयोग किए गए कोड को सरल बनाने का प्रयास किया।
इसलिए, डेटा संरेखण, सीमा बिंदुओं के प्रसंस्करण, और बहुत कुछ गुंजाइश से परे बने रहे।
मैं निम्नलिखित लेखों में इन सवालों का खुलासा करने की कोशिश करूंगा। अगर पाठक इस बात में दिलचस्पी रखता है कि इस कार्यक्षमता को लागू करने वाले युद्ध के कोड को किस तरह देखा जाएगा, तो वह इसे यहां ढूंढ सकेगा।