Aide-memoire: 2017

2017年12月11日月曜日

人間の耳で音の位相は聞き取れるのか？

「人間の耳で音の位相は聞き取れるのか？」問題のテスト。音痴な自分でも結構分かる気はする...（特にphase=0と0.5の差は相当顕著）。

import ddf.minim.*;    
import ddf.minim.ugens.*;    

Minim minim;
AudioOutput out;    
Oscil sin1, sin2;
float f=440;
float A=1.0;
float trig=0.5;
void setup() {

  size(512, 200);
  smooth();
  minim = new Minim(this);
  out = minim.getLineOut(Minim.MONO, 1024);

  sin1=new Oscil(f, A, Waves.SINE);
  sin2=new Oscil(f*2, A*0.5, Waves.SINE);
  sin1.setPhase(0);
  sin2.setPhase(0);

  sin1.patch(out);
  sin2.patch(out);
}


void draw() {
  background(255);
  stroke(0);
  float x=0, Amp_old=999;
  int flag=0;
  for (int i = 0; i < out.bufferSize()-1; i++) {
    float Amp=out.left.get(i);
    if ((Amp_old=trig)||flag==1) {
      point(x, height/2 + Amp*height/4);
      flag=1;
      x+=2;
    }
    Amp_old=Amp;
  }
}

void mouseMoved() {
  float phase=map(mouseX, 0, width, 0, 1);
  sin2.setPhase(phase);
}


void stop() {
  out.close();
  minim.stop();
  super.stop();
}

2017年12月7日木曜日

Pyorbitalで方位仰角、Rangeの計算

Pyorbitalで方位仰角、Rangeの計算。こんなのがさっくりできるのは、とてもありがたい...。どのくらい信用できるのかはわからないけど。

from pyorbital import orbital　＃これを追加すると動く（外すと動かない）のは何故
from pyorbital import orbital
import pyorbital
import numpy
import datetime

#Sat info
satname='ITUPSAT 1'
f=437.325
tlef='TLE.txt'

#Position Info
#静岡大学
lon=138.431704
lat=34.965720
alt=90

#Obs Info
start_jst=datetime.datetime(2018,02,16,15,25,00)
end_jst=datetime.datetime(2018,02,16,15,45,00)

#Conf
c=299792.458 # km/s
dt=datetime.timedelta(0,30)

def jst2utc(time) :
    return time-datetime.timedelta(0,9.*60*60,0)

def obsf_ICR6(obsf) :
    f0=obsf*1000.
    f1=int(obsf*100.)*10
    df=f0-f1
    if df <2.5 :
        return f1/1000.
    elif df <7.5 : 
        return (f1+5.)/1000. 
    else :
        return (f1+10.)/1000. 
    
tle=pyorbital.tlefile.read(satname,tlef)
orb=pyorbital.orbital.Orbital(satname,tlef)

print "TLE Epoch:",",",str(datetime.datetime(2000+int(tle.epoch_year),1,1)+datetime.timedelta(days=tle.epoch_day)).split(' ')[0]
print "Date:",",",str(start_jst).split(' ')[0]
print "Lon:",",",lon
print "Lat:",",",lat
print "Alt:",",",alt
print ""

jst=start_jst
while jst < end_jst :

    now=jst2utc(jst)
    obs_pos=numpy.array(pyorbital.astronomy.observer_position(now,lon,lat,alt)[0])
    sat_pos,sat_vel=orb.get_position(now,False)

    range_dis=numpy.linalg.norm(sat_pos-obs_pos)
    range_vec=(sat_pos-obs_pos)/range_dis
    rang_rate=numpy.dot(sat_vel,range_vec)

    Az,El=orb.get_observer_look(now,lon,lat,alt)
    obsf=f*(1-rang_rate/c)
    
    print str(jst).split(' ')[1],",",Az,",",El,",",obsf_ICR6(obsf),",",obsf
    jst+=dt

2017年6月26日月曜日

大学入試問題をモンテカルロで解く

情報処理の課題として、「大学入試問題をモンテカルロで解く」というのは面白いかと思い、こちらの問題について作成。しかし、これなら素直に数え上げた方が遥かに速い（汗）。

int a,b,c;
float n1=0;
float n2=0;
void setup(){
  size(1000,20);
}
void draw(){
   a=int(random(1,7));
   b=int(random(1,7));
   c=int(random(1,7));
   
   if(log(a+b)/log(0.25)>log(c)/log(0.5)) n1+=1;
   if((pow(2,a)+pow(2,b)+pow(2,c))%3==0) n2+=1;
   background(255);
   fill(0);
   text(str(frameCount),width*0.1,height/2);
   text(str(n1/frameCount),width*0.5,height/2);
   text(str(n2/frameCount),width*0.8,height/2);
}

数えた。コンピュータ便利。

int a, b, c;
float n1=0;
float n2=0;
float n=0;
size(1000, 20);
for (a=1; a<7; a+=1) {
  for (b=1; b<7; b+=1) {
    for (c=1; c<7; c+=1) {
      if (log(a+b)/log(0.25)>log(c)/log(0.5)) n1+=1;
      if ((pow(2, a)+pow(2, b)+pow(2, c))%3==0) n2+=1;
      n+=1;
    }
  }
}
background(255);
fill(0);
text(str(n), width*0.1, height/2);
text(str(n1/n), width*0.5, height/2);
text(str(n2/n), width*0.8, height/2);

2017年5月2日火曜日

画像データの値を取得

こちらを参考に作成。面白い。

PImage myPhoto;
int myPhotoWidth, myPhotoHeight;

void setup() {
  size(500,750);
  //Thanks for http://gahag.net/004881-high-school-girl/
  myPhoto = loadImage("Test.jpg"); 
  myPhotoWidth = myPhoto.width;
  myPhotoHeight = myPhoto.height;
  noStroke();
  smooth();
  background(0);
  myPhoto.loadPixels();
}

void draw() {
  background(0,0,100);

  int x;
  int y; 
  for(x=0;x<myPhoto.width;x+=5){
  for(y=0;y<myPhoto.height;y+=5){
  
  strokeWeight(random(1,3));
  float f=random(0.5,1);
  stroke(red(myPhoto.pixels[y*width + x])*f,green(myPhoto.pixels[y*width + x])*f,blue(myPhoto.pixels[y*width + x])*f);

  float r0=random(5, 15);
  float r1=random(5, 15);
  float th0=random(0, 2*PI);
  float th1=random(0, 2*PI);
  float dx=random(-5, 5);
  float dy=random(-5, 5);

  line(r0*cos(th0)+x+dx, r0*sin(th0)+y+dy, r1*cos(th1)+x+dx, r1*sin(th1)+y+dy);
   }}
}

2017年4月30日日曜日

色々統合した練習

もう少し賢く書けそう。こちら。

PImage img;

int R=200, G=200, B=200;

void setup() {
  size( 800, 800 );
  background( R, G, B );
  //Thanks for Yume Yume Iro TOWN http://www.poipoi.com/yakko/cgi-bin/sb/log/eid1737.html
  img = loadImage("siruetto20080205a.png");
  stroke(255, 0.0);
  noCursor();
}
void draw() {
  background(R, G, B);
  float x=mouseX-280;
  float y=mouseY-500;


  for (int i=0; i<10; i++) {
    fill(R, G, B, 30);
    noStroke();
    rect(0, 0, width, height);
    float xi=map(i, 0, 10, width/2, x);
    float yi=map(i, 0, 10, height/2, y);
    float f=map(i, 0, 10, 0.1, 1);
    image(img, xi, yi, 560*f, 1000*f);
  }
}
void mouseMoved() {
  int a=int(random(1, 4));
  int pm;
  if (R+G+B<50) {
    pm=1;
  } else if (random(0, 1.0)<0.5) {
    pm=-1;
  } else {
    pm=1;
  }

  if (a==1) {
    if (R==200) {
      R=199;
    } else if (R==1) {
      R=2;
    } else {
      R=R+pm;
    }
  }
  if (a==2) {
    if (G==200) {
      G=199;
    } else if (G==1) {
      G=2;
    } else {
      G=G+pm;
    }
  } 
  if (a==3) {
     if (B==200) {
      B=199;
    } else if (B==1) {
      B=2;
    } else {
      B=B+pm;
    }
  } 
}

2017年4月29日土曜日

マウスの使い方練習

マウスの使い方練習。IRHH 2

PImage img;

void setup() {
  size( 600, 600 );
  background( 0, 0, 255 );
  //Thanks for Yume Yume Iro TOWN http://www.poipoi.com/yakko/cgi-bin/sb/log/eid1737.html
  img = loadImage("siruetto20080205a.png");
  stroke(255, 0.0);
}
void draw() {
  background(0, 0, 255);
  image(img, mouseX-140, mouseY-250, 280, 500);
  //image(img, mouseX-120, mouseY-80, 280, 500);

  beginShape();
  for (int i=0; i<300; i++) {
    stroke(random(0, 255), 0, 0);
    strokeWeight( random(1, 3) );

    float r0=random(0, 50);
    float r1=random(0, 50);
    float th0=random(0, 2*PI);
    float th1=random(0, 2*PI);
    line(r0*cos(th0)+mouseX-25, r0*sin(th0)+mouseY-170, r1*cos(th1)+mouseX-25, r1*sin(th1)+mouseY-170);
  }
}

ペジエの使い方の練習

ペジエの使い方の練習。こちら

PImage img;


void setup() {
  size( 800, 800 );
  background( 255 );
  //Thanks for Yume Yume Iro TOWN http://www.poipoi.com/yakko/cgi-bin/sb/log/eid1737.html
  img = loadImage("siruetto20080205a.png");
  strokeWeight( 6 );
}

void draw() {
  if (frameCount%3==0) {
    float f=random(0.1, 0.8);
    image(img, random(0, 500), random(0, 500), 560*f, 1000*f);
  }
  stroke(  random(255), random(255), random(255) );
  fill( random(255), random(255), random(255));


  beginShape();
  int x0=radomxy();
  int y0=radomxy();

  vertex( x0, y0 );
  bezierVertex( radomxy(), radomxy(), radomxy(), radomxy(), radomxy(), radomxy() );
  int i=0;
  while (i==0) {
    if (random(6)>3) {
      i=1;
    }
    bezierVertex( radomxy(), radomxy(), radomxy(), radomxy(), radomxy(), radomxy() );
  }
  bezierVertex( radomxy(), radomxy(), radomxy(), radomxy(), x0, y0 );
  endShape();
}

void mousePressed() {
  //background( 255 );
  delay(1000);
}

int radomxy() {
  return int(random(-100, height+100));
}

2017年4月26日水曜日

なにかとむなしい

こちら。

int x, y, vx, vy;

void setup() {
  smooth();
  size(800, 800);
  x=width/2;
  y=height/2;
  background(255);
  int v=20;
  vx=int(random(-1*v, v));
  int a;
  if (random(1)>0.5) {
    a=1;
  } else {
    a=-1;
  }

  vy=a*int(sqrt(pow(v, 2)-pow(vx, 2)));
}

void draw() {
  float a=random(1, 4);
  int xnew=x+int(a*vx);
  int ynew=y+int(a*vy);
  if (xnew<0) {
    xnew=0;
    vx=-1*vx;
    vy=int(vy*random(0.9, 1.1));
  }
  if (xnew>width) {
    xnew=width;
    vx=-1*vx;
    vy=int(vy*random(0.9, 1.1));
  }
  if (ynew<0) {
    ynew=0;
    vy=-1*vy;
    vx=int(vx*random(0.9, 1.1));
  }
  if (ynew>height) {
    ynew=height;
    vy=-1*vy;
    vx=int(vx*random(0.8, 1.2));
  }
  if (vx*vx<30) {
    vx=int(random(-20, 20));
  }
  if (vy*vy<30) {
    vy=int(random(-20, 20));
  }  
  stroke(random(255), random(255), random(255));
  strokeWeight(random(30));
  line(x, y, xnew, ynew);
  x=xnew;
  y=ynew;
  if (frameCount%10==0) {
    fill(255, 20);
    noStroke();
    rect(0, 0, width, height);
  }
}

2017年4月19日水曜日

月の自転・公転周期デモ

月の自転・公転周期の説明デモ。動いているのはこちら。意図せずもドラoもんカラー。

float theta=0;
int n=0;

void setup()
{
  size(1000, 1000);
  background(255);
  noLoop();
}

void draw()
{
  float r=width/3;
  float re=width/8;
  float rm=width/10;

  background(255);
  stroke(100);
  fill(0, 0, 255);
  ellipse(width/2, height/2, re, re);
  stroke(0);

  fill(255, 255, 0);

  ellipse(width/2+r*cos(theta), height/2+r*sin(theta), rm, rm);
  fill(255, 0, 0);

  ellipse(width/2+r*cos(theta)+0.35*rm*cos(PI+1*theta), height/2+r*sin(theta)+0.35*rm*sin(PI+1*theta), rm*0.2, rm*0.2);
  theta+=PI/200;
}

void mousePressed() {
  if (n==0) {
    loop();
    n=1;
  } else {
    noLoop();
    n=0;
  }
}

2017年4月14日金曜日

マル9円周率その3

意味無しモンテカルロ。こちら。

float n=1;
float s=0;

void setup() {
  size(300, 350);
  background(255);  
  noFill();
  stroke(0, 2555, 0);
  ellipse(0, 0, 600, 600);
}

void draw() {

  float x=random(0.0, 1.0);
  float y=random(0.0, 1.0);
  if (pow(x, 2)+pow(y, 2)<1.0) {
    stroke(255, 0, 0);
    s+=1;
  } else {
    stroke(0, 0, 255);
  }
  ellipse(x*300, y*300, 1, 1);
  stroke(255);
  fill(255);
  rect(0, 310, 300, 50);
  fill(0);
  text(str(int(n)), 5, 330);  
  text(str(s/n*4.0), 150, 330);
  n++;
}

void mousePressed() {
  n=1;
  s=0;
  background(255);
  noFill();
  stroke(0, 2555, 0);
  ellipse(0, 0, 600, 600);
}

マル9な円周率その2

その2。動いているのはこちら。

float n=1;

void setup(){
  size(300,350);
  background(255);  
}

void draw(){
  float s=0;
  background(255);
  fill(255);
  noFill();
  stroke(255,0,0);
  ellipse(0,0,600,600);
  float x0=0;
  float y0=1.0;


  for(int m=1;m<n+1;m++){
    float x=m/n;
    float y=sqrt(1-pow((x),2.0));
    stroke(0);
    line(x*300,y*300,x0*300,y0*300);
    s=s+sqrt(pow((x-x0),2.0)+pow((y-y0),2.0));
    x0=x;
    y0=y;
}
  n+=int(n/10+1);
 if(n<300){
   delay(int(1000/sqrt(n)));
  }
  fill(0);
  textSize(15);
  text(str(int(n)),5,330);  
  text(str(s*2.0),150,330);
  if(n>10000000){
    noLoop();
  }    
}

void mousePressed(){
  n=1;
  loop();
}

2017年4月13日木曜日

マル9な円周率計算デモ

頭の悪い円周率計算デモ。Processing.jsで動かすとこちら。

float n=1;

void setup(){
  size(300,350);
  background(255);  
}

void draw(){
  float s=0;
  background(255);
  fill(255);
  for(int m=0;m<n;m++){
    float x=m/n;
    float f=sqrt(1-pow((x),2.0));
    stroke(0);
    rect(x*300,0,1/n*300,f*300);
    noFill();
    stroke(255,0,0);
    ellipse(0,0,600,600);

    s=s+f*1/n;
  }
  n+=int(n/10+1);
 if(n<300){
    delay(int(1000/sqrt(n)));
  }
  fill(0);
  textSize(15);
  text(str(int(n)),5,330);  
  text(str(s*4.0),150,330);
  if(n>10000000){
    noLoop();
  }    
}

void mousePressed(){
  n=1;
  loop();
}

2017年4月7日金曜日

地球の自転に対するアリストテレス or ニュートン的描像

天文学の授業で「地動説に対する反論」を説明する為にProcessingで作成したアニメーション。こういう教材を簡単に作れるので、本当はどこかの授業で学生さんにProcessingを触ってもらいたいところ。Processing.jsで動かしたアリストテレス or ニュートン的な描像はそれぞれ、こちらとこちら。

float glevel=0.8;

float x;

float xb, yb, vy=0, vxw;

float headsize;

int t=0;

int ball_start=150;

float g;

float[] xbo=new float[500];

float[] ybo=new float[500];



void setup() {

  x=100;

  size(1100, 700);

  colorMode(RGB, 256);

  background(255, 255, 255);

  headsize=width*0.03;

  frameRate(50);

  xb=0;

  yb=height*glevel*0.15;

  g=0.05;

  vxw=width*0.002;



  for (int i=0; i<500 i="" p="">    xbo[i]=0;

    ybo[i]=0;

  }

}

void draw() {

  background(255, 255, 255);

  //smooth() ;

  //World

  stroke(0, 0, 0);

  strokeWeight(1);

  line(0, height*glevel, width, height*glevel);

  noStroke();

  //point

  fill(255, 0, 0);

  ellipse(x, height*glevel, width*0.01, width*0.01);

  //man

  fill(0, 0, 0);



  //head

  ellipse(x+width*0.1, height*glevel*0.71, headsize, headsize);

  //body

  rect(x+width*(0.1-0.025), height*glevel*0.75, width*(0.05), height*0.1);

  //leg

  rect(x+width*(0.1-0.025), height*glevel*0.89, width*(0.02), height*0.085);

  rect(x+width*(0.1+0.005), height*glevel*0.89, width*(0.02), height*0.085);

  //arm

  rect(x+width*(0.1-0.045), height*glevel*0.75, width*(0.015), height*0.08);

  rect(x+width*(0.1+0.03), height*glevel*0.75, width*(0.015), height*0.08);



  //Bar

  strokeWeight(10);

  stroke(183, 65, 14);

  line(x+width*0.2, height*glevel, x+width*0.2, height*glevel*0.1);

  line(x, height*glevel*0.1, x+width*0.2, height*glevel*0.1);



  //Ball Shadow

  for (int i=1; i    if ( i%20==0) {

      stroke(150);

      strokeWeight(1);

      noFill();

      ellipse(xbo[i], ybo[i], width*0.03, width*0.03);

    }

  }



  //Ball

  xbo[t]=xb;

  ybo[t]=yb;

  if (t>ball_start & yb<=height*glevel) {

    vy=vy+g;

    yb=yb+vy;

    //下の1文を入れるか入れないかが両者の差。

    xb=x;

  } else {

    xb=x;

  }



  if (yb>=height*glevel) {

    vy=0;

    yb=height*glevel;

    xb=xb+vxw;

  }





  noStroke();

  fill(0, 0, 255);

  ellipse(xb, yb, width*0.03, width*0.03);

  x=x+vxw;

  print(t, "\n");

  if (t<499 p="">    t++;

  }

}

2017年4月5日水曜日

Processingで矩形波のフーリエ変換の可視・音声化デモ

こちらを参考にさせていただきました。学生さんの卒論に向けたProcessingでのフーリエ合成（というかMinimの使い方）練習。Processing、本当に楽しい...。

import ddf.minim.*;   
import ddf.minim.signals.*;   


Minim minim;
AudioOutput out;  
SineWave[] sine;
float waveH = 50;  
int n=30;



void setup()
{

  size(512, 200);
  smooth();
  minim = new Minim(this);
  sine=new SineWave[n];

  out = minim.getLineOut(Minim.STEREO);

  int i;

  for(i=0;i < n;i++){
    int m=i+1;
    sine[i] = new SineWave(440*m, 2*(1-cos(m*PI))/(m*PI*440)*100, out.sampleRate());
  }



  for(i=0;i < n;i++){
  out.addSignal(sine[i]);
  }

}



void draw()
{
  background(255);
  stroke(0);
  for(int i = 0; i < out.bufferSize()-1; i++)

  {

    point(i, 50 + out.left.get(i)*waveH); 
    point(i, 150 + out.right.get(i)*waveH);

  }

}

void mouseMoved() {
  int t = int(map(mouseX, 0, width, 1, n));

  for(int i=0;i < t;i++){
    int m=i+1;
    sine[i].setAmp(2*(1-cos(m*PI))/(m*PI*440.)*100); 

  }
  for(int i=t;i < n;i++){
    sine[i].setAmp(0); 
  }
}



void stop()
{
  out.close();
  minim.stop();
  super.stop();
}

登録: 投稿 (Atom)

Aide-memoire

2017年12月11日月曜日

人間の耳で音の位相は聞き取れるのか？

2017年12月7日木曜日

Pyorbitalで方位仰角、Rangeの計算

2017年6月26日月曜日

大学入試問題をモンテカルロで解く

2017年5月2日火曜日

画像データの値を取得

2017年4月30日日曜日

色々統合した練習

2017年4月29日土曜日

マウスの使い方練習

ペジエの使い方の練習

2017年4月26日水曜日

なにかとむなしい

2017年4月19日水曜日

月の自転・公転周期デモ

2017年4月14日金曜日

マル9円周率その3

マル9な円周率その2

2017年4月13日木曜日

マル9な円周率計算デモ

2017年4月7日金曜日

地球の自転に対するアリストテレス or ニュートン的描像

2017年4月5日水曜日

Processingで矩形波のフーリエ変換の可視・音声化デモ

プロフィール

ブログアーカイブ

ラベル

このブログを検索

2017年12月11日月曜日

2017年12月7日木曜日

2017年6月26日月曜日

2017年5月2日火曜日

2017年4月30日日曜日

2017年4月29日土曜日

2017年4月26日水曜日

2017年4月19日水曜日

2017年4月14日金曜日

2017年4月13日木曜日

2017年4月7日金曜日

2017年4月5日水曜日

プロフィール

ブログ アーカイブ

ラベル

このブログを検索

ブログアーカイブ