00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024 static char const Id[] = "$Id: writedatapng.c 1246 2008-01-30 18:50:24Z scholz $\n\n";
00025 static char const Td[] = "$Today: " __FILE__ " " __DATE__ " " __TIME__ " $\n\n";
00026
00027 #include "writepng.h"
00028 #include "util/util.h"
00029
00030 static int doinit=1;
00031 static int fieldon=0;
00032 static int pix[ND];
00033 static Mat Afe2sq;
00034 static Vec VIpol;
00037 int WritePNG2init(GridData *gdata)
00038 {
00039 MagparFunctionLogBegin;
00040
00041 int rank,size;
00042 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
00043 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
00044
00045 int nmax;
00046 PetscTruth flg;
00047
00048 PetscReal slice_n[ND];
00049 nmax=ND;
00050 ierr = PetscOptionsGetRealArray(PETSC_NULL,"-slice2_n",slice_n,(PetscInt*)&nmax,&flg);CHKERRQ(ierr);
00051 if (flg!=PETSC_TRUE) {
00052 slice_n[0]=slice_n[1]=0;slice_n[2]=1;
00053 PetscPrintf(PETSC_COMM_WORLD,
00054 "Option -slice2_n not found, using default value: (%g,%g,%g)\n",
00055 slice_n[0],slice_n[1],slice_n[2]
00056 );
00057 }
00058 else if (nmax!=ND) {
00059 SETERRQ(PETSC_ERR_ARG_INCOMP,
00060 "Wrong number of values for option -slice2_n!\n");
00061 }
00062
00063 PetscReal slice_p[ND];
00064 nmax=ND;
00065 ierr = PetscOptionsGetRealArray(PETSC_NULL,"-slice2_p",slice_p,(PetscInt*)&nmax,&flg);CHKERRQ(ierr);
00066 if (flg!=PETSC_TRUE) {
00067 slice_p[0]=slice_p[1]=slice_p[2]=PETSC_MAX;
00068 PetscPrintf(PETSC_COMM_WORLD,
00069 "Option -slice2_p not found, using default value: (%g,%g,%g)\n",
00070 slice_p[0],slice_p[1],slice_p[2]
00071 );
00072
00073 PetscPrintf(PETSC_COMM_WORLD,"PNG output off\n");
00074 MagparFunctionLogReturn(0);
00075 }
00076 else if (nmax!=ND) {
00077 SETERRQ(PETSC_ERR_ARG_INCOMP,
00078 "Wrong number of values for option -slice2_p!\n");
00079 }
00080
00081 int slice_g;
00082 ierr = PetscOptionsGetInt(PETSC_NULL,"-slice2_g",(PetscInt*)&slice_g,&flg);CHKERRQ(ierr);
00083 if (flg!=PETSC_TRUE) {
00084 slice_g=0;
00085 PetscPrintf(PETSC_COMM_WORLD,
00086 "Option -slice2_g not found, using default value: %i\n",
00087 slice_g
00088 );
00089 }
00090 if (PetscAbsInt(slice_g)>gdata->n_prop) {
00091 SETERRQ(PETSC_ERR_ARG_INCOMP,
00092 "slice2_g > number property ids\n");
00093 }
00094 if (slice_g==0) {
00095 PetscPrintf(PETSC_COMM_WORLD,
00096 "slice2_g == 0: drawing magnetization of any volume\n"
00097 );
00098 }
00099 else if (slice_g<0) {
00100 PetscPrintf(PETSC_COMM_WORLD,
00101 "slice2_g < 0: drawing magnetization of any volume except pid %i\n",
00102 slice_g
00103 );
00104 }
00105 else {
00106 PetscPrintf(PETSC_COMM_WORLD,
00107 "slice2_g > 0: drawing only magnetization of volumes with pid %i\n",
00108 slice_g
00109 );
00110 }
00111
00112 int res;
00113 ierr = PetscOptionsGetInt (PETSC_NULL,"-res",(PetscInt*)&res,&flg);CHKERRQ(ierr);
00114 if (flg!=PETSC_TRUE) {
00115 res=50;
00116 PetscPrintf(PETSC_COMM_WORLD,
00117 "Option -res not found, using default value: %i\n",
00118 res
00119 );
00120 }
00121
00122 PetscPrintf(PETSC_COMM_WORLD,"shift and tilt to avoid cutting 'between' elements\n");
00123 #define DIST 1e-12
00124 slice_p[0]+=DIST;slice_n[0]+=DIST;
00125 slice_p[1]+=DIST;slice_n[1]+=DIST;
00126 slice_p[2]+=DIST;slice_n[2]+=DIST;
00127 my_dscal(ND,1.0+DIST,slice_n,1);
00128 my_dscal(ND,1.0+DIST,slice_p,1);
00129
00130 PetscReal *vertxyz2;
00131 PetscReal axes2[ND*ND];
00132 ierr = PetscMalloc(ND*gdata->n_vert*sizeof(PetscReal),&vertxyz2);CHKERRQ(ierr);
00133 ierr = AxesRotation(slice_p,slice_n,axes2,0,gdata->n_vert,gdata->vertxyz,vertxyz2);CHKERRQ(ierr);
00134
00135
00136
00137
00138 slice_g--;
00139 PetscReal bbox[2*ND];
00140 pix[0]=res;
00141 ierr = bbox2(gdata,vertxyz2,1,&slice_g,bbox,pix);CHKERRQ(ierr);
00142
00143
00144 bbox[2]=0;bbox[ND+2]=D_EPS;
00145
00146 ierr = CalAfe2sq(gdata,vertxyz2,bbox,pix,-1,&Afe2sq);CHKERRQ(ierr);
00147 ierr = PetscFree(vertxyz2);CHKERRQ(ierr);
00148
00149 if (Afe2sq==PETSC_NULL) {
00150 MagparFunctionLogReturn(0);
00151 }
00152 fieldon=1;
00153
00154 ierr = VecCreate(PETSC_COMM_WORLD,&VIpol);CHKERRQ(ierr);
00155 if (!rank) {
00156 ierr = VecSetSizes(VIpol,
00157 ND*pix[0]*pix[1],
00158 ND*pix[0]*pix[1]
00159 );CHKERRQ(ierr);
00160 }
00161 else {
00162 ierr = VecSetSizes(VIpol,0,ND*pix[0]*pix[1]);CHKERRQ(ierr);
00163 }
00164 ierr = VecSetFromOptions(VIpol);CHKERRQ(ierr);
00165
00166 MagparFunctionLogReturn(0);
00167 }
00168
00169
00170 int WritePNG2(GridData *gdata)
00171 {
00172 MagparFunctionInfoBegin;
00173
00174 int rank,size;
00175 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
00176 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
00177
00178 if (doinit) {
00179 ierr = WritePNG2init(gdata);CHKERRQ(ierr);
00180 doinit=0;
00181 }
00182 if (!fieldon) {
00183 MagparFunctionInfoReturn(0);
00184 }
00185
00186 ierr = MatMult(Afe2sq,gdata->M,VIpol);CHKERRQ(ierr);
00187 if (rank) {
00188 MagparFunctionInfoReturn(0);
00189 }
00190
00191 png_bytepp image;
00192 ierr = PetscMalloc(pix[1]*sizeof(png_bytep),&image);CHKERRQ(ierr);
00193 ierr = PetscMalloc(3*pix[1]*pix[0]*sizeof(png_byte),&image[0]);CHKERRQ(ierr);
00194
00195 for (int i=1;i<pix[1];i++) {
00196 image[i] = (png_bytep) (image[0]+3*pix[0]*i);
00197 assert(image[i]-image[i-1]==3*pix[0]);
00198 }
00199
00200 PetscReal *a;
00201 VecGetArray(VIpol,&a);
00202
00203 for (int k=0;k<ND;k++) {
00204 int mcomp;
00205
00206 mcomp=k;
00207
00208 for (int i=0;i<pix[0];i++) {
00209 for (int j=0;j<pix[1];j++) {
00210 png_byte red,green,blue;
00211
00212
00213
00214 if (a[ND*(i*pix[1]+pix[1]-1-j)+0]*a[ND*(i*pix[1]+pix[1]-1-j)+0]+
00215 a[ND*(i*pix[1]+pix[1]-1-j)+1]*a[ND*(i*pix[1]+pix[1]-1-j)+1]+
00216 a[ND*(i*pix[1]+pix[1]-1-j)+2]*a[ND*(i*pix[1]+pix[1]-1-j)+2]
00217 < D_EPS
00218 ) {
00219 red=green=blue=255;
00220 }
00221 else {
00222 PetscReal t_val;
00223 switch(mcomp) {
00224 case 0:
00225
00226 t_val=a[ND*(i*pix[1]+pix[1]-1-j)+0];
00227 break;
00228 case 1:
00229
00230 t_val=a[ND*(i*pix[1]+pix[1]-1-j)+1];
00231 break;
00232 case 2:
00233
00234 t_val=a[ND*(i*pix[1]+pix[1]-1-j)+2];
00235 break;
00236 default:
00237 SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"mcomp has illegal value\n");
00238 break;
00239 }
00240
00241
00242
00243
00244
00245
00246
00247
00248
00249
00250 PetscReal t_r,t_g,t_b;
00251 PetscReal t_min=-1.0,t_max=+1.0,t_mid=0.5;
00252
00253 t_b=(t_val < (t_max+t_min)*t_mid) ?
00254 t_val/(-(t_max+t_min)*t_mid+t_min)-
00255 (t_max+t_min)*t_mid/(-(t_max+t_min)*t_mid+t_min) :
00256 0.0;
00257 t_g=(t_val < (t_max+t_min)*t_mid) ?
00258 t_val/(+(t_max+t_min)*t_mid-t_min)-
00259 t_min/(+(t_max+t_min)*t_mid-t_min) :
00260 t_val/(-t_max+(t_max+t_min)*t_mid)-
00261 t_max/(-t_max+(t_max+t_min)*t_mid);
00262 t_r=(t_val > (t_max+t_min)*t_mid) ?
00263 t_val/(t_max-(t_max+t_min)*t_mid)-
00264 (t_max+t_min)*t_mid/(t_max-(t_max+t_min)*t_mid) :
00265 0.0;
00266
00267
00268 if (t_b>1.0) t_b=1.0; if (t_b<0.0) t_b=0.0;
00269 if (t_g>1.0) t_g=1.0; if (t_g<0.0) t_g=0.0;
00270 if (t_r>1.0) t_r=1.0; if (t_r<0.0) t_r=0.0;
00271
00272 red= (png_byte)(t_r*254.0);
00273 green=(png_byte)(t_g*254.0);
00274 blue= (png_byte)(t_b*254.0);
00275 }
00276
00277
00278
00279
00280
00281
00282
00283 image[j][3*i+0]=red;
00284 image[j][3*i+1]=green;
00285 image[j][3*i+2]=blue;
00286 }
00287 }
00288
00289 char fmesh[256],desc[256];
00290 ierr = PetscSNPrintf(fmesh,255,"%s_%01i.%04i.%i.png",gdata->simname,2,gdata->inp,mcomp);CHKERRQ(ierr);
00291 ierr = PetscSNPrintf(desc,255,"t=%g ns, Etot=%g (J/m^3)",gdata->time*gdata->tscale*1e9,gdata->Etot*gdata->escale);CHKERRQ(ierr);
00292 WritePNGfile(fmesh,pix[0],pix[1],image,desc);
00293 }
00294
00295 VecRestoreArray(VIpol,&a);
00296
00297 ierr = PetscFree(image[0]);CHKERRQ(ierr);
00298 ierr = PetscFree(image);CHKERRQ(ierr);
00299
00300
00301 MagparFunctionInfoReturn(0);
00302 }
00303