+ case OP_MEDIAN:
+ stackunderflow(0);
+ {
+ int elements = (int) rpnstack->s[stptr--];
+ int final_elements = elements;
+ double *element_ptr = rpnstack->s + stptr - elements + 1;
+ double *goodvals = element_ptr;
+ double *badvals = element_ptr + elements - 1;
+
+ stackunderflow(elements - 1);
+
+ /* move values to consolidate the non-NANs for sorting, keeping
+ * track of how many NANs we encounter. */
+ while (goodvals < badvals) {
+ if (isnan(*goodvals)) {
+ *goodvals = *badvals--;
+ --final_elements;
+ } else {
+ ++goodvals;
+ }
+ }
+
+ stptr -= elements;
+ if (!final_elements) {
+ /* no non-NAN elements; push NAN */
+ rpnstack->s[++stptr] = DNAN;
+ } else {
+ /* when goodvals and badvals meet, they might have met on a
+ * NAN, which wouldn't decrease final_elements. so, check
+ * that now. */
+ if (isnan(*goodvals)) --final_elements;
+ /* and finally, take the median of the remaining non-NAN
+ * elements. */
+ qsort(element_ptr, final_elements, sizeof(double),
+ rpn_compare_double);
+ if (final_elements % 2 == 1){
+ rpnstack->s[++stptr] = element_ptr[ final_elements / 2 ];
+ }
+ else {
+ rpnstack->s[++stptr] = 0.5 * ( element_ptr[ final_elements / 2 ] + element_ptr[ final_elements / 2 - 1 ] );
+ }
+ }
+ }
+ break;
+