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 slice_g--;
00137 PetscReal bbox[2*ND];
00138 pix[0]=res;
00139 ierr = bbox2(gdata,vertxyz2,1,&slice_g,bbox,pix);CHKERRQ(ierr);
00140
00141
00142 bbox[2]=0;bbox[ND+2]=D_EPS;
00143
00144 ierr = CalAfe2sq(gdata,vertxyz2,bbox,pix,slice_g,&Afe2sq);CHKERRQ(ierr);
00145 ierr = PetscFree(vertxyz2);CHKERRQ(ierr);
00146
00147 if (Afe2sq==PETSC_NULL) {
00148 MagparFunctionLogReturn(0);
00149 }
00150 fieldon=1;
00151
00152 ierr = VecCreate(PETSC_COMM_WORLD,&VIpol);CHKERRQ(ierr);
00153 if (!rank) {
00154 ierr = VecSetSizes(VIpol,
00155 ND*pix[0]*pix[1],
00156 ND*pix[0]*pix[1]
00157 );CHKERRQ(ierr);
00158 }
00159 else {
00160 ierr = VecSetSizes(VIpol,0,ND*pix[0]*pix[1]);CHKERRQ(ierr);
00161 }
00162 ierr = VecSetFromOptions(VIpol);CHKERRQ(ierr);
00163
00164 MagparFunctionLogReturn(0);
00165 }
00166
00167
00168 int WritePNG2(GridData *gdata)
00169 {
00170 MagparFunctionInfoBegin;
00171
00172 int rank,size;
00173 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
00174 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
00175
00176 if (doinit) {
00177 ierr = WritePNG2init(gdata);CHKERRQ(ierr);
00178 doinit=0;
00179 }
00180 if (!fieldon) {
00181 MagparFunctionInfoReturn(0);
00182 }
00183
00184 ierr = MatMult(Afe2sq,gdata->M,VIpol);CHKERRQ(ierr);
00185 if (rank) {
00186 MagparFunctionInfoReturn(0);
00187 }
00188
00189 png_bytepp image;
00190 ierr = PetscMalloc(pix[1]*sizeof(png_bytep),&image);CHKERRQ(ierr);
00191 ierr = PetscMalloc(3*pix[1]*pix[0]*sizeof(png_byte),&image[0]);CHKERRQ(ierr);
00192
00193 for (int i=1;i<pix[1];i++) {
00194 image[i] = (png_bytep) (image[0]+3*pix[0]*i);
00195 assert(image[i]-image[i-1]==3*pix[0]);
00196 }
00197
00198 PetscReal *a;
00199 VecGetArray(VIpol,&a);
00200
00201 for (int k=0;k<ND;k++) {
00202 int mcomp;
00203
00204 mcomp=k;
00205
00206 for (int i=0;i<pix[0];i++) {
00207 for (int j=0;j<pix[1];j++) {
00208 png_byte red,green,blue;
00209
00210
00211
00212 if (a[ND*(i*pix[1]+pix[1]-1-j)+0]*a[ND*(i*pix[1]+pix[1]-1-j)+0]+
00213 a[ND*(i*pix[1]+pix[1]-1-j)+1]*a[ND*(i*pix[1]+pix[1]-1-j)+1]+
00214 a[ND*(i*pix[1]+pix[1]-1-j)+2]*a[ND*(i*pix[1]+pix[1]-1-j)+2]
00215 < D_EPS
00216 ) {
00217 red=green=blue=255;
00218 }
00219 else {
00220 PetscReal t_val;
00221 switch(mcomp) {
00222 case 0:
00223
00224 t_val=a[ND*(i*pix[1]+pix[1]-1-j)+0];
00225 break;
00226 case 1:
00227
00228 t_val=a[ND*(i*pix[1]+pix[1]-1-j)+1];
00229 break;
00230 case 2:
00231
00232 t_val=a[ND*(i*pix[1]+pix[1]-1-j)+2];
00233 break;
00234 default:
00235 SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"mcomp has illegal value\n");
00236 break;
00237 }
00238
00239
00240
00241
00242
00243
00244
00245
00246
00247
00248 PetscReal t_r,t_g,t_b;
00249 PetscReal t_min=-1.0,t_max=+1.0,t_mid=0.5;
00250
00251 t_b=(t_val < (t_max+t_min)*t_mid) ?
00252 t_val/(-(t_max+t_min)*t_mid+t_min)-
00253 (t_max+t_min)*t_mid/(-(t_max+t_min)*t_mid+t_min) :
00254 0.0;
00255 t_g=(t_val < (t_max+t_min)*t_mid) ?
00256 t_val/(+(t_max+t_min)*t_mid-t_min)-
00257 t_min/(+(t_max+t_min)*t_mid-t_min) :
00258 t_val/(-t_max+(t_max+t_min)*t_mid)-
00259 t_max/(-t_max+(t_max+t_min)*t_mid);
00260 t_r=(t_val > (t_max+t_min)*t_mid) ?
00261 t_val/(t_max-(t_max+t_min)*t_mid)-
00262 (t_max+t_min)*t_mid/(t_max-(t_max+t_min)*t_mid) :
00263 0.0;
00264
00265
00266 if (t_b>1.0) t_b=1.0; if (t_b<0.0) t_b=0.0;
00267 if (t_g>1.0) t_g=1.0; if (t_g<0.0) t_g=0.0;
00268 if (t_r>1.0) t_r=1.0; if (t_r<0.0) t_r=0.0;
00269
00270 red= (png_byte)(t_r*254.0);
00271 green=(png_byte)(t_g*254.0);
00272 blue= (png_byte)(t_b*254.0);
00273 }
00274
00275
00276
00277
00278
00279
00280
00281 image[j][3*i+0]=red;
00282 image[j][3*i+1]=green;
00283 image[j][3*i+2]=blue;
00284 }
00285 }
00286
00287 char fmesh[256],desc[256];
00288 ierr = PetscSNPrintf(fmesh,255,"%s_%01i.%04i.%i.png",gdata->simname,2,gdata->inp,mcomp);CHKERRQ(ierr);
00289 ierr = PetscSNPrintf(desc,255,"t=%g ns, Etot=%g (J/m^3)",gdata->time*gdata->tscale*1e9,gdata->Etot*gdata->escale);CHKERRQ(ierr);
00290 WritePNGfile(fmesh,pix[0],pix[1],image,desc);
00291 }
00292
00293 VecRestoreArray(VIpol,&a);
00294
00295 ierr = PetscFree(image[0]);CHKERRQ(ierr);
00296 ierr = PetscFree(image);CHKERRQ(ierr);
00297
00298
00299 MagparFunctionInfoReturn(0);
00300 }
00301