//#####################################################################
// Copyright 2011, Alexey Stomakhin.
// This file is part of PhysBAM whose distribution is governed by the license contained in the accompanying file PHYSBAM_COPYRIGHT.txt.
//#####################################################################
// Class BOWL
//##################################################################### 
#ifndef __BOWL__
#define __BOWL__
#include <Core/Basic_Geometry/AUTODIFF_LEVELSET.h>
#include <Core/Meshes/MESHES_FORWARD.h>
namespace PhysBAM{

template<class T>
class BOWL:public AUTODIFF_LEVELSET<VECTOR<T,3>,BOWL<T> >
{
    typedef VECTOR<T,3> TV;

public:
    typedef TV VECTOR_T;

    TV bottom,axis;
    T bottom_radius,inner_radius,thickness;

    BOWL()
        :axis(0,1,0),bottom_radius(2),inner_radius(1),thickness(.5)
    {}

    // The middle of the inside of the bowl is at bottom+axis*thickness.
    BOWL(const TV& bottom,const TV& axis,T bottom_radius,T inner_radius,
        T thickness)
        :bottom(bottom),axis(axis.Normalized()),bottom_radius(bottom_radius),
        inner_radius(inner_radius),thickness(thickness)
    {
    }

    template<class IN>
    auto Raw_Phi(const IN& X) const
    {
        auto Z=X-bottom;
        auto outer_radius=thickness+inner_radius;
        auto y=Z.Dot(axis);
        auto r=(Z-y*axis).Magnitude();
        auto cr=r-bottom_radius;
        auto top_rim=y-outer_radius;
        decltype(r) p;

        if(cr.x<0)
        {
            if(y.x<thickness/2)
            {
                auto out_bot=-y;
                p.Fill_From(out_bot);
            }
            else
            {
                auto in_bot=y-thickness;
                auto in_top_rim=hypot(top_rim,cr-inner_radius);
                p.Fill_From(min(in_bot,in_top_rim));
            }
        }
        else if(top_rim.x<0)
        {
            auto raw_curve=hypot(top_rim,cr);
            auto in_curve=inner_radius-raw_curve;
            auto out_curve=raw_curve-outer_radius;
            p.Fill_From(max(in_curve,out_curve));
        }
        else if(cr.x<inner_radius)
        {
            auto in_top_rim=hypot(top_rim,cr-inner_radius);
            p.Fill_From(in_top_rim);
        }
        else if(cr.x>outer_radius)
        {
            auto out_top_rim=hypot(top_rim,cr-outer_radius);
            p.Fill_From(out_top_rim);
        }
        else p.Fill_From(top_rim);
        return p;
    }

//#####################################################################
    RANGE<TV> Bounding_Box() const;
    static std::string Name();
    template<class RW> void Read(UNTYPED_ISTREAM<RW>& input);
    template<class RW> void Write(UNTYPED_OSTREAM<RW>& output) const;
    void Fill_Uniform(RANDOM_NUMBERS<T>& rand,const RANGE<TV>& box);
//#####################################################################
};   
template<class TV> class RAY;
template<class TV> class RAY_HIT;

template<class T,class TV> BOWL<T> Scaled(const BOWL<T>& b,const TV& center,T scale);
template<class T,class TV> BOWL<T> Transformed(const BOWL<T>& b,const FRAME<TV>& f);
template<class T> TRIANGULATED_SURFACE<T>* Generate_Triangles(const BOWL<T>& bowl,const int n_radial=32,const int n_vertical=8);
template<class T> RAY_HIT<T> Intersects(const RAY<VECTOR<T,3> >& ray,const BOWL<T>& b,T isocontour=0);
template<class T> VECTOR<RAY_HIT<T>,4>
Intersection_Interval(const RAY<VECTOR<T,3> >& ray,const BOWL<T>& b,T isocontour=0);
template<class TV,class T> void
All_Intersections(ARRAY<RAY_HIT<T> >& hits,const RAY<TV>& ray,const BOWL<T>& b,const T isocontour=0);
}
#endif
