ArDrone SDK 1.8 added
[mardrone] / mardrone / ARDrone_SDK_Version_1_8_20110726 / ARDroneLib / VP_SDK / Examples / win32 / mjpeg_decoder.c
1 #include <stdlib.h>
2 #include <ctype.h>
3
4 #include <VP_Api/vp_api.h>
5 #include <VP_Api/vp_api_thread_helper.h>
6 #include <VP_Api/vp_api_error.h>
7 #include <VP_Api/vp_api_picture.h>
8 #include <VP_Stages/vp_stages_configs.h>
9 #include <VP_Stages/vp_stages_io_console.h>
10 #include <VP_Stages/vp_stages_o_sdl.h>
11 #include <VP_Stages/vp_stages_io_file.h>
12 #include <VP_Os/vp_os_print.h>
13 #include <VP_Os/vp_os_malloc.h>
14 #include <VP_Os/vp_os_delay.h>
15
16 #include <MJPEG/mjpeg.h>
17 #include <MJPEG/dct.h>
18
19 // #define USE_EULER_ANGLES
20
21 // #define ACQ_WIDTH  (176+0*16)
22 // #define ACQ_HEIGHT (144+0*16)
23
24 // #define QVGA_WIDTH  352
25 // #define QVGA_HEIGHT 288
26
27 #define ACQ_WIDTH  (352)
28 #define ACQ_HEIGHT (288)
29
30 #define NB_STAGES 4
31
32
33 PIPELINE_HANDLE pipeline_handle;
34
35
36 PROTO_THREAD_ROUTINE(app, nomParams);
37 PROTO_THREAD_ROUTINE(dct, nomParams);
38
39 BEGIN_THREAD_TABLE
40   THREAD_TABLE_ENTRY(app, 20)
41   THREAD_TABLE_ENTRY(dct, 20)
42 END_THREAD_TABLE
43
44
45 ///*******************************************************************************************************************///
46
47
48 typedef struct _buffer_to_picture_config_t
49 {
50   vp_api_picture_t* picture;
51
52 } buffer_to_picture_config_t;
53
54 C_RESULT
55 buffer_to_picture_open(buffer_to_picture_config_t *cfg)
56 {
57   return (SUCCESS);
58 }
59
60 C_RESULT
61 buffer_to_picture_transform(buffer_to_picture_config_t *cfg, vp_api_io_data_t *in, vp_api_io_data_t *out)
62 {
63   vp_os_mutex_lock(&out->lock);
64
65
66   if(out->status == VP_API_STATUS_INIT)
67   {
68     out->numBuffers   = 1;
69     out->size         = (ACQ_WIDTH*ACQ_HEIGHT*3)/2;
70     out->buffers      = (int8_t **) cfg->picture;
71     out->indexBuffer  = 0;
72     out->status       = VP_API_STATUS_PROCESSING;
73   }
74
75   if(out->status == VP_API_STATUS_ENDED)
76   {
77   }
78
79   if(out->status == VP_API_STATUS_PROCESSING)
80   {
81     vp_os_memcpy( cfg->picture->y_buf, in->buffers[0], ACQ_WIDTH*ACQ_HEIGHT );
82     vp_os_memcpy( cfg->picture->cb_buf, in->buffers[0] + ACQ_WIDTH*ACQ_HEIGHT, ACQ_WIDTH*ACQ_HEIGHT/4 );
83     vp_os_memcpy( cfg->picture->cr_buf, in->buffers[0] + ACQ_WIDTH*ACQ_HEIGHT + ACQ_WIDTH*ACQ_HEIGHT/4, ACQ_WIDTH*ACQ_HEIGHT/4 );
84   }
85
86   out->status = in->status;
87
88   vp_os_mutex_unlock(&out->lock);
89
90   return (SUCCESS);
91 }
92
93 C_RESULT
94 buffer_to_picture_close(buffer_to_picture_config_t *cfg)
95 {
96   return (SUCCESS);
97 }
98
99 const vp_api_stage_funcs_t buffer_to_picture_funcs =
100 {
101   NULL,
102   (vp_api_stage_open_t)buffer_to_picture_open,
103   (vp_api_stage_transform_t)buffer_to_picture_transform,
104   (vp_api_stage_close_t)buffer_to_picture_close
105 };
106
107
108 ///*******************************************************************************************************************///
109
110
111 typedef struct _picture_to_buffer_config_t
112 {
113   vp_api_picture_t* picture;
114
115 } picture_to_buffer_config_t;
116
117 C_RESULT
118 picture_to_buffer_open(buffer_to_picture_config_t *cfg)
119 {
120   return C_OK;
121 }
122
123 C_RESULT
124 picture_to_buffer_transform(buffer_to_picture_config_t *cfg, vp_api_io_data_t *in, vp_api_io_data_t *out)
125 {
126   vp_os_mutex_lock(&out->lock);
127
128   if(out->status == VP_API_STATUS_INIT)
129   {
130     out->numBuffers   = 1;
131     out->size         = (ACQ_WIDTH*ACQ_HEIGHT*3)/2;
132     out->buffers      = (int8_t **) vp_os_malloc(out->size*sizeof(int8_t) + sizeof(int8_t*));
133     out->indexBuffer  = 0;
134     out->status       = VP_API_STATUS_PROCESSING;
135
136     out->buffers[0]   = (int8_t *)(out->buffers+1);
137   }
138
139   if(out->status == VP_API_STATUS_PROCESSING)
140   {
141     if( in->size == 1 )
142     {
143       // got a picture
144       vp_os_memcpy( out->buffers[0], cfg->picture->y_buf, ACQ_WIDTH*ACQ_HEIGHT );
145       vp_os_memcpy( out->buffers[0] + ACQ_WIDTH*ACQ_HEIGHT, cfg->picture->cb_buf, ACQ_WIDTH*ACQ_HEIGHT/4);
146       vp_os_memcpy( out->buffers[0] + ACQ_WIDTH*ACQ_HEIGHT + ACQ_WIDTH*ACQ_HEIGHT/4, cfg->picture->cr_buf, ACQ_WIDTH*ACQ_HEIGHT/4);
147     }
148   }
149
150   // out->status = in->status;
151
152   vp_os_mutex_unlock(&out->lock);
153
154   return (SUCCESS);
155 }
156
157 C_RESULT
158 picture_to_buffer_close(buffer_to_picture_config_t *cfg)
159 {
160   return (SUCCESS);
161 }
162
163 const vp_api_stage_funcs_t picture_to_buffer_funcs =
164 {
165   NULL,
166   (vp_api_stage_open_t) picture_to_buffer_open,
167   (vp_api_stage_transform_t) picture_to_buffer_transform,
168   (vp_api_stage_close_t) picture_to_buffer_close
169 };
170
171
172 ///*******************************************************************************************************************///
173
174
175 enum {
176   MJPEG_ENCODER,
177   MJPEG_DECODER
178 };
179
180 typedef struct _mjpeg_stage_encoding_config_t
181 {
182   stream_t          stream;
183   mjpeg_t           mjpeg;
184   vp_api_picture_t* picture;
185
186   uint32_t          out_buffer_size;
187
188 } mjpeg_stage_encoding_config_t;
189
190 C_RESULT mjpeg_stage_encoding_open(mjpeg_stage_encoding_config_t *cfg)
191 {
192   stream_new( &cfg->stream, INPUT_STREAM );
193
194   return mjpeg_init( &cfg->mjpeg, MJPEG_ENCODE, cfg->picture->width, cfg->picture->height, cfg->picture->format );
195 }
196
197 C_RESULT mjpeg_stage_encoding_transform(mjpeg_stage_encoding_config_t *cfg, vp_api_io_data_t *in, vp_api_io_data_t *out)
198 {
199   C_RESULT res;
200   uint32_t num_frames;
201   bool_t got_image;
202
203   res = C_OK;
204
205   vp_os_mutex_lock(&out->lock);
206
207   if( out->status == VP_API_STATUS_INIT )
208   {
209     out->numBuffers   = 1;
210     out->buffers      = (int8_t**) vp_os_malloc( sizeof(int8_t*) + cfg->out_buffer_size*sizeof(int8_t) );
211     out->buffers[0]   = (int8_t*) ( out->buffers + 1 );
212     out->indexBuffer  = 0;
213
214     out->status = VP_API_STATUS_PROCESSING;
215   }
216
217   if( out->status == VP_API_STATUS_PROCESSING )
218   {
219     stream_config( &cfg->stream, cfg->out_buffer_size, out->buffers[0] );
220
221     num_frames = cfg->mjpeg.num_frames;
222     res = mjpeg_encode( &cfg->mjpeg, cfg->picture, &cfg->stream, &got_image );
223     if( got_image )
224     {
225       PRINT("Frame complete. Size = %d bytes\n", cfg->stream.index);
226     }
227     out->size = cfg->stream.index;
228   }
229
230   if( out->status == VP_API_STATUS_ENDED )
231   {
232     PRINT("End of data\n");
233   }
234
235   vp_os_mutex_unlock( &out->lock );
236
237   return C_OK;
238 }
239
240 C_RESULT mjpeg_stage_encoding_close(mjpeg_stage_encoding_config_t *cfg)
241 {
242   return mjpeg_release( &cfg->mjpeg );
243 }
244
245
246 ///*******************************************************************************************************************///
247
248
249 typedef struct _mjpeg_stage_decoding_config_t
250 {
251   stream_t          stream;
252   mjpeg_t           mjpeg;
253   vp_api_picture_t* picture;
254
255   uint32_t          out_buffer_size;
256
257 } mjpeg_stage_decoding_config_t;
258
259 C_RESULT mjpeg_stage_decoding_open(mjpeg_stage_decoding_config_t *cfg)
260 {
261   stream_new( &cfg->stream, OUTPUT_STREAM );
262
263   return mjpeg_init( &cfg->mjpeg, MJPEG_DECODE, cfg->picture->width, cfg->picture->height, cfg->picture->format );
264 }
265
266 C_RESULT mjpeg_stage_decoding_transform(mjpeg_stage_decoding_config_t *cfg, vp_api_io_data_t *in, vp_api_io_data_t *out)
267 {
268   bool_t got_image;
269
270   vp_os_mutex_lock( &out->lock );
271
272   if(out->status == VP_API_STATUS_INIT)
273   {
274     out->numBuffers   = 1;
275     out->buffers      = (int8_t**)&cfg->picture;
276     out->indexBuffer  = 0;
277     out->lineSize     = 0;
278
279     out->status = VP_API_STATUS_PROCESSING;
280   }
281
282   if( in->status == VP_API_STATUS_ENDED )
283     out->status = in->status;
284
285   // Several cases must be handled in this stage
286   // 1st: Input buffer is too small to decode a complete picture
287   // 2nd: Input buffer is big enough to decode 1 frame
288   // 3rd: Input buffer is so big we can decode more than 1 frame
289
290   if( out->status == VP_API_STATUS_PROCESSING )
291   {
292     // Reinit stream with new data
293     stream_config( &cfg->stream, in->size, in->buffers[in->indexBuffer] );
294   }
295
296   if(out->status == VP_API_STATUS_PROCESSING || out->status == VP_API_STATUS_STILL_RUNNING)
297   {
298     // If out->size == 1 it means picture is ready
299     out->size = 0;
300     out->status = VP_API_STATUS_PROCESSING;
301
302     mjpeg_decode( &cfg->mjpeg, cfg->picture, &cfg->stream, &got_image );
303
304     if( got_image )
305     {
306       // we got one picture (handle case 1)
307       out->size = 1;
308
309       PRINT( "%d picture decoded\n", cfg->mjpeg.num_frames );
310
311       // handle case 2 & 3
312       if( FAILED(stream_is_empty( &cfg->stream )) )
313       {
314         // Some data are still in stream
315         // Next time we run this stage we don't want this data to be lost
316         // So flag it!
317         out->status = VP_API_STATUS_STILL_RUNNING;
318       }
319     }
320   }
321
322   vp_os_mutex_unlock( &out->lock );
323
324   return C_OK;
325 }
326
327 C_RESULT mjpeg_stage_decoding_close(mjpeg_stage_decoding_config_t *cfg)
328 {
329   stream_delete( &cfg->stream );
330
331   return mjpeg_release( &cfg->mjpeg );
332 }
333
334
335 ///*******************************************************************************************************************///
336
337
338 const vp_api_stage_funcs_t mjpeg_encoding_funcs = {
339   (vp_api_stage_handle_msg_t) NULL,
340   (vp_api_stage_open_t) mjpeg_stage_encoding_open,
341   (vp_api_stage_transform_t) mjpeg_stage_encoding_transform,
342   (vp_api_stage_close_t) mjpeg_stage_encoding_close
343 };
344
345
346 const vp_api_stage_funcs_t mjpeg_decoding_funcs = {
347   (vp_api_stage_handle_msg_t) NULL,
348   (vp_api_stage_open_t) mjpeg_stage_decoding_open,
349   (vp_api_stage_transform_t) mjpeg_stage_decoding_transform,
350   (vp_api_stage_close_t) mjpeg_stage_decoding_close
351 };
352
353
354 ///*******************************************************************************************************************///
355
356 int32_t codec = MJPEG_DECODER;
357
358 int
359 main(int argc, char **argv)
360 {
361   // START_THREAD(escaper, NO_PARAM);
362   START_THREAD(app, 0);
363
364   // JOIN_THREAD(escaper);
365   JOIN_THREAD(app);
366
367   return EXIT_SUCCESS;
368 }
369
370
371 PROTO_THREAD_ROUTINE(app,argv)
372 {
373   uint32_t num_stages = 0;
374   vp_api_picture_t picture;
375
376   vp_api_io_pipeline_t    pipeline;
377   vp_api_io_data_t        out;
378   vp_api_io_stage_t       stages[NB_STAGES];
379
380   vp_stages_input_file_config_t     ifc;
381   vp_stages_output_file_config_t    ofc;
382   // vp_stages_output_sdl_config_t     osc;
383
384   buffer_to_picture_config_t        bpc;
385   mjpeg_stage_encoding_config_t     mec;
386
387   picture_to_buffer_config_t        pbc;
388   mjpeg_stage_decoding_config_t     dec;
389
390   /// Picture configuration
391   picture.format              = PIX_FMT_YUV420P;
392
393   picture.width               = ACQ_WIDTH;
394   picture.height              = ACQ_HEIGHT;
395   picture.framerate           = 15;
396
397   picture.y_buf               = vp_os_malloc( ACQ_WIDTH*ACQ_HEIGHT );
398   picture.cr_buf              = vp_os_malloc( ACQ_WIDTH*ACQ_HEIGHT/4 );
399   picture.cb_buf              = vp_os_malloc( ACQ_WIDTH*ACQ_HEIGHT/4 );
400
401   picture.y_line_size         = ACQ_WIDTH;
402   picture.cb_line_size        = ACQ_WIDTH / 2;
403   picture.cr_line_size        = ACQ_WIDTH / 2;
404
405   picture.y_pad               = 0;
406   picture.c_pad               = 0;
407
408   vp_os_memset(&ifc,0,sizeof(vp_stages_input_file_config_t));
409
410   ifc.name                    = "../temp.mjpg";
411   ifc.buffer_size             = (ACQ_WIDTH*ACQ_HEIGHT*3)/2;
412
413   ofc.name                    = "../out.yuv";
414
415   stages[num_stages].type     = VP_API_INPUT_FILE;
416   stages[num_stages].cfg      = (void *)&ifc;
417   stages[num_stages].funcs    = vp_stages_input_file_funcs;
418
419   num_stages++;
420
421   if( codec == MJPEG_ENCODER )
422   {
423     bpc.picture         = &picture;
424
425     mec.picture         = &picture;
426     mec.out_buffer_size = 4096 * 4;
427
428     stages[num_stages].type      = VP_API_FILTER_DECODER;
429     stages[num_stages].cfg       = (void *)&bpc;
430     stages[num_stages].funcs     = buffer_to_picture_funcs;
431
432     num_stages++;
433
434     stages[num_stages].type      = MJPEG_ENCODER;
435     stages[num_stages].cfg       = (void*)&mec;
436     stages[num_stages].funcs     = mjpeg_encoding_funcs;
437   }
438   else if( codec == MJPEG_DECODER )
439   {
440     dec.picture         = &picture;
441     dec.out_buffer_size = 4096 * 4;
442
443     pbc.picture         = &picture;
444
445     stages[num_stages].type      = MJPEG_DECODER;
446     stages[num_stages].cfg       = (void*)&dec;
447     stages[num_stages].funcs     = mjpeg_decoding_funcs;
448
449     num_stages++;
450
451     stages[num_stages].type      = VP_API_FILTER_ENCODER;
452     stages[num_stages].cfg       = (void *)&pbc;
453     stages[num_stages].funcs     = picture_to_buffer_funcs;
454   }
455
456   num_stages++;
457
458   stages[num_stages].type      = VP_API_OUTPUT_FILE;
459   stages[num_stages].cfg       = (void*)&ofc;
460   stages[num_stages].funcs     = vp_stages_output_file_funcs;
461
462   num_stages++;
463
464   pipeline.nb_stages  = num_stages;
465   pipeline.stages     = &stages[0];
466
467   PRINT("Pipeline configured with %d stages\n", num_stages);
468
469   vp_api_open(&pipeline, &pipeline_handle);
470   out.status = VP_API_STATUS_PROCESSING;
471   while(SUCCEED(vp_api_run(&pipeline, &out)) && (out.status == VP_API_STATUS_PROCESSING || out.status == VP_API_STATUS_STILL_RUNNING));
472
473   vp_api_close(&pipeline, &pipeline_handle);
474
475   return EXIT_SUCCESS;
476 }
477
478
479 ///*******************************************************************************************************************///
480
481
482 // static THREAD_HANDLE dct_thread_handle;
483 static vp_os_mutex_t dct_start_mutex;
484 static vp_os_cond_t  dct_start_cond;
485 static vp_os_mutex_t critical_section;
486
487 static dct_io_buffer_t* current_io_buffer;
488 static dct_io_buffer_t* result_io_buffer;
489
490 static void fdct(const unsigned short* in, short* out);
491 static void idct(const short* in, unsigned short* out);
492
493
494 //-----------------------------------------------------------------------------
495 // DCT Thread
496 //-----------------------------------------------------------------------------
497
498 PROTO_THREAD_ROUTINE(dct, params)
499 {
500   uint32_t i;
501
502   PRINT("DCT thread start\n");
503
504   while(1)
505   {
506     if( current_io_buffer == NULL )
507     {
508       vp_os_mutex_lock(&dct_start_mutex);
509         vp_os_cond_wait(&dct_start_cond);
510       vp_os_mutex_unlock(&dct_start_mutex);
511     }
512
513     if( current_io_buffer->dct_mode == DCT_MODE_FDCT )
514     {
515       for( i = 0; i < current_io_buffer->num_total_blocks; i++ )
516       {
517         fdct(current_io_buffer->input[i], current_io_buffer->output[i]);
518       }
519     }
520     else if( current_io_buffer->dct_mode == DCT_MODE_IDCT )
521     {
522       for( i = 0; i < current_io_buffer->num_total_blocks; i++ )
523       {
524         idct(current_io_buffer->input[i], current_io_buffer->output[i]);
525       }
526     }
527
528     vp_os_mutex_lock(&critical_section);
529       result_io_buffer = current_io_buffer;
530       current_io_buffer = NULL;
531     vp_os_mutex_unlock(&critical_section);
532   }
533
534   return 0;
535 }
536
537
538 //-----------------------------------------------------------------------------
539 // DCT API
540 //-----------------------------------------------------------------------------
541
542
543 bool_t dct_init(void)
544 {
545   vp_os_mutex_init(&dct_start_mutex);
546   vp_os_cond_init(&dct_start_cond, &dct_start_mutex);
547
548   vp_os_mutex_init(&critical_section);
549
550   current_io_buffer = NULL;
551   result_io_buffer  = NULL;
552
553   return TRUE;
554 }
555
556 bool_t dct_compute( dct_io_buffer_t* io_buffer )
557 {
558   bool_t res = FALSE;
559
560   assert(io_buffer != NULL);
561
562   if( current_io_buffer == NULL && result_io_buffer == NULL )
563   {
564     current_io_buffer = io_buffer;
565
566     res = TRUE;
567   }
568
569   return res;
570 }
571
572 dct_io_buffer_t* dct_result( void )
573 {
574   uint32_t i;
575   dct_io_buffer_t* io_buffer;
576
577   io_buffer = NULL;
578
579   if( current_io_buffer != NULL)
580   {
581     if( current_io_buffer->dct_mode == DCT_MODE_FDCT )
582     {
583       for( i = 0; i < current_io_buffer->num_total_blocks; i++ )
584       {
585         fdct(current_io_buffer->input[i], current_io_buffer->output[i]);
586       }
587     }
588     else if( current_io_buffer->dct_mode == DCT_MODE_IDCT )
589     {
590       for( i = 0; i < current_io_buffer->num_total_blocks; i++ )
591       {
592         idct(current_io_buffer->input[i], current_io_buffer->output[i]);
593       }
594     }
595
596     io_buffer = current_io_buffer;
597     current_io_buffer = NULL;
598   }
599
600   return io_buffer;
601 }
602
603 //-----------------------------------------------------------------------------
604 // DCT Computation
605 //-----------------------------------------------------------------------------
606
607
608 #define FIX_0_298631336  ((INT32)  2446)        /* FIX(0.298631336) */
609 #define FIX_0_390180644  ((INT32)  3196)        /* FIX(0.390180644) */
610 #define FIX_0_541196100  ((INT32)  4433)        /* FIX(0.541196100) */
611 #define FIX_0_765366865  ((INT32)  6270)        /* FIX(0.765366865) */
612 #define FIX_0_899976223  ((INT32)  7373)        /* FIX(0.899976223) */
613 #define FIX_1_175875602  ((INT32)  9633)        /* FIX(1.175875602) */
614 #define FIX_1_501321110  ((INT32)  12299)       /* FIX(1.501321110) */
615 #define FIX_1_847759065  ((INT32)  15137)       /* FIX(1.847759065) */
616 #define FIX_1_961570560  ((INT32)  16069)       /* FIX(1.961570560) */
617 #define FIX_2_053119869  ((INT32)  16819)       /* FIX(2.053119869) */
618 #define FIX_2_562915447  ((INT32)  20995)       /* FIX(2.562915447) */
619 #define FIX_3_072711026  ((INT32)  25172)       /* FIX(3.072711026) */
620
621 #define INT32       int
622 #define DCTELEM     int
623 #define DCTSIZE     8
624 #define DCTSIZE2    64
625 #define CONST_BITS  13
626 #define PASS1_BITS  1
627 #define ONE     ((INT32) 1)
628 #define MULTIPLY(var,const)  ((var) * (const))
629 #define DESCALE(x,n)  RIGHT_SHIFT((x) + (ONE << ((n)-1)), n)
630 #define RIGHT_SHIFT(x,shft)     ((x) >> (shft))
631
632 static void fdct(const unsigned short* in, short* out)
633 {
634   INT32 tmp0, tmp1, tmp2, tmp3, tmp4, tmp5, tmp6, tmp7;
635   INT32 tmp10, tmp11, tmp12, tmp13;
636   INT32 z1, z2, z3, z4, z5;
637   int ctr;
638   // SHIFT_TEMPS
639
640   int data[DCTSIZE * DCTSIZE];
641   int i, j;
642   int* dataptr = data;
643
644   for( i = 0; i < DCTSIZE; i++ )
645   {
646     for( j = 0; j < DCTSIZE; j++ )
647     {
648       int temp;
649
650       temp = in[i*DCTSIZE + j];
651       dataptr[i*DCTSIZE + j] = temp;
652     }
653   }
654
655   /* Pass 1: process rows. */
656   /* Note results are scaled up by sqrt(8) compared to a true DCT; */
657   /* furthermore, we scale the results by 2**PASS1_BITS. */
658
659   dataptr = data;
660   for (ctr = DCTSIZE-1; ctr >= 0; ctr--) {
661     tmp0 = dataptr[0] + dataptr[7];
662     tmp7 = dataptr[0] - dataptr[7];
663     tmp1 = dataptr[1] + dataptr[6];
664     tmp6 = dataptr[1] - dataptr[6];
665     tmp2 = dataptr[2] + dataptr[5];
666     tmp5 = dataptr[2] - dataptr[5];
667     tmp3 = dataptr[3] + dataptr[4];
668     tmp4 = dataptr[3] - dataptr[4];
669
670     /* Even part per LL&M figure 1 --- note that published figure is faulty;
671      * rotator "sqrt(2)*c1" should be "sqrt(2)*c6".
672      */
673
674     tmp10 = tmp0 + tmp3;
675     tmp13 = tmp0 - tmp3;
676     tmp11 = tmp1 + tmp2;
677     tmp12 = tmp1 - tmp2;
678
679     dataptr[0] = (DCTELEM) ((tmp10 + tmp11) << PASS1_BITS);
680     dataptr[4] = (DCTELEM) ((tmp10 - tmp11) << PASS1_BITS);
681
682     z1 = MULTIPLY(tmp12 + tmp13, FIX_0_541196100);
683     dataptr[2] = (DCTELEM) DESCALE(z1 + MULTIPLY(tmp13, FIX_0_765366865), CONST_BITS-PASS1_BITS);
684     dataptr[6] = (DCTELEM) DESCALE(z1 + MULTIPLY(tmp12, - FIX_1_847759065), CONST_BITS-PASS1_BITS);
685
686     /* Odd part per figure 8 --- note paper omits factor of sqrt(2).
687      * cK represents cos(K*pi/16).
688      * i0..i3 in the paper are tmp4..tmp7 here.
689      */
690
691     z1 = tmp4 + tmp7;
692     z2 = tmp5 + tmp6;
693     z3 = tmp4 + tmp6;
694     z4 = tmp5 + tmp7;
695     z5 = MULTIPLY(z3 + z4, FIX_1_175875602); /* sqrt(2) * c3 */
696
697     tmp4 = MULTIPLY(tmp4, FIX_0_298631336); /* sqrt(2) * (-c1+c3+c5-c7) */
698     tmp5 = MULTIPLY(tmp5, FIX_2_053119869); /* sqrt(2) * ( c1+c3-c5+c7) */
699     tmp6 = MULTIPLY(tmp6, FIX_3_072711026); /* sqrt(2) * ( c1+c3+c5-c7) */
700     tmp7 = MULTIPLY(tmp7, FIX_1_501321110); /* sqrt(2) * ( c1+c3-c5-c7) */
701     z1 = MULTIPLY(z1, - FIX_0_899976223); /* sqrt(2) * (c7-c3) */
702     z2 = MULTIPLY(z2, - FIX_2_562915447); /* sqrt(2) * (-c1-c3) */
703     z3 = MULTIPLY(z3, - FIX_1_961570560); /* sqrt(2) * (-c3-c5) */
704     z4 = MULTIPLY(z4, - FIX_0_390180644); /* sqrt(2) * (c5-c3) */
705
706     z3 += z5;
707     z4 += z5;
708
709     dataptr[7] = (DCTELEM) DESCALE(tmp4 + z1 + z3, CONST_BITS-PASS1_BITS);
710     dataptr[5] = (DCTELEM) DESCALE(tmp5 + z2 + z4, CONST_BITS-PASS1_BITS);
711     dataptr[3] = (DCTELEM) DESCALE(tmp6 + z2 + z3, CONST_BITS-PASS1_BITS);
712     dataptr[1] = (DCTELEM) DESCALE(tmp7 + z1 + z4, CONST_BITS-PASS1_BITS);
713
714     dataptr += DCTSIZE;         /* advance pointer to next row */
715   }
716
717   /* Pass 2: process columns.
718    * We remove the PASS1_BITS scaling, but leave the results scaled up
719    * by an overall factor of 8.
720    */
721
722   dataptr = data;
723   for (ctr = DCTSIZE-1; ctr >= 0; ctr--) {
724     tmp0 = dataptr[DCTSIZE*0] + dataptr[DCTSIZE*7];
725     tmp7 = dataptr[DCTSIZE*0] - dataptr[DCTSIZE*7];
726     tmp1 = dataptr[DCTSIZE*1] + dataptr[DCTSIZE*6];
727     tmp6 = dataptr[DCTSIZE*1] - dataptr[DCTSIZE*6];
728     tmp2 = dataptr[DCTSIZE*2] + dataptr[DCTSIZE*5];
729     tmp5 = dataptr[DCTSIZE*2] - dataptr[DCTSIZE*5];
730     tmp3 = dataptr[DCTSIZE*3] + dataptr[DCTSIZE*4];
731     tmp4 = dataptr[DCTSIZE*3] - dataptr[DCTSIZE*4];
732
733     /* Even part per LL&M figure 1 --- note that published figure is faulty;
734      * rotator "sqrt(2)*c1" should be "sqrt(2)*c6".
735      */
736
737     tmp10 = tmp0 + tmp3;
738     tmp13 = tmp0 - tmp3;
739     tmp11 = tmp1 + tmp2;
740     tmp12 = tmp1 - tmp2;
741
742     dataptr[DCTSIZE*0] = (DCTELEM) DESCALE(tmp10 + tmp11, PASS1_BITS);
743     dataptr[DCTSIZE*4] = (DCTELEM) DESCALE(tmp10 - tmp11, PASS1_BITS);
744
745     z1 = MULTIPLY(tmp12 + tmp13, FIX_0_541196100);
746     dataptr[DCTSIZE*2] = (DCTELEM) DESCALE(z1 + MULTIPLY(tmp13, FIX_0_765366865), CONST_BITS+PASS1_BITS);
747     dataptr[DCTSIZE*6] = (DCTELEM) DESCALE(z1 + MULTIPLY(tmp12, - FIX_1_847759065), CONST_BITS+PASS1_BITS);
748
749     /* Odd part per figure 8 --- note paper omits factor of sqrt(2).
750      * cK represents cos(K*pi/16).
751      * i0..i3 in the paper are tmp4..tmp7 here.
752      */
753
754     z1 = tmp4 + tmp7;
755     z2 = tmp5 + tmp6;
756     z3 = tmp4 + tmp6;
757     z4 = tmp5 + tmp7;
758     z5 = MULTIPLY(z3 + z4, FIX_1_175875602); /* sqrt(2) * c3 */
759
760     tmp4 = MULTIPLY(tmp4, FIX_0_298631336); /* sqrt(2) * (-c1+c3+c5-c7) */
761     tmp5 = MULTIPLY(tmp5, FIX_2_053119869); /* sqrt(2) * ( c1+c3-c5+c7) */
762     tmp6 = MULTIPLY(tmp6, FIX_3_072711026); /* sqrt(2) * ( c1+c3+c5-c7) */
763     tmp7 = MULTIPLY(tmp7, FIX_1_501321110); /* sqrt(2) * ( c1+c3-c5-c7) */
764     z1 = MULTIPLY(z1, - FIX_0_899976223); /* sqrt(2) * (c7-c3) */
765     z2 = MULTIPLY(z2, - FIX_2_562915447); /* sqrt(2) * (-c1-c3) */
766     z3 = MULTIPLY(z3, - FIX_1_961570560); /* sqrt(2) * (-c3-c5) */
767     z4 = MULTIPLY(z4, - FIX_0_390180644); /* sqrt(2) * (c5-c3) */
768
769     z3 += z5;
770     z4 += z5;
771
772     dataptr[DCTSIZE*7] = (DCTELEM) DESCALE(tmp4 + z1 + z3, CONST_BITS+PASS1_BITS);
773     dataptr[DCTSIZE*5] = (DCTELEM) DESCALE(tmp5 + z2 + z4, CONST_BITS+PASS1_BITS);
774     dataptr[DCTSIZE*3] = (DCTELEM) DESCALE(tmp6 + z2 + z3, CONST_BITS+PASS1_BITS);
775     dataptr[DCTSIZE*1] = (DCTELEM) DESCALE(tmp7 + z1 + z4, CONST_BITS+PASS1_BITS);
776
777     dataptr++;  /* advance pointer to next column */
778   }
779
780   for( i = 0; i < DCTSIZE; i++ )
781     for( j = 0; j < DCTSIZE; j++ )
782       out[i*DCTSIZE + j] = data[i*DCTSIZE + j] >> 3;
783 }
784
785 static void idct(const short* in, unsigned short* out)
786 {
787   INT32 tmp0, tmp1, tmp2, tmp3;
788   INT32 tmp10, tmp11, tmp12, tmp13;
789   INT32 z1, z2, z3, z4, z5;
790   int* wsptr;
791   int* outptr;
792   const short* inptr;
793   int ctr;
794   int workspace[DCTSIZE2];      /* buffers data between passes */
795   int data[DCTSIZE2];
796   // SHIFT_TEMPS
797
798   /* Pass 1: process columns from input, store into work array. */
799   /* Note results are scaled up by sqrt(8) compared to a true IDCT; */
800   /* furthermore, we scale the results by 2**PASS1_BITS. */
801
802   inptr = in;
803   wsptr = workspace;
804   for (ctr = DCTSIZE; ctr > 0; ctr--) {
805     /* Due to quantization, we will usually find that many of the input
806      * coefficients are zero, especially the AC terms.  We can exploit this
807      * by short-circuiting the IDCT calculation for any column in which all
808      * the AC terms are zero.  In that case each output is equal to the
809      * DC coefficient (with scale factor as needed).
810      * With typical images and quantization tables, half or more of the
811      * column DCT calculations can be simplified this way.
812      */
813
814     if( inptr[DCTSIZE*1] == 0 && inptr[DCTSIZE*2] == 0 &&
815         inptr[DCTSIZE*3] == 0 && inptr[DCTSIZE*4] == 0 &&
816         inptr[DCTSIZE*5] == 0 && inptr[DCTSIZE*6] == 0 &&
817         inptr[DCTSIZE*7] == 0 ) {
818       /* AC terms all zero */
819       int dcval = inptr[DCTSIZE*0] << PASS1_BITS;
820
821       wsptr[DCTSIZE*0] = dcval;
822       wsptr[DCTSIZE*1] = dcval;
823       wsptr[DCTSIZE*2] = dcval;
824       wsptr[DCTSIZE*3] = dcval;
825       wsptr[DCTSIZE*4] = dcval;
826       wsptr[DCTSIZE*5] = dcval;
827       wsptr[DCTSIZE*6] = dcval;
828       wsptr[DCTSIZE*7] = dcval;
829
830       inptr++;  /* advance pointers to next column */
831       wsptr++;
832       continue;
833     }
834
835     /* Even part: reverse the even part of the forward DCT. */
836     /* The rotator is sqrt(2)*c(-6). */
837
838     z2 = inptr[DCTSIZE*2];
839     z3 = inptr[DCTSIZE*6];
840
841     z1 = MULTIPLY(z2 + z3, FIX_0_541196100);
842     tmp2 = z1 + MULTIPLY(z3, - FIX_1_847759065);
843     tmp3 = z1 + MULTIPLY(z2, FIX_0_765366865);
844
845     z2 = inptr[DCTSIZE*0];
846     z3 = inptr[DCTSIZE*4];
847
848     tmp0 = (z2 + z3) << CONST_BITS;
849     tmp1 = (z2 - z3) << CONST_BITS;
850
851     tmp10 = tmp0 + tmp3;
852     tmp13 = tmp0 - tmp3;
853     tmp11 = tmp1 + tmp2;
854     tmp12 = tmp1 - tmp2;
855
856     /* Odd part per figure 8; the matrix is unitary and hence its
857      * transpose is its inverse.  i0..i3 are y7,y5,y3,y1 respectively.
858      */
859
860     tmp0 = inptr[DCTSIZE*7];
861     tmp1 = inptr[DCTSIZE*5];
862     tmp2 = inptr[DCTSIZE*3];
863     tmp3 = inptr[DCTSIZE*1];
864
865     z1 = tmp0 + tmp3;
866     z2 = tmp1 + tmp2;
867     z3 = tmp0 + tmp2;
868     z4 = tmp1 + tmp3;
869     z5 = MULTIPLY(z3 + z4, FIX_1_175875602); /* sqrt(2) * c3 */
870
871     tmp0 = MULTIPLY(tmp0, FIX_0_298631336); /* sqrt(2) * (-c1+c3+c5-c7) */
872     tmp1 = MULTIPLY(tmp1, FIX_2_053119869); /* sqrt(2) * ( c1+c3-c5+c7) */
873     tmp2 = MULTIPLY(tmp2, FIX_3_072711026); /* sqrt(2) * ( c1+c3+c5-c7) */
874     tmp3 = MULTIPLY(tmp3, FIX_1_501321110); /* sqrt(2) * ( c1+c3-c5-c7) */
875     z1 = MULTIPLY(z1, - FIX_0_899976223); /* sqrt(2) * (c7-c3) */
876     z2 = MULTIPLY(z2, - FIX_2_562915447); /* sqrt(2) * (-c1-c3) */
877     z3 = MULTIPLY(z3, - FIX_1_961570560); /* sqrt(2) * (-c3-c5) */
878     z4 = MULTIPLY(z4, - FIX_0_390180644); /* sqrt(2) * (c5-c3) */
879
880     z3 += z5;
881     z4 += z5;
882
883     tmp0 += z1 + z3;
884     tmp1 += z2 + z4;
885     tmp2 += z2 + z3;
886     tmp3 += z1 + z4;
887
888     /* Final output stage: inputs are tmp10..tmp13, tmp0..tmp3 */
889
890     wsptr[DCTSIZE*0] = (int) DESCALE(tmp10 + tmp3, CONST_BITS-PASS1_BITS);
891     wsptr[DCTSIZE*7] = (int) DESCALE(tmp10 - tmp3, CONST_BITS-PASS1_BITS);
892     wsptr[DCTSIZE*1] = (int) DESCALE(tmp11 + tmp2, CONST_BITS-PASS1_BITS);
893     wsptr[DCTSIZE*6] = (int) DESCALE(tmp11 - tmp2, CONST_BITS-PASS1_BITS);
894     wsptr[DCTSIZE*2] = (int) DESCALE(tmp12 + tmp1, CONST_BITS-PASS1_BITS);
895     wsptr[DCTSIZE*5] = (int) DESCALE(tmp12 - tmp1, CONST_BITS-PASS1_BITS);
896     wsptr[DCTSIZE*3] = (int) DESCALE(tmp13 + tmp0, CONST_BITS-PASS1_BITS);
897     wsptr[DCTSIZE*4] = (int) DESCALE(tmp13 - tmp0, CONST_BITS-PASS1_BITS);
898
899     inptr++;  /* advance pointers to next column */
900     wsptr++;
901   }
902
903   /* Pass 2: process rows from work array, store into output array. */
904   /* Note that we must descale the results by a factor of 8 == 2**3, */
905   /* and also undo the PASS1_BITS scaling. */
906
907   wsptr = workspace;
908   outptr = data;
909   for (ctr = 0; ctr < DCTSIZE; ctr++) {
910     /* Even part: reverse the even part of the forward DCT. */
911     /* The rotator is sqrt(2)*c(-6). */
912
913     z2 = (INT32) wsptr[2];
914     z3 = (INT32) wsptr[6];
915
916     z1 = MULTIPLY(z2 + z3, FIX_0_541196100);
917     tmp2 = z1 + MULTIPLY(z3, - FIX_1_847759065);
918     tmp3 = z1 + MULTIPLY(z2, FIX_0_765366865);
919
920     tmp0 = ((INT32) wsptr[0] + (INT32) wsptr[4]) << CONST_BITS;
921     tmp1 = ((INT32) wsptr[0] - (INT32) wsptr[4]) << CONST_BITS;
922
923     tmp10 = tmp0 + tmp3;
924     tmp13 = tmp0 - tmp3;
925     tmp11 = tmp1 + tmp2;
926     tmp12 = tmp1 - tmp2;
927
928     /* Odd part per figure 8; the matrix is unitary and hence its
929      * transpose is its inverse.  i0..i3 are y7,y5,y3,y1 respectively.
930      */
931
932     tmp0 = (INT32) wsptr[7];
933     tmp1 = (INT32) wsptr[5];
934     tmp2 = (INT32) wsptr[3];
935     tmp3 = (INT32) wsptr[1];
936
937     z1 = tmp0 + tmp3;
938     z2 = tmp1 + tmp2;
939     z3 = tmp0 + tmp2;
940     z4 = tmp1 + tmp3;
941     z5 = MULTIPLY(z3 + z4, FIX_1_175875602); /* sqrt(2) * c3 */
942
943     tmp0 = MULTIPLY(tmp0, FIX_0_298631336); /* sqrt(2) * (-c1+c3+c5-c7) */
944     tmp1 = MULTIPLY(tmp1, FIX_2_053119869); /* sqrt(2) * ( c1+c3-c5+c7) */
945     tmp2 = MULTIPLY(tmp2, FIX_3_072711026); /* sqrt(2) * ( c1+c3+c5-c7) */
946     tmp3 = MULTIPLY(tmp3, FIX_1_501321110); /* sqrt(2) * ( c1+c3-c5-c7) */
947     z1 = MULTIPLY(z1, - FIX_0_899976223); /* sqrt(2) * (c7-c3) */
948     z2 = MULTIPLY(z2, - FIX_2_562915447); /* sqrt(2) * (-c1-c3) */
949     z3 = MULTIPLY(z3, - FIX_1_961570560); /* sqrt(2) * (-c3-c5) */
950     z4 = MULTIPLY(z4, - FIX_0_390180644); /* sqrt(2) * (c5-c3) */
951
952     z3 += z5;
953     z4 += z5;
954
955     tmp0 += z1 + z3;
956     tmp1 += z2 + z4;
957     tmp2 += z2 + z3;
958     tmp3 += z1 + z4;
959
960     /* Final output stage: inputs are tmp10..tmp13, tmp0..tmp3 */
961
962     outptr[0] = (tmp10 + tmp3) >> ( CONST_BITS+PASS1_BITS+3 );
963     outptr[7] = (tmp10 - tmp3) >> ( CONST_BITS+PASS1_BITS+3 );
964     outptr[1] = (tmp11 + tmp2) >> ( CONST_BITS+PASS1_BITS+3 );
965     outptr[6] = (tmp11 - tmp2) >> ( CONST_BITS+PASS1_BITS+3 );
966     outptr[2] = (tmp12 + tmp1) >> ( CONST_BITS+PASS1_BITS+3 );
967     outptr[5] = (tmp12 - tmp1) >> ( CONST_BITS+PASS1_BITS+3 );
968     outptr[3] = (tmp13 + tmp0) >> ( CONST_BITS+PASS1_BITS+3 );
969     outptr[4] = (tmp13 - tmp0) >> ( CONST_BITS+PASS1_BITS+3 );
970
971     wsptr += DCTSIZE; /* advance pointer to next row */
972     outptr += DCTSIZE;
973   }
974
975   for(ctr = 0; ctr < DCTSIZE2; ctr++)
976     out[ctr] = data[ctr];
977 }