turbine/flowsepcount.cpp at master · chunmingchen/turbine · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
//C++ standard Header Files
#include <iostream>
#include <cmath>
#include <ctime>
#include <vector>
#include <sstream>
#include <string>
#include <cstring>
#include <fstream>
#include <cstdlib>
#include <stdio.h>
#include <float.h>
#include <stdlib.h>
#include <algorithm>
#include <map>
#include <assert.h>
#include <jclib/system/path.h>
#include <jclib/file/nrrd.h>
#include <jclib/macros.h>

// VTK#include <vtkFloatArray.h>
#include <vtkPointData.h>
#include <vtkVersion.h>
#include <vtkCellArray.h>
#include <vtkPoints.h>
#include <vtkXMLStructuredGridWriter.h>
#include <vtkStructuredGrid.h>
#include <vtkSmartPointer.h>
#include <vtkImageData.h>
#include <vtkXMLImageDataReader.h>
#include <vtkPoints.h>
#include <vtkMultiPieceDataSet.h>
#include <vtkNew.h>
#include <vtkMultiBlockPLOT3DReader.h>
#include <vtkMultiBlockDataSet.h>
#include <vtkXMLPolyDataReader.h>
#include <vtkPolyData.h>
#include <vtkFloatArray.h>
#include <vtkCellArray.h>
#include <vtkQuad.h>
#include <vtkPoints.h>
#include <vtkXMLPolyDataWriter.h>

#include "VectorMatrix.h"


//run: ~/Project/turbine/build/flowsepcount surface.dvt.list
#define EXTRACT_SURFACE

int dim[3]={150, 70, 36}; // -1 from the original sizes
const int xstart=40, xend=105; // end not included
const int ystart=30, yend=61; // hub to tip (end not included)

inline int xyz_id_surface(int x, int y, int b)
{
    return (x-xstart)+(xend-xstart)*(y-ystart + (yend-ystart)*b);
}

void run(char *list_fname)
{
    //char list_filename[1024];
	//sprintf(list_filename, "%s/list", DATA_PATH);
	FILE *fp = fopen(list_fname, "rt");
	char s[1024];
	fgets(s, 1024, fp);

    int files = atoi(s);
    int i;

    FILE *fstat = fopen("flowsepcount.txt", "wt");

    for (i=0; i<files; i++)
    {
        fgets(s, 1024, fp);
        char fname[1024];
        sscanf(s, "%s", fname);

        // load with vtk
        vtkNew<vtkXMLPolyDataReader> reader;
        reader->SetFileName(fname);
        reader->UpdateInformation();
        reader->Update();

        // get data
        vtkSmartPointer<vtkPolyData> data = reader->GetOutput() ;

        // output
        vtkFloatArray* dvt_ary = vtkFloatArray::SafeDownCast(
                    data->GetPointData()->GetArray("d_tangent_vel")
                    );
        if (dvt_ary==NULL) {
            printf("cannot get d_tangent_vel\n");
            exit(1);
        }

        float *p = (float *)dvt_ary->GetVoidPointer(0);
        int n = data->GetNumberOfPoints();

        int x,y,b;
        for (b=0; b<dim[2]; b++)
        {
            vector<float> out_ary(xend-xstart, 0.f);
            int count=0;
            for (y=ystart ;y<yend; y++)
                for (x=xstart; x<xend; x++)
                {
                    int idx = x+dim[0]*(y+dim[1]*b);
                    if (p[idx] < 0)
                        out_ary[x-xstart]  = max(out_ary[x-xstart], -p[idx]);
                }

            for (x=0; x<xend-xstart; x++)
            {
                fprintf(fstat, "%g ", out_ary[x]);
            }
            fprintf(fstat, "\n");
        }

        printf("%d\n", i);
        fflush(fstat);
    }
    fclose(fp);
    fclose(fstat);

}

int main(int argc, char **argv)
{
    printf("usage: filename \n");
    char *filename = argv[1];

	// load data
    run(filename);


}