#include "netcdfcpp.h"
#include<math.h>
#include <iostream>
using namespace std;
#define  lon 1442
#define  lat 281
#define  lat2 562

int  main()
{    

    //读取现有矢量场.nc文件,长、宽
    NcFile dataReadFile1("global10.nc",NcFile::Replace);
    if (!dataReadFile1.is_valid())
    {
        std::cout<<"couldn‘t open file!"<<std::endl;
    }
    //添加维度
    dataReadFile1.add_dim("LON",lon);
    dataReadFile1.add_dim("LAT",lat);
    dataReadFile1.add_dim("Time",1);
    dataReadFile1.add_dim("bnds",2);
    dataReadFile1.add_dim("UX",1*lon*lat);
    dataReadFile1.add_dim("VY",1*lon*lat);
    dataReadFile1.add_dim("MAG",1*lon*lat);


    //添加变量信息
    dataReadFile1.add_var("LON",ncDouble,dataReadFile1.get_dim(0));
    dataReadFile1.add_var("LAT",ncDouble,dataReadFile1.get_dim(1));    
    dataReadFile1.add_var("Time",ncDouble,dataReadFile1.get_dim(2));
    dataReadFile1.add_var("Time_bnds",ncDouble,dataReadFile1.get_dim(2),dataReadFile1.get_dim(3));
    dataReadFile1.add_var("UX",ncFloat,dataReadFile1.get_dim(2),dataReadFile1.get_dim(1),dataReadFile1.get_dim(0));
    dataReadFile1.add_var("VY",ncFloat,dataReadFile1.get_dim(2),dataReadFile1.get_dim(1),dataReadFile1.get_dim(0));
    dataReadFile1.add_var("MAG",ncFloat,dataReadFile1.get_dim(2),dataReadFile1.get_dim(1),dataReadFile1.get_dim(0));        

    
    //添加属性
    dataReadFile1.add_att("history","FERRET V6.05   20-Oct-14");
    dataReadFile1.get_var(0)->add_att("long_name","Longitude");
    dataReadFile1.get_var(0)->add_att("point_spacing","even");
    dataReadFile1.get_var(0)->add_att("units","degrees");
    dataReadFile1.get_var(0)->add_att("axis","X");
    dataReadFile1.get_var(0)->add_att("modulo","360");

    dataReadFile1.get_var(1)->add_att("long_name","Latitude");
    dataReadFile1.get_var(1)->add_att("point_spacing","even");
    dataReadFile1.get_var(1)->add_att("units","degrees");
    dataReadFile1.get_var(1)->add_att("axis","Y");
    
    dataReadFile1.get_var(2)->add_att("units","days since 1800-01-01 00:00:00");
    dataReadFile1.get_var(2)->add_att("long_name","Time");
    dataReadFile1.get_var(2)->add_att("axis","T");
    dataReadFile1.get_var(2)->add_att("time_origin","01-JAN-1800 00:00:00");


    //dataReadFile1.get_var(4)->add_att("missing_value",-3.39999995214436E+38);
//    dataReadFile1.get_var(4)->add_att("_FillValue",-3.39999995214436E+38);
    dataReadFile1.get_var(6)->add_att("long_name","wind speed");
    dataReadFile1.get_var(6)->add_att("units","meter per second");
    dataReadFile1.get_var(6)->add_att("history","From monthly");

    //    dataReadFile1.add_att("UX","cFictional Model Output");
    dataReadFile1.get_var(4)->add_att("long_name","Zonal Wind Speed");
    dataReadFile1.get_var(4)->add_att("units","meter second-1");    
    //dataReadFile1.get_var(5)->add_att("missing_value",-3.39999995214436E+38);
    //dataReadFile1.get_var(5)->add_att("_FillValue",-3.39999995214436E+38);
    dataReadFile1.get_var(4)->add_att("history","From monthly");

    //dataReadFile1.add_att("UY","cFictional Model Output");
    dataReadFile1.get_var(5)->add_att("long_name","Meridional Wind Speed");
    dataReadFile1.get_var(5)->add_att("units","meter second-1");
//    dataReadFile1.get_var(6)->add_att("missing_value",-3.39999995214436E+38);
//    dataReadFile1.get_var(6)->add_att("_FillValue",-3.39999995214436E+38);
    dataReadFile1.get_var(5)->add_att("history","From monthly");

    //输出.nc文件
    

    double *Lon = new double[lon];
    double *Lat = new double[lat];
    double *Lat2 = new double[lat2];
    double temp = 0.125;
    Lon[0] = -180;
    Lat[0] = 70;
    
    for (int i = 1;i<lon;i++)
    {
        Lon[i] = Lon[i-1]+temp;
    }
    dataReadFile1.get_var("LON")->put(Lon,lon);

    for (int j = 1;j<lat;j++)
    {
        Lat[j] = Lat[j-1]-0.25;
    }
    dataReadFile1.get_var("LAT")->put(Lat,lat);

    //线性插值，向UX里写数据

    int LonNum=0;
    int LatNum=0;
    //读取现有矢量场.nc文件,长、宽
    NcFile dataReadFile("global_origin.nc",NcFile::ReadOnly);
    if (!dataReadFile.is_valid())
    {
        std::cout<<"couldn‘t open file!"<<std::endl;
    }

    for (int i = 0;i<=dataReadFile.num_dims()-1;i++)
    {
        if (i==0)
        {
            LonNum  =dataReadFile.get_dim(i)->size();
        }
        else if (i==1)
        {
            LatNum=dataReadFile.get_dim(i)->size();
        }
    }
    
    //    double *VY = new double[lon*lat];
    long count1[3];
    count1[0] = 1;
    count1[1] = LatNum;
    count1[2] = LonNum*2-1;



    static const int Time = 1;
    static const int TMP = Time*LonNum*LatNum;
    double *UXValue = http://www.mamicode.com/new double[lon*lat];
    double *VYValue = http://www.mamicode.com/new double[lon*lat];
    double *MAGValue = http://www.mamicode.com/new double[lon*lat];

    double *Tmp_UX = new double[TMP];
    double *Tmp_VY = new double[TMP];
    double *Tmp_LAT = new double[TMP];
    double *Tmp_LON = new double[TMP];

    NcVar *dataTmp_LAT = dataReadFile.get_var("LAT");    
    NcVar *dataTmp_LON = dataReadFile.get_var("LONN359_361");    
    NcVar *dataTmp_UX = dataReadFile.get_var("UX");    
    NcVar *dataTmp_VY = dataReadFile.get_var("VY");    
    dataTmp_LAT->get(Tmp_LAT,LatNum,LatNum);
    dataTmp_LON->get(Tmp_LON,LonNum,LonNum);
    dataTmp_UX->get(Tmp_UX,Time,LatNum,LonNum);
    dataTmp_VY->get(Tmp_VY,Time,LatNum,LonNum);

    long k = 0;

    for (int j =0;j<LatNum;j++)
    {

        for (int i=0;i<LonNum;i++)
        {        
            UXValue[k]=Tmp_UX[j*LonNum+i];
            UXValue[k+1] = (Tmp_UX[j*LonNum+i]+Tmp_UX[j*LonNum+i+1])/2;

            VYValue[k] = Tmp_VY[j*LonNum+i];
            VYValue[k+1] = (Tmp_VY[j*LonNum+i]+Tmp_VY[j*LonNum+i+1])/2;

            MAGValue[k] =sqrt(UXValue[k]*UXValue[k]+VYValue[k]*VYValue[k]);
            MAGValue[k+1] = sqrt(UXValue[k+1]*UXValue[k+1]+VYValue[k+1]*VYValue[k+1]);

            k=k+2;        
        }
    }    
    dataReadFile1.get_var("UX")->put(UXValue,count1);
    dataReadFile1.get_var("VY")->put(VYValue,count1);
    //dataReadFile1.get_var("MAG")->put(MAGValue,count1);

    //纵向扩大2倍
    NcFile dataReadFile2("global15.nc",NcFile::Replace);
    //添加维度
    dataReadFile2.add_dim("LON",lon);
    dataReadFile2.add_dim("LAT",lat2);
    dataReadFile2.add_dim("Time",1);
    dataReadFile2.add_dim("bnds",2);
    dataReadFile2.add_dim("UX",1*lon*lat2);
    dataReadFile2.add_dim("VY",1*lon*lat2);
    dataReadFile2.add_dim("MAG",1*lon*lat2);
    
    //添加变量信息
    dataReadFile2.add_var("LON",ncDouble,dataReadFile2.get_dim(0));
    dataReadFile2.add_var("LAT",ncDouble,dataReadFile2.get_dim(1));    
    dataReadFile2.add_var("Time",ncDouble,dataReadFile2.get_dim(2));
    dataReadFile2.add_var("Time_bnds",ncDouble,dataReadFile2.get_dim(2),dataReadFile2.get_dim(3));
    dataReadFile2.add_var("UX",ncFloat,dataReadFile2.get_dim(2),dataReadFile2.get_dim(1),dataReadFile2.get_dim(0));
    dataReadFile2.add_var("VY",ncFloat,dataReadFile2.get_dim(2),dataReadFile2.get_dim(1),dataReadFile2.get_dim(0));
    dataReadFile2.add_var("MAG",ncFloat,dataReadFile2.get_dim(2),dataReadFile2.get_dim(1),dataReadFile2.get_dim(0));        


    //添加属性
    dataReadFile2.add_att("history","FERRET V6.05   20-Oct-14");
    dataReadFile2.get_var(0)->add_att("long_name","Longitude");
    dataReadFile2.get_var(0)->add_att("point_spacing","even");
    dataReadFile2.get_var(0)->add_att("units","degrees");
    dataReadFile2.get_var(0)->add_att("axis","X");
    dataReadFile2.get_var(0)->add_att("modulo","360");

    dataReadFile2.get_var(1)->add_att("long_name","Latitude");
    dataReadFile2.get_var(1)->add_att("point_spacing","even");
    dataReadFile2.get_var(1)->add_att("units","degrees");
    dataReadFile2.get_var(1)->add_att("axis","Y");

    dataReadFile2.get_var(2)->add_att("units","days since 1800-01-01 00:00:00");
    dataReadFile2.get_var(2)->add_att("long_name","Time");
    dataReadFile2.get_var(2)->add_att("axis","T");
    dataReadFile2.get_var(2)->add_att("time_origin","01-JAN-1800 00:00:00");


    //dataReadFile1.get_var(4)->add_att("missing_value",-3.39999995214436E+38);
    //    dataReadFile1.get_var(4)->add_att("_FillValue",-3.39999995214436E+38);
    dataReadFile2.get_var(6)->add_att("long_name","wind speed");
    dataReadFile2.get_var(6)->add_att("units","meter per second");
    dataReadFile2.get_var(6)->add_att("history","From monthly");

    //    dataReadFile1.add_att("UX","cFictional Model Output");
    dataReadFile2.get_var(4)->add_att("long_name","Zonal Wind Speed");
    dataReadFile2.get_var(4)->add_att("units","meter second-1");    
    //dataReadFile1.get_var(5)->add_att("missing_value",-3.39999995214436E+38);
    //dataReadFile1.get_var(5)->add_att("_FillValue",-3.39999995214436E+38);
    dataReadFile2.get_var(4)->add_att("history","From monthly");

    //dataReadFile1.add_att("UY","cFictional Model Output");
    dataReadFile2.get_var(5)->add_att("long_name","Meridional Wind Speed");
    dataReadFile2.get_var(5)->add_att("units","meter second-1");
    //    dataReadFile1.get_var(6)->add_att("missing_value",-3.39999995214436E+38);
    //    dataReadFile1.get_var(6)->add_att("_FillValue",-3.39999995214436E+38);
    dataReadFile2.get_var(5)->add_att("history","From monthly");
    Lat2[0] = 70;
    for (int i = 1;i<lon;i++)
    {
        Lon[i] = Lon[i-1]+temp;
    }
    dataReadFile2.get_var("LON")->put(Lon,lon);

    for (int j = 1;j<lat2;j++)
    {
        Lat2[j] = Lat2[j-1]-0.125;
        cout<<j<<"  "<<Lat2[j]<<endl;
    }
    dataReadFile2.get_var("LAT")->put(Lat2,lat2);

    int LonNum2 = 0;
    int LatNum2 = 0;
    for (int i = 0;i<=dataReadFile1.num_dims()-1;i++)
    {
        if (i==0)
        {
             LonNum2  =dataReadFile1.get_dim(i)->size();
        }
        else if (i==1)
        {
            LatNum2 =dataReadFile1.get_dim(i)->size();
        }
    }
//     cout<<LonNum2<<endl;
//     cout<<LatNum2<<endl;
    long count2[3];
    count2[0] = 1;
    count2[1] = LatNum*2;
    count2[2] = LonNum*2;

    static const int TMP2 = Time*LonNum2*LatNum2;
    double *UXValue2= new double[lon*lat2];
    double *VYValue2= new double[lon*lat2];
    double *MAGValue2 = new double[lon*lat2];

    //cout<<lon*lat2<<endl;
    double *Tmp_UX2 = new double[TMP2];
    double *Tmp_VY2 = new double[TMP2];
    double *Tmp_LAT2 = new double[TMP2];
    double *Tmp_LON2 = new double[TMP2];


    NcVar *dataTmp_LAT2= dataReadFile1.get_var("LAT");    
    NcVar *dataTmp_LON2 = dataReadFile1.get_var("LON");    
    NcVar *dataTmp_UX2 = dataReadFile1.get_var("UX");    
    NcVar *dataTmp_VY2 = dataReadFile1.get_var("VY");    
    dataTmp_LAT2->get(Tmp_LAT2,LatNum2,LatNum2);
    dataTmp_LON2->get(Tmp_LON2,LonNum2,LonNum2);
    dataTmp_UX2->get(Tmp_UX2,Time,LatNum2,LonNum2);
    dataTmp_VY2->get(Tmp_VY2,Time,LatNum2,LonNum2);

    int k2= 0;
    for (int m = 0;m<LatNum2;m++)
    {
        for (int n = 0;n<LonNum2;n++)
        {
            
             if (k2<=558)
             {
                for (int p = 0;p<LonNum2;p++)
                {
                    UXValue2[k2*LonNum2+p] = Tmp_UX2[m*LonNum2+p];
                }
                UXValue2[(k2+1)*(LonNum2-1)+n] = (Tmp_UX2[m*LonNum2+n]+Tmp_UX2[(m+1)*LonNum2+n])/2;

            }
            else{
                for (int p = 0;p<LonNum2;p++)
                {
                    UXValue2[(k2)*LonNum2+p] = Tmp_UX2[m*LonNum2+p];
                }
            }
            
//cout<<(k2+1)*LonNum2+n<<endl;
            //VYValue[k] = Tmp_VY[j*LonNum+i];
        //    VYValue[k+1] = (Tmp_VY[j*LonNum+i]+Tmp_VY[j*LonNum+i+1])/2;

        //    MAGValue[k] =sqrt(UXValue[k]*UXValue[k]+VYValue[k]*VYValue[k]);
        //    MAGValue[k+1] = sqrt(UXValue[k+1]*UXValue[k+1]+VYValue[k+1]*VYValue[k+1]);

        }
        k2=k2+2;        
        
    }


//     for (int j2 =0;j2<LatNum2;j2++)
//     {
// 
//         for (int i2=0;i2<LonNum2;i2++)
//         {        
//             UXValue2[m*LonNum2+i2]=Tmp_UX2[j2*LonNum2+i2];
//             UXValue2[(m+1)*LonNum2+i2] = (Tmp_UX2[j2*LonNum2+i2]+Tmp_UX2[(j2+1)*LonNum2+i2])/2;
// 
//             //VYValue2[m]=Tmp_VY2[j2*LonNum2+i2];
//         //    VYValue2[(m+1)*LonNum2] = (Tmp_VY2[j2*LonNum2+i2]+Tmp_VY2[(j2+1)*LonNum2+i2])/2;
// 
//             //MAGValue2[k] =sqrt(UXValue[k]*UXValue[k]+VYValue[k]*VYValue[k]);
//             //MAGValue2[k+1] = sqrt(UXValue[k+1]*UXValue[k+1]+VYValue[k+1]*VYValue[k+1]);
// 
//         //        
//         }
//         m=m++;    
//     }    
    dataReadFile2.get_var("UX")->put(UXValue2,count2);
    //dataReadFile1.get_var("VY")->put(VYValue2,count2);
    //dataReadFile1.get_var("MAG")->put(MAGValue,count2);


    return 0;
}