%% $RCSfile: sdspmdn2.tlc,v $
%% $Revision: 1.2 $ 
%% $Date: 2000/06/15 20:37:44 $
%%
%% Copyright 1995-2000 The MathWorks, Inc. 
%%
%% Abstract: Output the median value of an input vector

%implements sdspmdn2 "C"

%%============================================================================
%openfile buffer
/* DSP Blockset Median Block dependencies: */
extern void qid_RecSort(const real_T *qid_array, int_T *qid_index, int_T i, int_T j);

%closefile buffer
%<LibCacheFunctionPrototype(buffer)> 
%%============================================================================


%% Function: InitializeConditions =============================================
%%
%function InitializeConditions(block, system) Output
  /* DSP Blockset Median (%<ParamSettings.FunctionName>) - %<Name> */
  %%
  %assign INPORT_DATA = 0
  %assign width       = LibBlockInputSignalWidth(INPORT_DATA)
  %assign numDims     = LibBlockInputSignalNumDimensions(INPORT_DATA)
  %assign dims        = LibBlockInputSignalDimensions(INPORT_DATA)  
  %%
  %assign isFullMatrix = LibBlockInputSignalIsFullMatrix(INPORT_DATA)
  %assign frame        = LibBlockInputSignalIsFrameData(INPORT_DATA)
  %assign nChans       = (frame || isFullMatrix) ? dims[1] : 1
  %assign nSamps       = (frame || isFullMatrix) ? dims[0] : width
  %%
  {
	int32_T *index = %<LibBlockDWorkAddr(Index,"","",0)>;
    int_T    i;
    
	for (i=0; i < %<nChans>; i++) {
		int j;
		for (j=0; j < %<nSamps>; j++) {
			*index++ = j;
		}
    }
  }
%endfunction %% InitializeConditions


%% Function: Outputs ==========================================================
%%
%function Outputs(block, system) Output
  /* DSP Blockset Median (%<ParamSettings.FunctionName>) - %<Name> */
  %%
  %assign INPORT_DATA  = 0
  %assign OUTPORT_DATA = 0
  %%
  %assign cplx   = (LibBlockInputSignalIsComplex(INPORT_DATA) != 0)
  %assign dtype  =  cplx ? "creal_T" : "real_T"
  %%
  %assign width   = LibBlockInputSignalWidth(INPORT_DATA)
  %assign numDims = LibBlockInputSignalNumDimensions(INPORT_DATA)
  %assign dims    = LibBlockInputSignalDimensions(INPORT_DATA)  
  %%
  %assign isFullMatrix = LibBlockInputSignalIsFullMatrix(INPORT_DATA)
  %assign frame        = LibBlockInputSignalIsFrameData(INPORT_DATA)
  %assign nChans       = (frame || isFullMatrix) ? dims[1] : 1
  %assign nSamps       = (frame || isFullMatrix) ? dims[0] : width
  %%
  %if width != 1
	{
		const %<dtype>*u       = %<LibBlockInputSignalAddr(INPORT_DATA, "", "", 0)>;
		real_T        *sort_in = %<LibBlockDWorkAddr(Sort_In,"","",0)>;
		int32_T       *index   = %<LibBlockDWorkAddr(Index,"","",0)>;
		int_T          i;
		for(i=0; i<%<nChans>; i++) {
			int j;
			for(j=0;j<%<nSamps>;j++) {
				%if cplx
					sort_in[j] = CMAGSQ(*u);				
					u++;
				%else
					sort_in[j] = *u++;
				%endif
			}
			/* Sort input values */
			qid_RecSort(sort_in, index, 0, %<nSamps>-1);

			sort_in += %<nSamps>;
			index   += %<nSamps>;
		}
	}
  %endif
  %%
  %if nChans > 1
	{
		int32_T  *index = %<LibBlockDWorkAddr(Index,"","",0)>;
		%<dtype> *y     = %<LibBlockOutputSignalAddr(OUTPORT_DATA, "","", 0)>;
		%<dtype> *u     = %<LibBlockInputSignalAddr(INPORT_DATA, "","", 0)>;
		int i; 

		for (i=0;i< %<nChans>;i++) {
			%if (nSamps % 2 == 0) 
				/* Even number of elements - interpolate: */
				%if !cplx
					*y++ = 0.5 * (u[index[%<nSamps>/2-1]] + u[index[%<nSamps>/2]]);
				%else
					const creal_T *u0 = &u[index[%<nSamps>/2-1]];
					const creal_T *u1 = &u[index[%<nSamps>/2]];

					y->re     = 0.5 * (u0->re + u1->re);
					(y++)->im = 0.5 * (u0->im + u1->im);
				%endif
			%else 
				%if !cplx
					*y++ = u[index[(%<nSamps>-1)/2]];
				%else
					*y++ = u[index[(%<nSamps>-1)/2]];
				%endif
			%endif
			u += %<nSamps>;
		}
	}
  %else
	  %if (nSamps % 2 == 0)      
			/* Even number of elements - interpolate: */
			%%
			%if !cplx
				  %% REAL:
				  %assign u0 = LibBlockInputSignal(INPORT_DATA, "%<LibBlockDWork(Index, "%<nSamps/2-1>", "", 0)>", "", "")
				  %assign u1 = LibBlockInputSignal(INPORT_DATA, "%<LibBlockDWork(Index, "%<nSamps/2>", "", 0)>", "", "")
				  %%
				  %<LibBlockOutputSignal(OUTPORT_DATA, "","", 0)> = 0.5 * (%<u0> + %<u1>);
			%else
				  %% COMPLEX:	
				  %assign u0re = LibBlockInputSignal(INPORT_DATA, "%<LibBlockDWork(Index, "%<nSamps/2-1>", "", 0)>", "", "%<tRealPart>")
				  %assign u0im = LibBlockInputSignal(INPORT_DATA, "%<LibBlockDWork(Index, "%<nSamps/2-1>", "", 0)>", "", "%<tImagPart>")
				  %assign u1re = LibBlockInputSignal(INPORT_DATA, "%<LibBlockDWork(Index, "%<nSamps/2>", "", 0)>", "", "%<tRealPart>")
				  %assign u1im = LibBlockInputSignal(INPORT_DATA, "%<LibBlockDWork(Index, "%<nSamps/2>", "", 0)>", "", "%<tImagPart>")
				  %%
				  %<LibBlockOutputSignal(OUTPORT_DATA, "","", "%<tRealPart>0")> = 0.5 * (%<u0re> + %<u1re>);
				  %<LibBlockOutputSignal(OUTPORT_DATA, "","", "%<tImagPart>0")> = 0.5 * (%<u0im> + %<u1im>);
			%endif
			%%      
	  %else
			%%
			%% REAL AND COMPLEX uses same code in this case
			%%
			%<LibBlockOutputSignal(OUTPORT_DATA, "","", 0)> = \
			%<LibBlockInputSignal(INPORT_DATA, "%<LibBlockDWork(Index, "%<(nSamps-1)/2>", "", 0)>","", 0)>; 
	  %endif    
  %endif  
%endfunction

%% [EOF] sdspmdn2.tlc