/*!
\file   yof_sort.hpp
\brief  ソートアルゴリズム
\author Kenta Yoshimura

$Id$
*/
/*
Copyright (C) 2006 Kenta Yoshimura All Rights Reserved.

Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions
are met:

1. Redistributions of source code must retain the above copyright
   notice, this list of conditions and the following disclaimer.
2. Redistributions in binary form must reproduce the above copyright
   notice, this list of conditions and the following disclaimer in the
   documentation and/or other materials provided with the distribution.

THIS SOFTWARE IS PROVIDED BY AUTHOR PROJECT ``AS IS'' AND
ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
ARE DISCLAIMED. IN NO EVENT SHALL AUTHOR OR CONTRIBUTORS BE LIABLE
FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
SUCH DAMAGE.
*/

#ifndef YOF_SORT_HPP
#define YOF_SORT_HPP

#if !defined(_MSC_VER) | _MSC_VER > 1000
#pragma once
#endif

#include "yof_types.h"
#include "yof_integer.hpp"
#include <cassert>
#include <cstdlib>
#include <memory>

namespace yoffy {

//! バケットソート
template<typename T>
inline void bucket_sort_u8( const T * first, const T * last, T * result, unsigned offset )
{
	unsigned bucket[256];

	memset( bucket, 0, sizeof(bucket) );

	// 分布を調べる(bucket の中身はヒストグラム)
	for ( uint8_t * itr = offset + (uint8_t *)first; itr < offset + (uint8_t *)last; itr += sizeof(T) )
		bucket[*itr]++;

	// bucket の出力先をテーブル化する(bucket の中身は result の添字テーブル)
	for ( int i = 0, p = 0; i < 256; i++ )
	{
		int tmp = bucket[i];
		bucket[i] = p;
		p += tmp;
	}

	// 出力
	for ( const T * itr = first; itr < last; itr++ )
	{
		const uint8_t index = ((uint8_t *)itr)[offset];

		result[bucket[index]] = *itr;
		bucket[index]++;
	}
}

//! バケットソート(signed の最上位バイト用)
template<typename T>
inline void bucket_sort_s8( const T * first, const T * last, T * result, unsigned offset )
{
	// マイナス方向へ 128 だけシフトしている以外は bucket_sort_u8 と同じ
	unsigned   b[256];
	unsigned * bucket = b;

	memset( b, 0, sizeof(b) );

	bucket += 128;
	for ( int8_t * itr = offset + (int8_t *)first; itr < offset + (int8_t *)last; itr += sizeof(T) )
		bucket[*itr]++;

	for ( int i = -128, p = 0; i < 128; i++ )
	{
		int tmp = bucket[i];
		bucket[i] = p;
		p += tmp;
	}

	for ( const T * itr = first; itr < last; itr++ )
	{
		const int8_t index = ((int8_t *)itr)[offset];

		result[bucket[index]] = *itr;
		bucket[index]++;
	}
}

//! バケットソート(floating の最上位バイト用)
template<typename T>
inline void bucket_sort_float_sign( const T * first, const T * last, T * result, unsigned offset )
{
	unsigned bucket[256];

	memset( bucket, 0, sizeof(bucket) );

	for ( uint8_t * itr = offset + (uint8_t *)first; itr < offset + (uint8_t *)last; itr += sizeof(T) )
		bucket[*itr]++;

	// bucket の出力先をテーブル化する
	// IEEE 754 の構造上、ヒストグラムは 0..FLT_MAX..-0..-FLT_MAX の順に並んでいる為、
	// -FLT_MAX..-0..0..FLT_MAX の順に出力されるようにテーブルを作らなくてはならない
	int p = bucket[255];
	for ( int i = 254; 127 < i; i-- )
	{
		p += bucket[i];
		bucket[i] += bucket[i+1];
	}
	for ( int i = 0; i < 128; i++ )
	{
		int tmp = bucket[i];
		bucket[i] = p;
		p += tmp;
	}

	for ( const T * itr = first; itr < last; itr++ )
	{
		const uint8_t index = ((uint8_t *)itr)[offset];

		if ( 127 < index )
			// マイナス符号が立っているので -0..-FLT_MAX の順に並んでいる
			// なので、後ろから詰めて行かなくてはならない
			result[--bucket[index]] = *itr;
		else
			// プラス符号なので通常通り
			result[bucket[index]++] = *itr;
	}
}

//! 基数ソート(integer)
template<typename T>
inline void radix_sort_helper( T * first, T * last, int_to_type<false> )
{
	const int      size = last - first;
	std::auto_ptr<T> r( new T[size] );
	T *            result = &*r;
	T              *f1 = first,      *l1 = last;
	T              *f2 = &result[0], *l2 = (&result[size-1]) + 1;

#if yof_RT_BIG_ENDIAN
	for ( int i = sizeof(*first) - 1; 0 < i; i-- )
	{
		bucket_sort_u8( f1, l1, f2, i );
		T * tmp = f1; f1 = f2; f2 = tmp;
		tmp = l1; l1 = l2; l2 = tmp;
	}

	if ( 0 < T(-1) )
		bucket_sort_u8( f1, l1, f2, 0 );
	else
		bucket_sort_s8( f1, l1, f2, 0 );
#elif yof_RT_LITTLE_ENDIAN
	for ( int i = 0; i < sizeof(*first) - 1; i++ )
	{
		bucket_sort_u8( f1, l1, f2, i );
		T * tmp = f1; f1 = f2; f2 = tmp;
		tmp = l1; l1 = l2; l2 = tmp;
	}

	if ( 0 < T(-1) )
		bucket_sort_u8( f1, l1, f2, sizeof(*first) - 1 );
	else
		bucket_sort_s8( f1, l1, f2, sizeof(*first) - 1 );
#endif
}

//! 基数ソート(floating)
template<typename T>
inline void radix_sort_helper( T * first, T * last, int_to_type<true> )
{
	const int      size = last - first;
	std::auto_ptr<T> r( new T[size] );
	T *            result = &*r;
	T              *f1 = first,      *l1 = last;
	T              *f2 = &result[0], *l2 = (&result[size-1]) + 1;
	
#if yof_RT_BIG_ENDIAN
	for ( int i = sizeof(*first) - 1; 0 <= i; i++ )
	{
		bucket_sort_u8( f1, l1, f2, i );
		T * tmp = f1; f1 = f2; f2 = tmp;
		tmp = l1; l1 = l2; l2 = tmp;
	}

	bucket_sort_float_sign( f1, l1, f2, 0 );
#elif yof_RT_LITTLE_ENDIAN
	for ( int i = 0; i < sizeof(*first) - 1; i++ )
	{
		bucket_sort_u8( f1, l1, f2, i );
		T * tmp = f1; f1 = f2; f2 = tmp;
		tmp = l1; l1 = l2; l2 = tmp;
	}

	bucket_sort_float_sign( f1, l1, f2, sizeof(*first) - 1 );
#endif
}

/*!
\brief 基数ソート
\param first int や float などの整数型または浮動小数点数型ポインタ
\param last  int や float などの整数型または浮動小数点数型ポインタ
\warning sizeof(*first) が 1 を超え、且つ奇数になる型では使用できません! (sizeof(*first)==1 は可)
\warning 浮動小数点数では、IEEE 754 に準拠していない環境では使用できません!

要素によって計算量が変わらない、安定したソートを行います。\n
計算量は O(kn) であり、中規模〜大規模な配列では qsort や std::sort より高速です。\n
入力は数値のポインタで、ソートは昇順に限定されます。\n
std::vector のような、ポインタを持たないが連続性を保証しているクラスでは、
\code
std::vector<T> v(n);

yoffy::radix_sort( &v[0], (&v[v.size()-1]) + 1 );
\endcode
のようにしてソートを行うことが出来ます。

初期化のコストとして、unsigned[256] の初期化を sizeof(T) 回行います。\n
演算のコストとして、配列全体のサーチとコピーを sizeof(T) 回行います。\n
よって、要素数が少ない場合は他のアルゴリズムの使用を検討してください。
*/
template<typename T>
inline void radix_sort8( T * first, T * last )
{
	assert( first <= last );
	radix_sort_helper( first, last, int_to_type<int_t<T>::is_floating>() );
}

} // namespace yoffy

#endif // #ifndef YOF_SORT_HPP